the question is in the fileComputer Project 3: Numerical Integration Note: Make sure to use the...

Question

the question is in the fileComputer Project 3: Numerical Integration Note: Make sure to use the matlab ’format long e’ to print your results. Write a report that contains for each problem the matlab program, results, and discussion of results Problem I: 1. Write a matlab m-file program trapez.m as trapez(’f’,a,b,n) and returns the trapezoidal approximation T (f, h), h = (b−a)/n. Define your func- tion in a m-file as function f=myf(x) f=.... Use n = 15 and run trapez to compute the integral:∫ 1 −0.2 exp(x2)dx Print the result. 2. Write a matlab m-file program romberg.m which should be executed by calling romberg(’f’,a,b,delta,nmax) which in turn uses trapez.m, with n = 2m, to compute Romberg approximations. Your program should stop when either |R(m,m) − R(m − 1,m − 1)|  nmax and return the approximation R(m,m). Note: your function should be defined in a separate m-file f.m function f = f(x) f=.......... Inside trapez.m you should call to obtain f(x) as feval(f,x) Use delta = 10−8 and nmax = 30 and test you program on the integral∫ 3 1 cos(2x2)exp(x)dx 1 3. Print the Romberg table. Problem II: 1. Write a function as an m.file mygauss.m that should be called as my- gauss(’f’,a,b,ng) to approximate the integrate ∫ b a f(x)dx using ng Gauss points. The mfile should start with: function I=mygauss(f,a,b,ng) ... ... I= define ξ and their corresponding weights w from the file gauss.m posted in Scholar under Resources for ng=5. Shift them to [a,b] and use the Gauss-Legendre quadrature. Hint: inside mygauss.m use feval(f,s) to compute f(s) for a given value s. 2. Test your program on I1 = ∫ 1 0 x 4dx to obtain a zero error. 3. Run your program to compute the following integral and print your result I2 = ∫ 2 1 5cos(2x2)exp(x)dx. 4. Consider the curve described by the parametric equations (x(t) = exp(t)cos(t2), y(t) = sin(5t)), t ∈ [0, 2] and define its length as l =∫ 2 0 √ x′2 + y′2dt. Divide the interval [0, 2] into m = 2, 4, 8, 16, 32, 64 uni- form subintervals and use the composite Gauss quadrature with ng=4 to compute the length of the curve. Compute the errors for each case by assuming the true length using m = 50, ng = 10. Print a table with the first column contains the values of m, the second contains the lengths of the curve and the third column contains the absolute value of the errors. 5. Print and discuss all results. Problem III: 2 1. Write a function as an m.file mygauss2d.m that should be called as my- gauss2d(’f2d’,ng) to approximate the double integral ∫ b a ∫ d c f(x, y)dydx. on the domain 0

Rahul · Accepted Answer

Math_280420/c_2.m
% Matlab Code to change string format function to function handle
function fun = c_2(f)
str = strcat('@','(x)',f)
fun = str2func(str)
end
Math_280420/d_2.m
% Matlab Code to change string format function to function handle
function fun = c_2(f)
str = strcat('@','(x)',f)
fun = str2func(str)
end
Math_280420/f2d.m
% Matlab Code to change 2 d string format function to function handle
function fun = f2d(f)
str = strcat('@','(x,y)',f)
fun = str2func(str)
end
Math_280420/gauss.m
function [xi,w] = gauss(ng)
%  for ng=1,2,3,4,5,6,7 and 10.
%  call as 
% >[xi,w]=gauss(5) 
% for 5 point-quadrature
%
switch ng  
case 1
  xi = [0.0];
  w  = [2.0];
case 2
  xi = [-1.0/sqrt(3.0),1.0/sqrt(3.0)];
  w  = [1.0,1.0];
 
case 3
  xi = [-sqrt(3.0/5.0), 0.0, sqrt(3.0/5.0)];
  w  = [ 5.0/9.0, 8.0/9.0,5.0/9.0];
case 4
  xi = [-0.861136311594053,-0.339981043584856];
  xi = [xi,-xi];
  w  = [ 0.347854845137454, 0.652145154862546];
  w  = [w,w];
case 5
  xi = [-0.906179845938664,-0.538469310105683];
  xi = [xi,0,-xi];
  w  = [0.236926885056189,0.478628670499366];
  w  = [w,0.568888888888889,w];
case 6
  xi = [-0.932469514203152, -0.661209386466265,- 0.238619186083197];
  xi = [xi,-xi];    
  w  = [0.171324492379170,0.360761573048139,0.467913934572691];
  w  = [w,w];
