MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn moreOpportunities for recent engineering grads.

Apply Today**New to MATLAB?**

Asked by Tim
on 28 Aug 2014

Hi,

I am trying to implement Romberg's integration process in MAT LAB (i.e. recursive-trapezoidal to simpson's 1/3 rule(or Richardsons's extrapolation) to finally Boole's rule. After computing my answers for the trapezoidal rule, I then call on the results from my trapezoidal code to then apply Romberg intergration for Simpson's 1/3 rule. And then call on the results from the Rom berg simpson code to finally apply Romberg integration using Boole's rule.

My question is simply around the inability to call the results effectively. My code for Romberg Simpson's 1/3 is able to call on the trapezoidal results correctly, but for some reason my Romberg boole code keeps calling the trapzoidal data instead of the newly extrapolated Romberg Simpson's 1/3 data and applies the equation. And I dont understand why?

TRAPEZOIDAL

n=1; %first value of area (i.e. trapezoidal rule) result(n)=((b-a)./2*(func(b)-func(a))) %start error error_result=1; %while loop that creates error bounds. Will cycle until error is larger %than toler while error_result(n)>toler; it=2.^(n-1); tnm=it; h=(b-a)./tnm; %first x values is from startig from a +0.5 h interval due to split x=a+0.5.*h; %first sum=0 sum=0.0; for j=1:it %sum is sum of sum=sum+feval(func,x); %finding next x value x=x+h; end %Answer/Result of integration at given interval split result(n+1)=0.5.*(result(n)+h.*sum); %error analysis for each interval split error_result(n+1)=abs((result(n+1)-result(n))); %counter n=n+1; end

ROMBERG SIMPSON

function[result]=Rombergsimp(func,a,b,toler) %starting counter n=1; %error starter error_result=1; while error_result(n)>toler; %calling result values from trapezoidal code result = trapz(func,a,b,toler); %defining new result values by applying simpsons 1/3 romberg integration. result_romsim(n+1)=result(n+1)+1/3.*(result(n+1)-result(n)); %error analysis of each step error_result(n+1)=abs(((result_romsim(n+1)-result_romsim(n)))); %counter+1 n=n+1; %printing results end

ROMBERG BOOLE

function[result]=Rombergboole(func,a,b,toler) %counter n=1; %error starter error_result=1; %while lopp defining error analysis while error_result(n)>toler; %calling result values from Rombergsimp m-function result = Rombergsimp(func,a,b,toler) %Applying Boole's rule to the already extrapolated simpson data result_romboole(n+1)=result(n)+16/15.*(result(n+1)-result(n)); %displaying results disp(result_romboole); %error check between each new Romber_boole value. error_result(n+1)=abs(((result_romboole(n+1)-result_romboole(n)))); %counter+1 n=n+1; end

I am sorry for the lengthiness of the code and question but I didn't know how else to ask it.

thanks

*No products are associated with this question.*

Answer by Star Strider
on 28 Aug 2014

I don’t completely understand how you’re calling your routines, but it seems if you call `‘Rombergboole’` alone, it calls all the others, first `‘Rombergsimp’` that itself calls `‘trapz’` and returns the eventual result.

I don’t understand what `‘func’` is. The ** trapz** function takes a vector (or array) argument — not a function argument — and produces a scalar (or vector) result. There would seem to me to be nothing left to integrate. In other places, it seems that you are calling

Tim
on 28 Aug 2014

'func' is the specified function that we were asked to integrate. My lecturer explained it to me that;

result=trapz(func,a,b,toler)

simply calls all previous answers that trapz came up with into a vector, in which the Rombergsimp code is now able to apply its formula to. The exact same process should be done for Rombergboole, this time calling Rombergsimp values.

I forgot to add the function part for trapz at the very top it is;

function [result]=trapz(func,a,b,toler) %Start counter n=1; %first value of area (i.e. trapezoidal rule)

and ill add the first trapezoidal m-file with the data;

%ID:25118900 DATE:26/8/2014 %RECURSIVE TRAPEZOIDAL RULE FOR INTERGRATION-DATA FILE %==================================================== %Function handle of current equation where x=t func=@(x) (10.*exp((-2.5.*x)./1).*sin((3.5.*pi.*x))./1).^2 %starting intergrad bound a=0; %end bound b=0.5; %tolerance toler=10^-4; %Linkage to m-function. [result]=trapz(func,a,b,toler);

and unfortunately we can just use quad to integrate as they don't let us.

Star Strider
on 28 Aug 2014

The problem is that:

result=trapz(func,a,b,toler)

isn’t going to work. See the documentation on ** trapz**.

Creating your own function and calling it ** trapz** is quite definitely bad

The added script lines *do* help me understand your code.

Still, your description of how you want your code to work doesn’t match what it appears to me to do. (I couldn’t run it without the code you just now added.) Your (unfortunately named) `‘trapz’` function doesn’t seem to me to call any other functions. However your `‘Rombergboole’` function calls all of them.

I’ll work with your code in the morning. It’s ** GMT-6** here, and the end of my day.

## 0 Comments