The code below shall be used to calculate numerically the integral of f(x), using rectangles. The while loop ' error> er_user_max' seems not to evaluate the logical operation properly and changes the row dimension of matrix 'save_area'.

6 views (last 30 days)
clear all
%define function
f= @(x)x.^2;
%length of function domain
L=10;
%initial integration step
dx=1;
%introducing a parameter for the total sum underneath the
%curve.
sum_area=0;
%max refinment passes within one dx element
ref_max=10;
%define the allowed error, calculated by the equidistantly
%split rectangle
%in two parts to calculate each area of both split parts
%separatly.
% Used formular is written below @error calculation
%For comparision, the two parts are summed up and divided by the
%former,
%unsplit rectangle.
er_user_min=0.995;
er_user_max=1.2499999;
%____________________________________
for i=1:1 %for each dx
a=(i-1)*dx;
b=i*dx;
%introduce a running number for counting the partitions
%within a
%initial bar element with the length of dx.
pn=1;
rn=pn+1;
%max refinment passes within one dx element
%ref_max=10;
%initialize storage matrix of created partitions
save_area= zeros(ref_max,3);
%introducing a parameter for the partial areas within a dx
%element.
area_part= 0;
%introducing a parameter for the sum of rectangular pieces
%within dx
area_dx=0;
%calculate initial area of specific dx element.
area=((b-a))*f((i-1)*dx+(b-a)/2);
%intial error calculation
b_new=(i-1)*dx+(b-a)/2;
area_l = ((b_new-a))*f((i-1)*dx+(b_new-a)/2);
area_r = ((b-b_new))*f(b_new+(b-b_new)/2);
error = (area_l+area_r)/area;
%save area of dx to add to sum_area, if error< er_user.
area_comp= area_l+area_r;
%calculate the total area of all partitions within one
%initial bar
%element with the lengt of dx, with respect to the error
%evaluation.
% while rn>1
%split dx in partitions to grant a certain grade of accuracy
%due to
%limitations of maximum possible refinment passes.
while error> er_user_max
%calculate new boundaries of the left side rectangle.
b_new=(i-1)*dx+(b-a)/2;
area_l = ((b_new-a))*f((i-1)*dx+(b_new-a)/2);
area_r = ((b-b_new))*f(b_new+(b-b_new)/2);
%store boundaries and area of right side rectangle.
a_r= b_new;
b_r= b;
%save_area(pn,:)=[a_r,b_r,area_r];
save_area(pn,1)=a_r;
save_area(pn,2)=b_r;
save_area(pn,3)=area_r;
%update the latest considered area
area_part= area_l;
%Without if condition, the row dimension as defined
%of ref_max
%of matrix 'save_area' will be largly overexceeded!
%if pn < 10
%add information that one refinment has been done.
%pn= pn+1;
%rn=rn+1;
%adjust new boundaries
%b= b_new;
%a=a;
%else
% break
%end
%calculate new error
b_new=(i-1)*dx+(b-a)/2;
area_l = ((b_new-a))*f((i-1)*dx+(b_new-a)/2);
area_r = ((b-b_new))*f(b_new+(b-b_new)/2);
area_comp=area_l+area_r;
error = area_comp/area_part;
end
%add error evaluated left-side part to the sum of the
%area of dx.
area_dx=area_dx+ area_comp;
%...further coding, does not display the mentioned
%problem with the coding%
%...continuing relevant coding
%Sum up the calculated and error evaluated partial rectangles
%underneath
%the curve
%sum_area=sum_area+area_dx;
end
disp(sum_area)
disp(save_area)

Answers (1)

dpb
dpb on 19 Aug 2015
Way too much code instead of concentrating on the concern and put the question in the body of the text and make the title just that...but, that aside--
...
error = (area_l+area_r)/area;
...
%limitations of maximum possible refinment passes.
while error> er_user_max
...
As defined your error term is a ratio, not a relative nor absolute and the test as written requires the estimates always be high in comparison to the reference to fail; if the initial guess is low instead, nothing will happen.
Probably you intend something like
error = 1-(area_l+area_r)/area;
and
while abs(error)>er_user_max
Also note that error is a Matlab builtin function; suggest using some other name for the variable ( err or similar isn't taken).

Categories

Find more on Loops and Conditional Statements in Help Center and File Exchange

Products

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!