اطلاعیه

**manian** · 2009-01-21T22:23:17

پاسخ : انتگرال به روش ذوزنقه ای

با سلام
بهتر است برای یادگیری دروس، تکالیف دانشگاه را خودتان انجام دهید.

**jafarzadeh85** · 2009-01-22T20:06:52

پاسخ : انتگرال به روش ذوزنقه ای

کاربر گرامی manian
شما که در حد اکابر درس خوندی چرا در مورد چیزی که خبر نداری اینطوری نظر می دی؟؟؟
من درسم تموم شده، اینم پروژه یکی از دوستان منه که فوریم هست.

**manian** · 2009-01-22T20:18:19

پاسخ : انتگرال به روش ذوزنقه ای

نوشته اصلی توسط jafarzadeh85

کاربر گرامی manian
شما که در حد اکابر درس خوندی چرا در مورد چیزی که خبر نداری اینطوری نظر می دی؟؟؟
من درسم تموم شده، اینم پروژه یکی از دوستان منه که فوریم هست.

با سلام
شرمنده من فکر کردم که شما تکالیف دانشگاهی خودتان را میآ‌خواهید انجام دهید نه برای دوستانتان!
من سوادم اون قدر نیست به خصوص که الان فهمیدم در اکابر به ما به اشتباه یاد دادهآ‌اند که Matlab به فارسی متلب نوشته میآ‌شود نه مطلب!

**رضا سپاس یار** · 2009-01-22T20:35:34

پاسخ : انتگرال به روش ذوزنقه ای

سلام،

من قدیما یه سری M. فایل برای محاسبه ی انتگرال به روش های عددی آماده کرده بودم.

روش ذوزنقه ای:

کد:

function z = trapint(x, y)
%
% ------------------------------------------------------
%   Reza Sepas Yar
%   http://www.avr.ir
% ------------------------------------------------------
%
% Trapezoidal integration of (x,y) data
% UNEQUALLY SPACED DATA
%
% USAGE: z = trapint(x,y)
%
% input: x = vector of x values
%     y = vector of y values
%
% output: z = value of integral

n = length(x);

isum = 0;

for i = 1:n-1
 isum = isum + (y(i) + y(i+1))*(x(i+1) - x(i))/2;
end

z = isum;

مثال، محاسبه ی انتگرال 1./sqrt(1+x.^2) به روش ذوزنقه ای:

کد:

&gt;&gt; x = [0:0.25:1];
&gt;&gt; y = 1./sqrt(1+x.^2);
&gt;&gt; z = trapint(x, y)

z =

  0.8795

فایل های پیوست شده

trapint.m (353 بایت, 13 مشاهدات)

**رضا سپاس یار** · 2009-01-22T20:37:18

پاسخ : انتگرال به روش ذوزنقه ای

روش گاوس:

کد:

%
% ------------------------------------------------------
%   Reza Sepas Yar
%   http://www.avr.ir
% ------------------------------------------------------
%
% f is a function to integrating
% a, the lower limit of integration
% b, the upper limit of integration
% n, the maximum number of Gauss points, 1-10

% Example1:
% f= @(x) (x.*log(x))./(1+x);
% a=0;
% b=1;
% n=2;

% Example2:
% f = @(x) 1./sqrt(2+(cos(x)).^2);
% a=0;
% b=1;
% n=2;

% Example3
% f = @(x) sin(x)./x;
% a = 0 + 1e-10;
% b = pi/2 + 1e-10;
% n = 3;

% Example3
% f = @(x) sin(x)./x;
% a = 0 + 1e-10;
% b = pi/2 + 1e-10;
% n = 5;
 
% Example4
% f = @(x) 1./sqrt(1+(sin(x)).^2);
% a = 0;
% b = 1;
% n = 5;

% Example4
% f = @(x) (x.*log(x))./(1+x);
% a = 0;
% b = 1;
% n = 5;

% Example4
% f = @(x) sin(x).^2;
% a = 0;
% b = pi/2;
% n = 5;
 
% Direct Way: quad(f, a, b)
 
%**********************************************************************