case 7
  xi = [-0.949107912342759,-0.741531185599394,-0.405845151377397];
  xi = [ xi,0.0,-xi];
  w  = [0.129484966168870,0.279705391489277, 0.381830050505119];
  w  = [w,0.417959183673469,w];
otherwise
  xi = [ -0.973906528517172,-0.865063366688985,-0.679409568299024,-0.433395394129247,-0.148874338981631];
  xi = [ xi,-xi];
  w  = [ 0.066671344308688,0.149451349150581,0.219086362515982,0.269266719309996,0.295524224714753];
  w  = [ w ,w ];
end;
Math_280420/myf.m
% Matlab function to change string to function in this independent variable
% is x
% Input f in string format
% Ouput function_handle
function fun = myf(f)
str = strcat('@','(x)',f)
fun = str2func(str)
end
Math_280420/myf_1.m
% Matlab function to change string to function in this independent variable
% is x
% Input f in string format
% Ouput function_handle
function fun = myf_1(f)
str = strcat('@','(t)',f)
fun = str2func(str)
end
Math_280420/mygauss.m
% Matlab code to find the value of integration by gauss method
%% Input:
% f : function in the string format
% a = Lower Limit
% b = Upper limit
% ng = Number of gauss points
%%Function Run
function I = mygauss(f,a,b,ng)
[xi,w] = gauss(ng); % Find the xi and weight from the gauss function
f = myf(f);% Change string function to function_handle
I = 0; % Initialise the integration with zero value
% Shifting of the x = (a+b)/2 + (b-a)*xi/2 and dx = (b-a)/2*dxi 
for i = 1:ng % Loop to find the integration value
    I = I + w(i)*(b-a)/2*feval(f,((b-a)*xi(i)/2)+(b+a)/2);
end
%% output: I is the value of the integration
disp(['Value of the integration is = ',num2str(I)]);
end
Math_280420/mygauss_1.m
% Matlab code to find the value of integration by gauss method
%% Input:
% f : function in the string format
% a = Lower Limit
% b = Upper limit
% ng = Number of gauss points
%%Function Run
function I = mygauss_1(f,a,b,ng)
[xi,w] = gauss(ng); % Find the xi and weight from the gauss function
f = myf_1(f);% Change string function to function_handle
I = 0; % Initialise the integration with zero value
% Shifting of the x = (a+b)/2 + (b-a)*xi/2 and dx = (b-a)/2*dxi 
for i = 1:ng % Loop to find the integration value
    I = I + w(i)*(b-a)/2*feval(f,((b-a)*xi(i)/2)+(b+a)/2);
