MATLAB07:数值计算

╰+攻爆jí腚メ 2022-12-04 02:22 343阅读 0赞

转自：[https://blog.csdn.net/ncepu\_Chen/article/details/103165266][https_blog.csdn.net_ncepu_Chen_article_details_103165266]

**目录**

多项式的数值运算

使用MATLAB表示多项式

使用向量表示多项式

多项式求值:polyval()

多项式的乘法:conv()

多项式的数值运算

多项式的因式分解:roots()

多项式的微分:polyder()

多项式的积分:polyint()

非线性表达式的数值运算

方程(组)求根fsolve()

数值微分

求差分:diff()

求导数:diff(y)./diff(x)

数值积分

数值积分原理

数值积分函数:integral()

--------------------

学习一门技术最好的方式就是阅读官方文档,可以查看[MATLAB官方文档][MATLAB]

# 多项式的数值运算 #

## 使用MATLAB表示多项式 ##

### 使用向量表示多项式 ###

在MATLAB中,多项式可以用向量表示,向量中的元素为多项式的系数(降幂排序):第一位为多项式最高次项系数,最后一位为常数项.

例如多项式:

![f(x) = x^\{3\}-2x-5][f_x_ _ x_3_-2x-5]

可以用向量`p = [1 0 -2 -5]`表示.

### 多项式求值:`polyval()` ###

使用`polyval(p, x)`可以计算多项式`p`在`x`的每个点处的值.

a = [9,-5,3,7]; x = -2:0.01:5;
    f = polyval(a,x);
    plot(x,f);

![在这里插入图片描述][20191120155802665.jpg]

### 多项式的乘法:`conv()` ###

使用`conv(p1, p2)`函数可以对两个向量`p1`和`p2`进行卷积相乘,用于计算多项式的乘法.

例如多项式:  
                                                       ![f(x)=(x^2+1) (2x+7)][f_x_x_2_1_ _2x_7]

可以使用`conv()`函数得到展开后的多项式:

p = conv([1 0 1], [2 7])

得到`p = [2 7 2 7]`.

## 多项式的数值运算 ##

### 多项式的因式分解:`roots()` ###

使用`roots(p)`函数可以对多项式`p`进行因式分解,即求表达式值为0的根.

p = roots([1 -3.5 2.75 2.125 -3.875 1.25])

得到`p = [2 -1, 1+0.5i, 1-0.5i, 0.5]`,表示

![x^\{5\}-3.5x^\{4\}+2.75x^\{3\}+2.125x^\{2\}-3.875x+1.25\\\\ =(x-2)(x+1)(x-1-0.5i)(x-1+0.5i)(x-0.5)][x_5_-3.5x_4_2.75x_3_2.125x_2_-3.875x_1.25_ _x-2_x_1_x-1-0.5i_x-1_0.5i_x-0.5]

### 多项式的微分:`polyder()` ###

使用`polyder(p)`函数可以计算多项式的导数.

例如对下面多项式求导:

![f(x)= 5x^4 - 2x^2 + 1][f_x_ 5x_4 - 2x_2 _ 1]

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

得到`p = [20 0 -4 0]`,表示计算得到导数![f^\{'\}(x)= 20 x^3 - 4x][f_x_ 20 x_3 - 4x]

### 多项式的积分:`polyint()` ###

使用`polyint(p, k)`函数可以计算多项式`p`的积分,积分结果的常数项设为`k`.

例如对下面多项式求导:

![f(x)= 5x^4 - 2x^2 + 1][f_x_ 5x_4 - 2x_2 _ 1]

p = polyint([5 0 -2 0 1], 3)

得到`p = [1 0 -0.6667 0 1 3]`,表示计算得到积分

![\\int f(x) dx = x^5 -0.6667x^3 + x + 3][int f_x_ dx _ x_5 -0.6667x_3 _ x _ 3]

# 非线性表达式的数值运算 #

## 方程(组)求根`fsolve()` ##

使用`fsolve(fun, x0)`求非线性方程组的根,`fun`为待求方程的函数句柄,`x0`为初值.

1.  求方程![1.2x+x\\sin(x)+0.3=0][1.2x_x_sin_x_0.3_0]在![x=0][x_0]附近的解.
    
        f2 = @(x) (1.2*x+x*sin(x)+0.3);
        fsolve(f2,0)	% 得到 -0.3500
2.  解方程组

![\\begin\{aligned\}\\left\\\{\\begin\{aligned\}e^\{-e^\{-(x\_1+x\_2)\}\} - x\_2(1+x\_1^2) = 0 \\\\x\_1 \\cos x\_2 + x\_2 \\sin x\_1 - \\frac\{1\}\{2\} = 0\\end\{aligned\}\\right.\\end\{aligned\}][begin_aligned_left_begin_aligned_e_-e_-_x_1_x_2_ - x_2_1_x_1_2_ _ 0 _x_1 _cos x_2 _ x_2 _sin x_1 - _frac_1_2_ _ 0_end_aligned_right._end_aligned]

设定初值为\[0,0\]

fun = @(x) [exp(-exp(-(x(1)+x(2))))-x(2)*(1+x(1)^2)...
        x(1)*cos(x(2)) + x(2)*sin(x(1)) - 0.5]
    x0 = [0,0];
    x = fsolve(fun,x0)	% 得到[0.3532 0.6061]

## 数值微分 ##

### 求差分:`diff()` ###

使用`diff(X, n)`计算向量`X`的`n`阶差分,`n`默认为`1`.

x = [1 2 5 2 1];
    diff(x);	% 得到 [1 3 -3 1]
    diff(x,1);	% 得到 [1 3 -3 1]
    diff(x,2);	% 得到 [2 -6 2]

### 求导数:`diff(y)./diff(x)` ###

使用导数的定义

![f^\{'\}(x\_0) = \\lim\_\{h\\rightarrow 0\}\\frac\{f(x\_0+h)-f(x\_0)\}\{h\}][f_x_0_ _ _lim_h_rightarrow 0_frac_f_x_0_h_-f_x_0_h]

可以计算函数在某点的近似导数.

x0 = pi/2; h = 0.1;
    x = [x0 x0+h];
    y = [sin(x0) sin(x0+h)];
    m = diff(y)./diff(x)		% 得到 m = -0.005

`h`的取值越小,得到的近似导数越精确.

![在这里插入图片描述][20191120163442489.png]

下面程序计算![f(x) = x^3][f_x_ _ x_3]的一阶和二阶导数的值.

x = -2:0.005:2; y = x.^3;
    m = diff(y)./diff(x);			% 计算一阶导数
    m2 = diff(m)./diff(x(1:end-1));	% 计算二阶导数
    
    plot(x,y,x(1:end-1),m,x(1:end-2),m2);
    xlabel('x'); ylabel('y');
    legend('f(x) =x^3','f''(x)','f''''(x)', 4);

![在这里插入图片描述][20191120163919851.jpg]

## 数值积分 ##

### 数值积分原理 ###

有三种常见算法用于计算数值积分: **矩形法**,**梯形法**,**抛物线法**,它们分别把微分区间的图形视为矩形,梯形,抛物线以计算面积.

<table> 
 <thead> 
  <tr> 
   <th style="width:230px;">算法</th> 
   <th style="width:158px;">图示</th> 
   <th style="width:516px;">表达式</th> 
  </tr> 
 </thead> 
 <tbody> 
  <tr> 
   <td style="width:230px;">矩形法(Midpoint Rule)</td> 
   <td style="width:158px;"><img alt="在这里插入图片描述" src="https://img-blog.csdnimg.cn/20191120165104999.png"></td> 
   <td style="width:516px;"><img alt="\int_a^b f(x) dx = 2h \sum_{i=0}^{(n/2)-1} f(x_{2i+1})" src="https://private.codecogs.com/gif.latex?%5Cint_a%5Eb%20f%28x%29%20dx%20%3D%202h%20%5Csum_%7Bi%3D0%7D%5E%7B%28n/2%29-1%7D%20f%28x_%7B2i&amp;plus;1%7D%29"></td> 
  </tr> 
  <tr> 
   <td style="width:230px;">梯形法(Trapezoid Rule)</td> 
   <td style="width:158px;"><img alt="在这里插入图片描述" src="https://img-blog.csdnimg.cn/20191120165114624.png"></td> 
   <td style="width:516px;"><img alt="\int_a^b f(x) dx = \frac{h}{2} \left[ f(a) + 2\sum_{i=1}^{n-1} f(x_i) +f(b) \right]" src="https://private.codecogs.com/gif.latex?%5Cint_a%5Eb%20f%28x%29%20dx%20%3D%20%5Cfrac%7Bh%7D%7B2%7D%20%5Cleft%5B%20f%28a%29%20&amp;plus;%202%5Csum_%7Bi%3D1%7D%5E%7Bn-1%7D%20f%28x_i%29%20&amp;plus;f%28b%29%20%5Cright%5D"></td> 
  </tr> 
  <tr> 
   <td style="width:230px;">抛物线法(Simpson’s Rule)</td> 
   <td style="width:158px;"><img alt="在这里插入图片描述" src="https://img-blog.csdnimg.cn/20191120165124685.png"></td> 
   <td style="width:516px;"><img alt="\int_a^b f(x) dx = \frac{h}{3} \left[ f(a) + 2\sum_{i=1}^{(n/2)-1} f(x_{2i}) + 4\sum_{i=1}^{n/2} f(x_{2i-1}) + f(b) \right]" src="https://private.codecogs.com/gif.latex?%5Cint_a%5Eb%20f%28x%29%20dx%20%3D%20%5Cfrac%7Bh%7D%7B3%7D%20%5Cleft%5B%20f%28a%29%20&amp;plus;%202%5Csum_%7Bi%3D1%7D%5E%7B%28n/2%29-1%7D%20f%28x_%7B2i%7D%29%20&amp;plus;%204%5Csum_%7Bi%3D1%7D%5E%7Bn/2%7D%20f%28x_%7B2i-1%7D%29%20&amp;plus;%20f%28b%29%20%5Cright%5D"></td> 
  </tr> 
 </tbody> 
</table>

![在这里插入图片描述][20191120180550884.png]

下面分别使用三种方法计算![f(x) = 4x^3][f_x_ _ 4x_3]在区间(0,2)内的积分.

h = 0.05; x = 0:h:2;
    y = 4*midpoint.^3;
    
    % 使用矩形法计算近似积分
    midpoint = (x(1:end-1)+x(2:end))./2;
    s = sum(h*y)			% 得到 15.9950
    
    % 使用梯形法计算近似积分
    trapezoid = (y(1:end-1)+y(2:end))/2;
    s = h*sum(trapezoid)	% 得到 15.2246
    
    % 使用抛物线法计算数值积分
    s = h/3*(y(1)+2*sum(y(3:2:end-2))+4*sum(y(2:2:end))+y(end))	% 得到 15.8240

### 数值积分函数:`integral()` ###

`integral()`,`integral2()`,`integral3()`分别对函数在`xmin`至`xmax`间进行一重,二重,三重积分.

它们的第一个参数都应该是一个函数句柄,下面例子演示他们的用法:

1.  计算![\\int\_0^2 \\frac\{1\}\{x^3-2x-5\}][int_0_2 _frac_1_x_3-2x-5]
    
        f = @(x) 1./(x.^3-2*x-5);
        integral(f,0,2)		% 得到 -0.4605
2.  计算![\\int\_0^\\pi \\int\_\\pi^\{2\\pi\} (y\\sin(x) + x \\cos(y)) dx dy][int_0_pi _int_pi_2_pi_ _y_sin_x_ _ x _cos_y_ dx dy]
    
        f = @(x,y) y.*sin(x)+x.*cos(y);
        integral2(f,pi,2*pi,0,pi)	% 得到 -9.8696
3.  计算![\\int\_\{-1\}^1 \\int\_0^1 \\int\_0^\\pi (y\\sin(x) + z \\cos(y)) dx dy dz][int_-1_1 _int_0_1 _int_0_pi _y_sin_x_ _ z _cos_y_ dx dy dz]
    
        f = @(x,y,z) y.*sin(x)+z.*cos(y);
        integral3(f,0,pi,0,1,-1,1)

[https_blog.csdn.net_ncepu_Chen_article_details_103165266]: https://blog.csdn.net/ncepu_Chen/article/details/103165266
[MATLAB]: https://ww2.mathworks.cn/help/matlab/getting-started-with-matlab.html?lang=en
[f_x_ _ x_3_-2x-5]: https://private.codecogs.com/gif.latex?f%28x%29%20%3D%20x%5E%7B3%7D-2x-5
[20191120155802665.jpg]: /images/20221123/15af0aec26b4402e88dc560de4d7648c.png
[f_x_x_2_1_ _2x_7]: https://private.codecogs.com/gif.latex?f%28x%29%3D%28x%5E2&plus;1%29%20%282x&plus;7%29
[x_5_-3.5x_4_2.75x_3_2.125x_2_-3.875x_1.25_ _x-2_x_1_x-1-0.5i_x-1_0.5i_x-0.5]: https://private.codecogs.com/gif.latex?x%5E%7B5%7D-3.5x%5E%7B4%7D&plus;2.75x%5E%7B3%7D&plus;2.125x%5E%7B2%7D-3.875x&plus;1.25%5C%5C%20%3D%28x-2%29%28x&plus;1%29%28x-1-0.5i%29%28x-1&plus;0.5i%29%28x-0.5%29
[f_x_ 5x_4 - 2x_2 _ 1]: https://private.codecogs.com/gif.latex?f%28x%29%3D%205x%5E4%20-%202x%5E2%20&plus;%201
[f_x_ 20 x_3 - 4x]: https://private.codecogs.com/gif.latex?f%5E%7B%27%7D%28x%29%3D%2020%20x%5E3%20-%204x
[int f_x_ dx _ x_5 -0.6667x_3 _ x _ 3]: https://private.codecogs.com/gif.latex?%5Cint%20f%28x%29%20dx%20%3D%20x%5E5%20-0.6667x%5E3%20&plus;%20x%20&plus;%203
[1.2x_x_sin_x_0.3_0]: https://private.codecogs.com/gif.latex?1.2x&plus;x%5Csin%28x%29&plus;0.3%3D0
[x_0]: https://private.codecogs.com/gif.latex?x%3D0
[begin_aligned_left_begin_aligned_e_-e_-_x_1_x_2_ - x_2_1_x_1_2_ _ 0 _x_1 _cos x_2 _ x_2 _sin x_1 - _frac_1_2_ _ 0_end_aligned_right._end_aligned]: https://private.codecogs.com/gif.latex?%5Cbegin%7Baligned%7D%5Cleft%5C%7B%5Cbegin%7Baligned%7De%5E%7B-e%5E%7B-%28x_1&plus;x_2%29%7D%7D%20-%20x_2%281&plus;x_1%5E2%29%20%3D%200%20%5C%5Cx_1%20%5Ccos%20x_2%20&plus;%20x_2%20%5Csin%20x_1%20-%20%5Cfrac%7B1%7D%7B2%7D%20%3D%200%5Cend%7Baligned%7D%5Cright.%5Cend%7Baligned%7D
[f_x_0_ _ _lim_h_rightarrow 0_frac_f_x_0_h_-f_x_0_h]: https://private.codecogs.com/gif.latex?f%5E%7B%27%7D%28x_0%29%20%3D%20%5Clim_%7Bh%5Crightarrow%200%7D%5Cfrac%7Bf%28x_0&plus;h%29-f%28x_0%29%7D%7Bh%7D
[20191120163442489.png]: /images/20221123/89ecb3f4a9ff492aabbf1cecb6f20e79.png
[f_x_ _ x_3]: https://private.codecogs.com/gif.latex?f%28x%29%20%3D%20x%5E3
[20191120163919851.jpg]: /images/20221123/97c270106ec54819824caa92411fc954.png
[20191120180550884.png]: /images/20221123/4e3a9b16dca84e3296137b165e5327b8.png
[f_x_ _ 4x_3]: https://private.codecogs.com/gif.latex?f%28x%29%20%3D%204x%5E3
[int_0_2 _frac_1_x_3-2x-5]: https://private.codecogs.com/gif.latex?%5Cint_0%5E2%20%5Cfrac%7B1%7D%7Bx%5E3-2x-5%7D
[int_0_pi _int_pi_2_pi_ _y_sin_x_ _ x _cos_y_ dx dy]: https://private.codecogs.com/gif.latex?%5Cint_0%5E%5Cpi%20%5Cint_%5Cpi%5E%7B2%5Cpi%7D%20%28y%5Csin%28x%29%20&plus;%20x%20%5Ccos%28y%29%29%20dx%20dy
[int_-1_1 _int_0_1 _int_0_pi _y_sin_x_ _ z _cos_y_ dx dy dz]: https://private.codecogs.com/gif.latex?%5Cint_%7B-1%7D%5E1%20%5Cint_0%5E1%20%5Cint_0%5E%5Cpi%20%28y%5Csin%28x%29%20&plus;%20z%20%5Ccos%28y%29%29%20dx%20dy%20dz