matlab_practice(1)

MATLAB practice(1) – POLYNOMIAL AND NUMERIC DIFFERENTIAL AND INTEGRATION

polynomial differential

polyval(a,x)

To show a polynomial $y=x^3-2x-5$, we need to provide a vector p such that each entry represent the coefficient of corresponding order of $x$

>> %differential and integration
>> % polynomial 
>> p = [1,0,-2,-5];% 多项式的系数
>> x = -2:0.1:5;
>> f = polyval(p,x);% polyval(p,x)通过系数p和变量x得到值域f
>> plot(x,f,'linewidth',2);
>> xlabel('x');
>> ylabel('y');
>> ylabel('f(x)');
>> set(gca,'fontsize',14);

polyder(p)

To get the differential of polynomial, use polyder§, where p is the coefficent of $f(x)$

>> % polyder(p)
>> p = [5,0,-2,0,1];
>> polyder(p)

ans =

    20     0    -4     0
% know f'(7)
>> ans = polyval(polyder(p),7)

ans = 

    6832

exercise: plot $f(x)$ and $f^{\prime }(x)$ of
$f(x)=(20x^3-7x^2+5x+10)(4x^2+12x-3)$

answer: conv(u, v)

>> %f = conv(u,v)返回两个向量u，v的卷积
>> %如果这两个向量代表多项式大的系数，那么返回的结果就是多项式相乘得到的系数向量
>> %代表了卷积后的多项式
>> f = conv([20,-7,5,10],[4,12,-3]);
>> x = linspace(-2,1,100);
>> y1 = polyval(f,x);
>> y2 = polyval(polyder(f),x);

>> hold on;
>> h1 = plot(x,y1,'b--');
>> h2 = plot(x,y2,'r');

>> legend('f(x)','f"(x)');

polynomial integration

polyint(p,k)

polynomial integration is

>> % 要提供积分一个积分常数
>> % polyint(p,k), where k is a constant
>> p = [5,0,-2,0,1];
>> f = polyint(p,3);% k = 3
>> x = linspace(-1,2,100);
>> y = polyval(f,x);
>> plot(x,y,'r--');
>> polyval(f,7)

ans =

   1.6588e+04

numeric differential

diff(x)

>> %diff(x) x is a vector 
>> x = [1,2,5,2,1];
>> diff(x)

ans =

     1     3    -3    -1

>> clear
>> x = [1 2];
>> clear
>> x0 = pi/2;
>> h = 0.1;
>> x = [x0,x0+h];
>> y = [sin(x0),sin(x0+h)];
>> m = diff(y)./diff(x);
>> m

m =

   -0.0500
% if x is a vector, then we can get the differential of each x if h is small enough

different slide h can result in different answers:

>> x1 = 0:0.1:2*pi;
>> x2 = 0:0.01:2*pi;
>> x3 = 0:0.001:2*pi;

>> fx1 = exp(-x1).*sin((x1.^2)./2);
>> fx2 = exp(-x2).*sin((x2.^2)./2);
>> fx3 = exp(-x3).*sin((x3.^2)./2);

>> fprimex1 = diff(fx1)./diff(x1);
>> fprimex2 = diff(fx2)./diff(x2);
>> fprimex3 = diff(fx3)./diff(x3);

>> hold on;
>> plot(x1(1:62),fprimex1,'--');
>> plot(x2(1:628),fprimex2,'.');
>> plot(x3(1:6283),fprimex3,'.-');
>> plot(x,fx);
>> legend('h=0.1','h=0.01','h=0.001','f(x)');

numeric integration

midpoint rule

>> h = 0.01; %slide
>> x = 0:h:2;
>> midpoint = (x(1:end-1)+x(2:end))./2; % mdipoint rectangle rule
>> y = 4*midpoint.^3;
>> s = sum(h*y);
>> s

s =

   15.9998

trapz(y)

>> % trapz()梯形预测
>> h = 0.05;
>> x = 0:h:2;
>> y = 4*x.^3;
>> s = h*trapz(y);
>> s

s =

   16.0100

trapezoid rule

>> h = 0.05;
>> x = 0:h:2;
>> y = 4*x.^3;
>> trapezoid = (y(1:end-1)+y(2:end))/2;
>> s = h*sum(trapezoid);
>> s

s =

   16.0100
% this is an alternative method of trapz(y)