end
%% output: I is the value of the integration
disp(['Value of the integration is = ',num2str(I)]);
end
Math_280420/mygauss2D.asv
% Matlab Code to find double integration value by using gauss quadrature
function I = mygauss(f,a,b,ng)
%ng = 2;
[xi,w1] = gauss(ng);% find the weight by using ng point
[yi,w2] = gauss(ng);
%f = '(x^2 + y^2)*cos(x)';
f = f2d(f);% change string format function to function handle
% a = 0;
% b = 2.5;
% c = @(x)cos(x);
% d = @(x)2*exp(x);
% d_1 = @(x)2*exp(x);
% c_1 = @(x)-1*sin(x);
Iy = 0;% Intial valye of integration is zero
for i = 1:ng % loop to find the value of the inegration
    %shifted quadrature are obtained using the mapping from 
    Iy = Iy + w1(i)*w2(i)*((d(a + ((b-a)*(1+xi(i))/2))-c(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        ((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2)) %+   (1+yi(i))*w1(i)*w1(i)*((d_1(a + ((b-a)*(1+xi(i))/2))-c_1(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        %((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2))     
end
disp(['Value of the integration is = ',num2str(Iy)]);
end
Math_280420/mygauss2D.m
% Matlab Code to find double integration value by using gauss quadrature
function Iy = mygauss2D(f,a,b,c,d,ng)
%ng = 2;
[xi,w1] = gauss(ng);% find the weight by using ng point
[yi,w2] = gauss(ng);
%f = '(x^2 + y^2)*cos(x)';
f = f2d(f);% change string format function to function handle
% a = 0;
% b = 2.5;
% c = @(x)cos(x);
% d = @(x)2*exp(x);
% d_1 = @(x)2*exp(x);
% c_1 = @(x)-1*sin(x);
Iy = 0;% Intial valye of integration is zero
for i = 1:ng % loop to find the value of the inegration
    %shifted quadrature are obtained using the mapping from [-1,1] to
    %[a,b]x [c,d]
    % x = a + (b-a)*(1 + xi)/2
    % y = c(xi) + (d(a + (b-a)*(1 + xi)/2) - c(a + (b-a)*(1 + xi)/2)*(1 + yi)/2
    Iy = Iy + w1(i)*w2(i)*((d(a + ((b-a)*(1+xi(i))/2))-c(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        ((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2)) %+   (1+yi(i))*w1(i)*w1(i)*((d_1(a + ((b-a)*(1+xi(i))/2))-c_1(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        %((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2))     
end
disp(['Value of the integration is = ',num2str(Iy)]);
end
Math_280420/mygauss2D_m.asv
% Matlab Code to find double integration value by using gauss quadrature
function Iy = mygauss2D_m(f,a,b,c_1,d,ng)
%ng = 2;
[xi,w1] = gauss(ng);% find the weight by using ng point
[yi,w2] = gauss(ng);
%f = '(x^2 + y^2)*cos(x)';
f = f2d(f);% change string format function to function handle
c = c(c);
d = d(d);
% a = 0;
% b = 2.5;
% c = @(x)cos(x);
% d = @(x)2*exp(x);
% d_1 = @(x)2*exp(x);
% c_1 = @(x)-1*sin(x);
Iy = 0;% Intial valye of integration is zero
for i = 1:ng % loop to find the value of the inegration
    %shifted quadrature are obtained using the mapping from [-1,1] to
    %[a,b]x [c,d]
    % x = a + (b-a)*(1 + xi)/2
    % y = c(xi) + (d(a + (b-a)*(1 + xi)/2) - c(a + (b-a)*(1 + xi)/2)*(1 + yi)/2
    Iy = Iy + w1(i)*w2(i)*((d(a + ((b-a)*(1+xi(i))/2))-c(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        ((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2)) %+   (1+yi(i))*w1(i)*w1(i)*((d_1(a + ((b-a)*(1+xi(i))/2))-c_1(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        %((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2))     
end
disp(['Value of the integration is = ',num2str(Iy)]);
end
Math_280420/mygauss2D_m.m
% Matlab Code to find double integration value by using gauss quadrature
function Iy = mygauss2D_m(f,a,b,c_1,d_1,ng)
%ng = 2;
[xi,w1] = gauss(ng);% find the weight by using ng point
[yi,w2] = gauss(ng);
%f = '(x^2 + y^2)*cos(x)';
f = f2d(f);% change string format function to function handle
c = c_2(c_1);
d = d_2(d_1);
% a = 0;
% b = 2.5;
% c = @(x)cos(x);
% d = @(x)2*exp(x);
% d_1 = @(x)2*exp(x);
% c_1 = @(x)-1*sin(x);
Iy = 0;% Intial valye of integration is zero
for i = 1:ng % loop to find the value of the inegration
    %shifted quadrature are obtained using the mapping from [-1,1] to
    %[a,b]x [c,d]
    % x = a + (b-a)*(1 + xi)/2
    % y = c(xi) + (d(a + (b-a)*(1 + xi)/2) - c(a + (b-a)*(1 + xi)/2)*(1 + yi)/2
    Iy = Iy + w1(i)*w2(i)*((d(a + ((b-a)*(1+xi(i))/2))-c(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        ((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2)) %+   (1+yi(i))*w1(i)*w1(i)*((d_1(a + ((b-a)*(1+xi(i))/2))-c_1(a + ((b-a)*(1+xi(i))/2)))/2)*((b-a)/2)*((b-a)/2)*f(a + ((b-a)*(1+xi(i))/2),c(xi(i)) + ...
        %((d(a + ((b-a)*(1+xi(i))/2))-c(d(a + ((b-a)*(1+xi(i))/2))))*(1 + yi(i))/2))     
end
disp(['Value of the integration is = ',num2str(Iy)]);
end
Math_280420/problem_2_4.m
format long e
% Matlab Code to find the length of the curve for given function
% Input is m  and in this finding the length of the for different value of
% m
m = [2 4 8 16 32 64];
s = zeros(length(m),1); % length of the curve for corresponding to particular value of m
for k = 1:length(m) % for loop to take different m from the array
    a = 0; % lower lomit
    b = 2; % Upper Limit
    % dividing the domain of the integration to number of intervals and m is
    % number of intervel in domain range is divided
    % for eg if 0 to 2 is domain and m is 2.

Computer Project 3: Numerical Integration Note: Make sure to use the matlab ’format long e’ to print your results. Write a report that contains for each problem the matlab program, results, and...

Answer To: Computer Project 3: Numerical Integration Note: Make sure to use the matlab ’format long e’ to print...

Answer To This Question Is Available To Download

Related Questions & Answers

Submit New Assignment