% Displays what inputs are being used
disp(sprintf(&#039;\n\n********************************Input Data**********************************\n&#039;))
disp(sprintf(&#039;  f(x), integrand function&#039;)) 
disp(sprintf(&#039;  a = %g, lower limit of integration &#039;,a)) 
disp(sprintf(&#039;  b = %g, upper limit of integration &#039;,b)) 
disp(sprintf(&#039;  n = %g, number of Gauss Points, 1-10&#039;,n)) 
format long

% Setup Gaussian Coefficients and Evaluation Points
x_values = zeros(10,10) ;
c_values = zeros(10,10) ;

i=1 ;
x_values(i,1) = 0.0 ;
c_values(i,1) = 2.0 ;

i=2 ;
x_values(i,1) = -0.5773502691896260 ; 
x_values(i,2) = -x_values(i,1) ;
c_values(i,1) = 1.0 ;
c_values(i,2) = 1.0 ;

i=3 ;
x_values(i,1) = -0.7745966692414830 ;
x_values(i,2) = 0.0 ;
x_values(i,3) = -x_values(i,1) ;
c_values(i,1) = 0.5555555555555560 ;
c_values(i,2) = 0.8888888888888890 ;
c_values(i,3) =c_values(i,1) ;

i=4 ;
x_values(i,1) = -0.8611363115940530 ;
x_values(i,2) = -0.3399810435848560 ;
x_values(i,3) = -x_values(i,2) ;
x_values(i,4) = -x_values(i,1) ;
c_values(i,1) = 0.3478548451374540 ;
c_values(i,2) = 0.6521451548625460 ;
c_values(i,3) =c_values(i,2) ;
c_values(i,4) =c_values(i,1) ;

i=5 ;
x_values(i,1) = 0.9061798459386640 ;
x_values(i,2) = 0.5384693101056830 ;
x_values(i,3) = 0.0 ;
x_values(i,4) = -x_values(i,2) ;
x_values(i,5) = -x_values(i,1) ;
c_values(i,1) = 0.2369268850561890 ;
c_values(i,2) = 0.4786386704993660 ;
c_values(i,3) = 0.5688888888890000 ;
c_values(i,4) =c_values(i,2) ;
c_values(i,5) =c_values(i,1) ;

i=6 ;
x_values(i,1) = -.9324695142032520 ;
x_values(i,2) = -.6612093864662650 ;
x_values(i,3) = -.2386191860831970 ;
x_values(i,4) = -x_values(i,3) ;
x_values(i,5) = -x_values(i,2) ;
x_values(i,6) = -x_values(i,1) ;
c_values(i,1) = 0.1713244923791700 ;
c_values(i,2) = 0.3607615730481390 ;
c_values(i,3) = 0.4679139345726910 ;
c_values(i,4) =c_values(i,3) ;
c_values(i,5) =c_values(i,2) ;
c_values(i,6) =c_values(i,1) ;

i=7 ;
x_values(i,1) = -0.9491079123427590 ;
x_values(i,2) = -0.7415311855993940 ;
x_values(i,3) = -0.4058451513773970 ;
x_values(i,4) = 0.0 ;
x_values(i,5) = -x_values(i,3) ;
x_values(i,6) = -x_values(i,2) ;
x_values(i,7) = -x_values(i,1) ;
c_values(i,1) = 0.1294849661688700 ;
c_values(i,2) = 0.2797053914892770 ;
c_values(i,3) = 0.3818300505051190 ;
c_values(i,4) = 0.4179591836734690 ;
c_values(i,5) =c_values(i,3) ;
c_values(i,6) =c_values(i,2) ;
c_values(i,7) =c_values(i,1) ;

i=8 ;
x_values(i,1) = -0.9602898564975360 ;
x_values(i,2) = -0.7966664774136270 ;
x_values(i,3) = -0.5255324099163290 ;
x_values(i,4) = -0.1834346424956500 ;
x_values(i,5) = -x_values(i,4) ;
x_values(i,6) = -x_values(i,3) ;
x_values(i,7) = -x_values(i,2) ;
x_values(i,8) = -x_values(i,1) ;
c_values(i,1) = 0.1012285362903760 ;
c_values(i,2) = 0.2223810344533740 ;
c_values(i,3) = 0.3137066458778870 ;
c_values(i,4) = 0.3626837833783620 ;
c_values(i,5) =c_values(i,4) ;
c_values(i,6) =c_values(i,3) ;
c_values(i,7) =c_values(i,2) ;
c_values(i,8) =c_values(i,1) ;

i=9 ;
x_values(i,1) = -0.9681602395076260 ;
x_values(i,2) = -0.8360311073266360 ;
x_values(i,4) = -0.6133714327005900 ;
x_values(i,4) = -0.3242534234038090 ;
x_values(i,5) = 0.0 ;
x_values(i,6) = -x_values(i,4) ;
x_values(i,7) = -x_values(i,3) ;
x_values(i,8) = -x_values(i,2) ;
x_values(i,9) = -x_values(i,1) ;
c_values(i,1) = 0.0812743883615740 ;
c_values(i,2) = 0.1806481606948570 ;
c_values(i,3) = 0.2606106964029350 ;
c_values(i,4) = 0.3123470770400030 ;
c_values(i,5) = 0.3302393550012600 ;
c_values(i,6) =c_values(i,4) ;
c_values(i,7) =c_values(i,3) ;
c_values(i,8) =c_values(i,2) ;
c_values(i,9) =c_values(i,1) ;

i=10 ;
x_values(i,1) = -0.9739065285171720 ;
x_values(i,2) = -0.8650633666889850 ;
x_values(i,3) = -0.6794095682990240 ;
x_values(i,4) = -0.4333953941292470 ;
x_values(i,5) = -0.1488743389816310 ;
x_values(i,6) = -x_values(i,5) ;
x_values(i,7) = -x_values(i,4) ;
x_values(i,8) = -x_values(i,3) ;
x_values(i,9) = -x_values(i,2) ;
x_values(i,10) = -x_values(i,1) ;
c_values(i,1) = 0.0666713443086880 ;
c_values(i,2) = 0.1494513491505810 ;
c_values(i,3) = 0.2190863625159820 ;
c_values(i,4) = 0.2692667193099960 ;
c_values(i,5) = 0.2955242247147530 ;
c_values(i,6) = c_values(i,5) ;
c_values(i,7) = c_values(i,4) ;
c_values(i,8) = c_values(i,3) ;
c_values(i,9) = c_values(i,2) ;
c_values(i,10) = c_values(i,1) ;


% simulation of the steps needed to perform the method begins here
% disp(sprintf(&#039;\n*******************************Simulation*********************************\n&#039;))

% First do a lookup of the weights and evaluation locations
disp(&#039;1) Lookup Gauss point constants and evaluation points from the table.&#039;)
disp(&#039; &#039;)
disp(sprintf(&#039;     Evaluation Locations&#039;))

% Used to display the gauss points and their weights
for i=1:n
 disp(sprintf(&#039;      x%i=%1.15f&#039;,i,x_values(n,i)))
end
disp(&#039; &#039;)
disp(sprintf(&#039;     Evaluation Constants&#039;))
for i=1:n
 disp(sprintf(&#039;      c%i=%1.15f&#039;,i,c_values(n,i)))
end

% Now transform the evaluation points to fit your interval
% disp(&#039; &#039;)
% disp(&#039;2) The evaluation points are definted for the interval [-1,1]. Transform them&#039;)
% disp(&#039; to fit [a,b]. This is done by the following formula&#039;)
% disp(&#039; &#039;)
% disp(&#039;     xnew = (b-a)/2 * x + (b+a)/2&#039;)
% disp(&#039; &#039;)
% disp(sprintf(&#039;     Transformed Evaluation Locations&#039;))
% Transform the locations
for i=1:n
 xv(i)=(b-a)/2*x_values(n,i)+(b+a)/2 ;
 %disp(sprintf(&#039;      x%i=%1.15f&#039;,i,xv(i)));
end

% Apply the Gaussian procedure
disp(&#039; &#039;)
disp(&#039;2) The approximation is given as&#039;)
disp(&#039; &#039;)
disp(&#039;       c(i)     *     f(x(i))&#039;)

approx = 0 ;
disp(&#039; &#039;)

% approximation via Gaussian integration
for i=1:n-1
 approx = approx + c_values(n,i)*f(xv(i)) ;
 disp(sprintf(&#039;     %1.15f * %1.15f = %1.15f&#039;,c_values(n,i),f(xv(i)), c_values(n,i)*f(xv(i) )))
end
approx = approx + c_values(n,n)*f(xv(n)) ;
disp(sprintf(&#039;    + %1.15f * %1.15f = %1.15f&#039;,c_values(n,n),f(xv(n)), c_values(n,n)*f(xv(n) )))
disp(sprintf(&#039;    ----------------------------&#039;))
disp(sprintf(&#039;     %1.15f&#039;,approx))

% the result must be scaled to fit the interval also
disp(&#039; &#039;)
disp(&#039;3) The approximation must also be transformed to fit the [a,b] interval.&#039;)
disp(&#039; &#039;)
disp(&#039; Answer =   (b-a)/2 *  sum( c(i)* f(x(i)) ) &#039;)
disp(&#039; &#039;)
disp(sprintf(&#039;    = %1.15f * %1.15f&#039;,(b-a)/2,approx))

% also transform the approximation
approx = (b-a)/2 * approx ;
disp(sprintf(&#039;        = %1.15f&#039;,approx))

% Display results section
disp(sprintf(&#039;\n\n**************************Results****************************&#039;))

% The following finds what is called the &#039;Exact&#039; solution
exact = quad(f,a,b) ; 

disp(sprintf(&#039;\n Approximate = %1.15f&#039;,approx)) 
disp(sprintf(&#039; Exact   = %1.15f&#039;,exact))

کد:

**رضا سپاس یار** · 2009-01-22T20:40:16

پاسخ : انتگرال به روش ذوزنقه ای

روش کوادراتور:

کد:


% ------------------------------------------------------
%   Reza Sepas Yar
%   http://www.avr.ir
% ------------------------------------------------------

 f= x.*log(x)./(1+x);

% a, the lower limit of integration

 a=1e-10 ;

% b, the upper limit of integration

 b=1 ;

% n, the maximum number of Gauss points, 1-10

 n=2 ;

%**********************************************************************

% Displays what inputs are being used
disp(sprintf(&#039;\n\n********************************Input Data**********************************\n&#039;))
disp(sprintf(&#039;  f(x), integrand function&#039;)) 
disp(sprintf(&#039;  a = %g, lower limit of integration &#039;,a)) 
disp(sprintf(&#039;  b = %g, upper limit of integration &#039;,b)) 
disp(sprintf(&#039;  n = %g, number of Gauss Points, 1-10&#039;,n)) 
format short g

% Setup Gaussian Coefficients and Evaluation Points
x_values = zeros(10,10) ;
c_values = zeros(10,10) ;

i=1 ;
x_values(i,1) = 0.0 ;
c_values(i,1) = 2.0 ;

i=2 ;
x_values(i,1) = -0.5773502691896260 ; 
x_values(i,2) = -x_values(i,1) ;
c_values(i,1) = 1.0 ;
c_values(i,2) = 1.0 ;

i=3 ;
x_values(i,1) = -0.7745966692414830 ;
x_values(i,2) = 0.0 ;
x_values(i,3) = -x_values(i,1) ;
c_values(i,1) = 0.5555555555555560 ;
c_values(i,2) = 0.8888888888888890 ;
c_values(i,3) =c_values(i,1) ;

i=4 ;
x_values(i,1) = -0.8611363115940530 ;
x_values(i,2) = -0.3399810435848560 ;
x_values(i,3) = -x_values(i,2) ;
x_values(i,4) = -x_values(i,1) ;
c_values(i,1) = 0.3478548451374540 ;
c_values(i,2) = 0.6521451548625460 ;
c_values(i,3) =c_values(i,2) ;
c_values(i,4) =c_values(i,1) ;

i=5 ;
x_values(i,1) = -0.9061798459386640 ;
x_values(i,2) = -0.5384693101056830 ;
x_values(i,3) = 0.0 ;
x_values(i,4) = -x_values(i,2) ;
x_values(i,5) = -x_values(i,1) ;
c_values(i,1) = 0.2369368850561890 ;
c_values(i,2) = 0.4786386704993660 ;
c_values(i,3) = 0.5688888888888890 ;
c_values(i,4) =c_values(i,2) ;
c_values(i,5) =c_values(i,1) ;

i=6 ;
x_values(i,1) = -.9324695142032520 ;
x_values(i,2) = -.6612093864662650 ;
x_values(i,3) = -.2386191860831970 ;
x_values(i,4) = -x_values(i,3) ;
x_values(i,5) = -x_values(i,2) ;
x_values(i,6) = -x_values(i,1) ;
c_values(i,1) = 0.1713244923791700 ;
c_values(i,2) = 0.3607615730481390 ;
c_values(i,3) = 0.4679139345726910 ;
c_values(i,4) =c_values(i,3) ;
c_values(i,5) =c_values(i,2) ;
c_values(i,6) =c_values(i,1) ;

i=7 ;
x_values(i,1) = -0.9491079123427590 ;
x_values(i,2) = -0.7415311855993940 ;
x_values(i,3) = -0.4058451513773970 ;
x_values(i,4) = 0.0 ;
x_values(i,5) = -x_values(i,3) ;
x_values(i,6) = -x_values(i,2) ;
x_values(i,7) = -x_values(i,1) ;
c_values(i,1) = 0.1294849661688700 ;
c_values(i,2) = 0.2797053914892770 ;
c_values(i,3) = 0.3818300505051190 ;
c_values(i,4) = 0.4179591836734690 ;
c_values(i,5) =c_values(i,3) ;
c_values(i,6) =c_values(i,2) ;
c_values(i,7) =c_values(i,1) ;

i=8 ;
x_values(i,1) = -0.9602898564975360 ;
x_values(i,2) = -0.7966664774136270 ;
x_values(i,3) = -0.5255324099163290 ;
x_values(i,4) = -0.1834346424956500 ;
x_values(i,5) = -x_values(i,4) ;
x_values(i,6) = -x_values(i,3) ;
x_values(i,7) = -x_values(i,2) ;
x_values(i,8) = -x_values(i,1) ;
c_values(i,1) = 0.1012285362903760 ;
c_values(i,2) = 0.2223810344533740 ;
c_values(i,3) = 0.3137066458778870 ;
c_values(i,4) = 0.3626837833783620 ;
c_values(i,5) =c_values(i,4) ;
c_values(i,6) =c_values(i,3) ;
c_values(i,7) =c_values(i,2) ;
c_values(i,8) =c_values(i,1) ;

i=9 ;
x_values(i,1) = -0.9681602395076260 ;
x_values(i,2) = -0.8360311073266360 ;
x_values(i,4) = -0.6133714327005900 ;
x_values(i,4) = -0.3242534234038090 ;
x_values(i,5) = 0.0 ;
x_values(i,6) = -x_values(i,4) ;
x_values(i,7) = -x_values(i,3) ;
x_values(i,8) = -x_values(i,2) ;
x_values(i,9) = -x_values(i,1) ;
c_values(i,1) = 0.0812743883615740 ;
c_values(i,2) = 0.1806481606948570 ;
c_values(i,3) = 0.2606106964029350 ;
c_values(i,4) = 0.3123470770400030 ;
c_values(i,5) = 0.3302393550012600 ;
c_values(i,6) =c_values(i,4) ;
c_values(i,7) =c_values(i,3) ;
c_values(i,8) =c_values(i,2) ;
c_values(i,9) =c_values(i,1) ;

i=10 ;
x_values(i,1) = -0.9739065285171720 ;
x_values(i,2) = -0.8650633666889850 ;
x_values(i,3) = -0.6794095682990240 ;
x_values(i,4) = -0.4333953941292470 ;
x_values(i,5) = -0.1488743389816310 ;
x_values(i,6) = -x_values(i,5) ;
x_values(i,7) = -x_values(i,4) ;
x_values(i,8) = -x_values(i,3) ;
x_values(i,9) = -x_values(i,2) ;
x_values(i,10) = -x_values(i,1) ;
c_values(i,1) = 0.0666713443086880 ;
c_values(i,2) = 0.1494513491505810 ;
c_values(i,3) = 0.2190863625159820 ;
c_values(i,4) = 0.2692667193099960 ;
c_values(i,5) = 0.2955242247147530 ;
c_values(i,6) = c_values(i,5) ;
c_values(i,7) = c_values(i,4) ;
c_values(i,8) = c_values(i,3) ;
c_values(i,9) = c_values(i,2) ;
c_values(i,10) = c_values(i,1) ;


% simulation of the steps needed to perform the method begins here
disp(sprintf(&#039;\n*******************************Simulation*********************************\n&#039;))

% First do a lookup of the weights and evaluation locations
disp(&#039;1) Lookup Gauss point constants and evaluation points from the table.&#039;)
disp(&#039; &#039;)
disp(sprintf(&#039;     Evaluation Locations&#039;))

% Used to display the gauss points and their weights
for i=1:n
 disp(sprintf(&#039;      x%i=%1.15f&#039;,i,x_values(n,i)))
end
disp(&#039; &#039;)
disp(sprintf(&#039;     Evaluation Constants&#039;))
for i=1:n
 disp(sprintf(&#039;      c%i=%1.15f&#039;,i,c_values(n,i)))
end

% Now transform the evaluation points to fit your interval
disp(&#039; &#039;)
disp(&#039;2) The evaluation points are definted for the interval [-1,1]. Transform them&#039;)
disp(&#039; to fit [a,b]. This is done by the following formula&#039;)
disp(&#039; &#039;)
disp(&#039;     xnew = (b-a)/2 * x + (b+a)/2&#039;)
disp(&#039; &#039;)
disp(sprintf(&#039;     Transformed Evaluation Locations&#039;))
% Transform the locations
for i=1:n
 xv(i)=(b-a)/2*x_values(n,i)+(b+a)/2 ;
 disp(sprintf(&#039;      x%i=%1.15f&#039;,i,xv(i)))
end

% Apply the Gaussian procedure
disp(&#039; &#039;)
disp(&#039;3) The approximation is given as&#039;)
disp(&#039; &#039;)
disp(&#039;     approx = sum(i=1,n) c(i)*f(x(i))&#039;)

approx = 0 ;
disp(&#039; &#039;)

% approximation via Gaussian integration
for i=1:n-1
 approx = approx + c_values(n,i)*f(xv(i)) ;
 disp(sprintf(&#039;     %1.15f * %g&#039;,c_values(n,i),f(xv(i))))
end
approx = approx + c_values(n,n)*f(xv(n)) ;
disp(sprintf(&#039;    + %1.15f * %g&#039;,c_values(n,n),f(xv(n))))
disp(sprintf(&#039;    ----------------------------&#039;))
disp(sprintf(&#039;     %g&#039;,approx))

% the result must be scaled to fit the interval also
disp(&#039; &#039;)
disp(&#039;4) The approximation must also be transformed to fit the [a,b] interval.&#039;)
disp(&#039; &#039;)
disp(&#039;    approx = (b-a)/2 * approx&#039;)
disp(sprintf(&#039;        = %g * %g&#039;,(b-a)/2,approx))

% also transform the approximation
approx = (b-a)/2 * approx ;
disp(sprintf(&#039;        = %g&#039;,approx))

% Display results section
disp(sprintf(&#039;\n\n**************************Results****************************&#039;))

% The following finds what is called the &#039;Exact&#039; solution
exact = quad(f,a,b) ; 

disp(sprintf(&#039;\n Approximate = %g&#039;,approx)) 
disp(sprintf(&#039; Exact   = %g&#039;,exact))
disp(sprintf(&#039;\n True Error = Exact - Approximate&#039;)) 
disp(sprintf(&#039;       = %g - %g&#039;,exact,approx))
disp(sprintf(&#039;       = %g&#039;,exact-approx))
disp(sprintf(&#039;\n Absolute Relative True Error Percentage&#039;))
disp(sprintf(&#039;       = | ( Exact - Approximate ) / Exact | * 100&#039;))
disp(sprintf(&#039;       = | %g / %g | * 100&#039;,exact-approx,exact))
disp(sprintf(&#039;       = %g\n\n&#039;,abs( (exact-approx)/exact )*100))

**رضا سپاس یار** · 2009-01-22T20:41:36

پاسخ : انتگرال به روش ذوزنقه ای

روش سیمپسون:

کد:

function area = simpson(h,y)
%
% ------------------------------------------------------
%    Reza Sepas Yar
%    http://www.avr.ir
% ------------------------------------------------------
%
% area = simpson(h,y)
%
% Example:
% x = [0:0.25:1]
% y = 1/sqrt(1+x.^2)
% simpson(0.25,y)
% aiyi =
%  1.000000000000000  3.880570000581328  1.788854381999832  3.200000000000000  0.707106781186547
% ans =
%  0.881377596980642
%
% Function Simpson computes the area under a function defined by
% an odd number of ordinates in vector y spaced s units apart
% by Simpson&#039;s First Rule.
%

n=length(y);

aiyi(1) = y(1);
aiyi(n) = y(n);

for i=2:(n-1)
  if (mod(i,2) == 0)
    aiyi(i) = 4 * y(i);
  else
    aiyi(i) = 2 * y(i);
  end
end

aiyi

if (mod(n,2) == 0)
  disp(&#039;Y must have odd element!&#039;);
  area = 0;
else
  
  area = y(1)+y(n);
  for i=2:(n-1)
    if (mod(i,2) == 0)
      area = area + 4 * y(i);
    else
      area = area + 2 *y(i);
    end
  end
  area = h*area/3;
end

end

% Other Examples:

% x = 0 + 1e-10:pi/16:pi/2+1e-10
% x =
%  Columns 1 through 6
%  0.000000000100000  0.196349540949362  0.392699081798724  0.589048622648086  0.785398163497448  0.981747704346810
%  Columns 7 through 9
%  1.178097245196172  1.374446786045535  1.570796326894897
% y = sin(x)./x
% y =
%  Columns 1 through 6
%  1.000000000000000  0.993586851137686  0.974495358391543  0.943165320725420  0.900316316132506  0.846927992473694
%  Columns 7 through 9
%  0.784213303542454  0.713585487907166  0.636619772327053
% simpson(pi/16,y)
% aiyi =
%  Columns 1 through 6
%  1.000000000000000  3.974347404550744  1.948990716783087  3.772661282901681  1.800632632265012  3.387711969894778
%  Columns 7 through 9
%  1.568426607084908  2.854341951628665  0.636619772327053
% ans =
%  1.370764076042494

%====================================================================

% x = [0+1e-10:1/8:1+1e-10]
% x =
%  Columns 1 through 6
%  0.000000000100000  0.125000000100000  0.250000000100000  0.375000000100000  0.500000000100000  0.625000000100000
%  Columns 7 through 9
%  0.750000000100000  0.875000000100000  1.000000000100000
% y = (x.*log(x))./(1+x)
% y =
%  Columns 1 through 6
%  -0.000000002302585 -0.231049060262061 -0.277258872232701 -0.267498887164168 -0.231049060150788 -0.180770626589236
%  Columns 7 through 9
%  -0.123292316717300 -0.062314649841909  0.000000000050000
% simpson(1/8,y)
% aiyi =
%  Columns 1 through 6
%  -0.000000002302585 -0.924196241048244 -0.554517744465402 -1.069995548656670 -0.462098120301577 -0.723082506356943
%  Columns 7 through 9
%  -0.246584633434600 -0.249258599367635  0.000000000050000
% ans =
%  -0.176238891495152

%====================================================================

% x = [0:pi/16:pi/2]
% x =
%  Columns 1 through 6
%          0  0.196349540849362  0.392699081698724  0.589048622548086  0.785398163397448  0.981747704246810
%  Columns 7 through 9
%  1.178097245096172  1.374446785945535  1.570796326794897
% y = sin(x.^2)
% y =
%  Columns 1 through 6
%          0  0.038543592357938  0.153602060372138  0.340057724740335  0.578468789354558  0.821381311850125
%  Columns 7 through 9
%  0.983323424736072  0.949766418160832  0.624265952639699
% simpson(pi/16,y)
% aiyi =
%  Columns 1 through 6
%          0  0.154174369431752  0.307204120744276  1.360230898961340  1.156937578709116  3.285525247400499
%  Columns 7 through 9
%  1.966646849472143  3.799065672643329  0.624265952639699
% ans =
%  0.828205680955492
 
%====================================================================

% x = 0:1/8:1;
% y = 1./sqrt(2+(cos(x)).^2);
% simpson(1/8,y)
% aiyi =
%    0.57735    2.3154    1.1667    2.3628    1.2017    2.4536    1.2561    2.5762   0.66054
% ans =
%    0.6071
%

**رضا سپاس یار** · 2009-01-22T20:43:47

پاسخ : انتگرال به روش ذوزنقه ای

تمام این برنامه ها عملکردشون تست شده و صحیح کار می کنند. البته باید مانیان جان من رو ببخشند که در حضور ایشون جسارت کردم، به هر حال ایشون اکابر رو دارن و ما حالا حالا ها مونده به اکابر برسیم :smile:

**manian** · 2009-01-22T20:47:45

پاسخ : انتگرال به روش ذوزنقه ای

با سلام
ادعایی بر سواد ندارم.
اگر هدف فقط محاسبه با روش ذوزنقهآ‌ای است (نه پیادهآ‌سازی خود الگوریتم) میآ‌توانید به این صورت عمل کنید (خود کد گویا است):
(من کدآ‌های زیر را با octave در لینوکس تست کردم)

xmin=-1;xmax=3;
step = (xmax-xmin)/100;
X=xmin:step:xmax;
Y = sin(X) + X.^2;
I = trapz(X,Y)

I = 10.864

با تشکر از آقای سپاسیار، اگر علاقه داشته باشید مکتب ما برای ترم بعد دانشجو میآ‌گیره

اگر باعث ناراحتی آقای/خانم jafarzadeh85 شدهآ‌ام، عذر خواهی میآ‌کنم.

**رضا سپاس یار** · 2009-01-22T20:52:07

پاسخ : انتگرال به روش ذوزنقه ای

نوشته اصلی توسط manian

با تشکر از آقای سپاسیار، اگر علاقه داشته باشید مکتب ما برای ترم بعد دانشجو میآ‌گیره

مانیان عزیز و بزرگوار، با کمال میل حاضر در مکتب شما تلمذ کنم :nerd:

**رضا حاجی زاده** · 2009-01-22T20:59:45

پاسخ : انتگرال به روش ذوزنقه ای

من سوادم اون قدر نیست به خصوص که الان فهمیدم در اکابر به ما به اشتباه یاد دادهآ‌اند که Matlab به فارسی متلب نوشته میآ‌شود نه مطلب!

سلام
این نوشتن مت لب و مطلب ! هم واسه خودش داستانی شده
این ترم استاد مهندسی نرم افزار ما هم می گفت که خواهشا کسی به فارسی Matlab رو مطلب ننویسه که من بهش آلرژی دارم

**رضا سپاس یار** · 2009-01-22T21:08:39

پاسخ : انتگرال به روش ذوزنقه ای

سلام
این نوشتن مت لب و مطلب ! هم واسه خودش داستانی شده
این ترم استاد مهندسی نرم افزار ما هم می گفت که خواهشا کسی به فارسی Matlab رو مطلب ننویسه که من بهش آلرژی دارم

به نظر من بهتره این واژه رو باید با همون رسم الخط لاتین اش نوشت: MATLAB

**jafarzadeh85** · 2009-01-23T10:16:49

پاسخ : انتگرال به روش ذوزنقه ای

ضمن تشکر از همگی دوستان
جناب آقای سپاس یار عزیز با توجه به اینکه وقتی متن لاتینی که paste می شه به هم میریزه اگر براتون امکان داره m-file روش انتگرال گیری ذوزنقه ای رو برای من ایمیل کنید. :job:

jafarzadeh85@yahoo.com

باز هم از زحمتی که کشیدید ممنونم :wow:

**رضا حاجی زاده** · 2009-01-23T11:20:51

پاسخ : انتگرال به روش ذوزنقه ای

سلام
میتونی متن رو توی یه نت پدر بریزی از اون کپی بکنی توی محیط ادیتور MATLAB این طوری فکر نکنم مشکی پیش بیاد
موفق باشی

اطلاعیه

انتگرال به روش ذوزنقه ای

انتگرال به روش ذوزنقه ای

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه

دیدگاه