matlab从无到有系列文章汇总：

matlab从无到有系列（一）：认识matlab

matlab从无到有系列（二）：矩阵运算基础

matlab从无到有系列（三）：数值计算基础

全文共17248个字，码字总结不易，老铁们来个三连：点赞、关注、评论

作者：左手の明天

原创不易，转载请联系作者并注明出处

敬请期待下一篇：

matlab从无到有系列（五）：基本图形处理功能

目录

1、符号对象的创建

1.1 sym函数

1.2 syms函数

2、符号表达式的化简和替换

2.1 创建符号表达式

2.1.1 符号表达式定义：

2.1.2 建立代数式的方法：

2.1.3 建立符号方程

2.2 符号表达式的化简

2.2.1 因式分解

2.2.2 展开

2.2.3 同类项合并

2.2.4 化简

2.2.5 分式通分

2.2.6 嵌套

2.3 符号表达式的替换

2.3.1 函数subs

2.3.2 函数subexpr

3、符号微积分

3.1 符号积分

3.1.1 limit函数

3.1.2 int函数

3.2 符号微分

3.3 积分变换

3.3.1 Fourier变换

3.3.2 Laplace变换

3.3.3 Z变换

3.4 符号求和

3.5 Taylor级数展开

4、符号方程的求解

4.1 线性方程组的符号求解

4.2 非线性方程组的符号求解

4.3 常微分方程组的符号求解

4.3.1 dsolve 函数

4.3.2 ode函数

5、简易绘图函数

5.1 二维绘图函数

5.2 三维绘图函数

5.3 等高线绘图函数

5.3.1 函数ezcontour

5.3.2 函数contour

1、符号对象的创建

建立符号变量的方法有两种：sym与syms函数，通常应用sym建立符号表达式，应用syms同时定义多个符号变量。符号变量可以是字符、字符串、表达式或字符表达式。

1.1 sym函数

功能：用来建立单个符号量

语法形式：sym（‘表达式或变量’）

在符号运算中是以表达式形式呈现结果的，而在数值运算中是以数值表示结果的。

举例1：符号标量和数值标量的区别

a = sym('a');     % 定义符号变量a,下同
b = sym('b')
c = sym('c')
d = sym('d')
w = 4;x = 3;y = -8;z = 9;A = [a,b;c,d]     % 建立符号矩阵A
B = [w,x;y,z]     % 建立数值矩阵B
C = det(A)        % 计算符号矩阵A的行列式
D = det(B)        % 建立数值矩阵B的行列式

结果：

A =[ a, b]
[ c, d]B =4     3-8     9C =a*d - b*cD =60

在命令窗口的显示中，数值矩阵只显示元素的数值，而符号矩阵的每行元素放在一对方括号内。

举例2：比较符号常数与数值在代数运算时的差别

r1=8;
r2=2;
r3=3;
k1=sym('8');
k2=sym('2');
k3=sym('3');
sqrt(k1)
sqrt(r1)
sqrt(k3+sqrt(k2))
sqrt(r3+sqrt(r2))

结果：

ans =2*2^(1/2)ans =2.8284ans =(2^(1/2) + 3)^(1/2)ans =2.1010

总结：

①用符号常数进行计算更象进行数学 演算，所得列的结果是精确的数学表达式， 而数值计算将结果近似为一个有限小数；

②可以看到，sym一项只能定义一个符号变 量，使用不方便。

1.2 syms函数

功能：定义多个符号变量

语法：syms Var1 Var2 ……Varn

说明：函数定义符号变量var1，var2，……，Varn等。在用这种格式定义符号变量时无需在变量名上加字符分解符（‘’），变量间用空格而不要用逗号分隔，要用空格来分隔。

举例：在命令窗口创建符号变量a、b、n、x、y和t，建立函数​，然后求f的默认自变量。

>> syms a b n t x y
f=a*x^n+b*y+c*t
findsym(f,6)
ans =x,y,t,n,c,b   ％按最接近x顺序排序
>> findsym(f)
ans =
a, b, n, t, x, y ％按字母排序

2、符号表达式的化简和替换

2.1 创建符号表达式

2.1.1 符号表达式定义：

含有符号变量的表达式。 符号数学工具箱提供了符号表达式的因式分解、展开、合并、化简、通分等操作

2.1.2 建立代数式的方法：

①直接用sym函数建立符号表达式； ②使用已定义的符号变量组成符号表达式

举例：用两种方法建立符号表达式： 3x^2+5y+2xy+6

>> U=sym('3*x^2+5*y+2*x*y+6')
syms x y
V=3*x^2+5*y+2*x*y+6
2*U-V+6U =3*x^2 + 2*y*x + 5*y + 6V =3*x^2 + 2*y*x + 5*y + 6ans =3*x^2 + 2*y*x + 5*y + 12>>

2.1.3 建立符号方程

y=sym('a*x^2+b*x+c＝0')

举例1：计算3阶范得蒙德矩阵行列式的值​

>> syms a b c;
U=[a b c];
A=[[1 1 1];U;U.^2];
det(A)ans =- a^2*b + a^2*c + a*b^2 - a*c^2 - b^2*c + b*c^2>>

举例2：建立x、y的一般二元函数

>> syms x y;
>> f=sym('f(x,y)')f =f(x, y)>>

2.2 符号表达式的化简

2.2.1 因式分解

因式分解函数：factor(s)

举例：对表达式f=x^12-1进行因式分解

>> syms x;
f=x^12-1;
factor(f)ans =[ x - 1, x + 1, x^2 + x + 1, x^2 + 1, x^2 - x + 1, x^4 - x^2 + 1]>>

2.2.2 展开

展开函数：expand(s)

举例：展开(x-1)^12

>> syms x;
f=sym('(x-1)^12');
expand(f)ans =x^12 - 12*x^11 + 66*x^10 - 220*x^9 + 495*x^8 - 792*x^7 + 924*x^6 - 792*x^5 + 495*x^4 - 220*x^3 + 66*x^2 - 12*x + 1

2.2.3 同类项合并

同类项合并函数：collect(s,n) % 按指定变量 n 进行合并

举例：对于表达式f=x(x(x-6)+12)t, 分别将自变量x和t的同类项合并

>> syms x t;
f=x*(x*(x-6)+12)*t;
collect(f)     % 缺省按x合并同类项
collect(f,t)   % 按t合并同类项 ans =t*x^3 - 6*t*x^2 + 12*t*xans =(x*(x*(x - 6) + 12))*t

2.2.4 化简

化简函数：simplify(s)或simple(s)

举例：对表达式(x-1)^3/(x-1)进行化简

>> syms x;f=sym('(x-1)^3/(x-1)');simplify(f)ans =(x - 1)^2>>

函数simple通过对表达式尝试不同的算法进行化简，以寻求符号表达式S的最简形式。

  • y=simple(f): 对 f 尝试多种不同的算法进行简化，返回其中最简短的形式
  • [y,HOW]=simple(f): y 为 f 的最简短形式，How 中记录的为简化过程中使用的方法

综合实验：

>> syms x;s = (-7*x^2-8*y^2)*(-x^2+3*y^2);expand(s)collect(s,x)factor(ans)g = simplify(ans)ans =7*x^4 - 13*x^2*y^2 - 24*y^4ans =7*x^4 - 13*x^2*y^2 - 24*y^4ans =[ 7*x^2 + 8*y^2, x^2 - 3*y^2]g =[ 7*x^2 + 8*y^2, x^2 - 3*y^2]>>

2.2.5 分式通分

分式通分函数：[N,D]=numden(s) %N 为通分后的分子，D 为通分后的分母

举例：对表达式y=x/y+y/x,进行通分

>> syms x y;f=x/y+y/x;[n,d]=numden(f)n =x^2 + y^2d =x*y>>

2.2.6 嵌套

嵌套函数：horner(s)

​

>> syms x;
>> f=x^4+2*x^3+4*x^2+x+1;
>> g=horner(f)g =x*(x*(x*(x + 2) + 4) + 1) + 1>>

2.3 符号表达式的替换

2.3.1 函数subs

  • R=subs(S)：用工作区中的变量值替换符号表达式中的某一特定符号。
  • R=subs(S,New)：用新符号变量New来替换符号表达式S中的默认变量。
  • R=subs(S,Old,New)
>> syms x y
>> f=sym('x^2+x*y+y^2')f =x^2 + x*y + y^2>> x=2;
>> subs(f)ans =y^2 + 2*y + 4>> y=2;
>> subs(f)ans =12>> subs(f,'t^2')ans =t^4 + t^2*y + y^2>> subs(f,{'x','y'},{3,4})ans =37

2.3.2 函数subexpr

  • [Y,SIGMA]=subexpr(S,SIGMA)：指定用变量SIGMA的值（必须为符号对象）来代替符号表达式（可以是矩阵）中重复出现的字符串。替换后的结果由Y返回，被替换的字符串由SIGMA返回。
  • [Y,SIGMA]=subexpr(S,'SIGMA')：与上一种形式的不同在于第二个输入参数是字符或者字符串。
>> syms x a
>> f=sym('2*x^2+a*x+3+x^-1')f =a*x + 1/x + 2*x^2 + 3>> s=solve(f)s =(a^2/36 - 1/2)/(a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3) - a/6 + (a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3)- a/6 - (a^2/36 - 1/2)/(2*(a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3)) - (a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3)/2 + (3^(1/2)*((a^2/36 - 1/2)/(a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3) - (a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3))*i)/2- a/6 - (a^2/36 - 1/2)/(2*(a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3)) - (a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3)/2 - (3^(1/2)*((a^2/36 - 1/2)/(a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3) - (a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4)^(1/3))*i)/2>> r=subexpr(s)sigma = a/8 + ((a^3/216 - a/8 + 1/4)^2 - (a^2/36 - 1/2)^3)^(1/2) - a^3/216 - 1/4r =(a^2/36 - 1/2)/sigma^(1/3) - a/6 + sigma^(1/3)- a/6 + (3^(1/2)*((a^2/36 - 1/2)/sigma^(1/3) - sigma^(1/3))*i)/2 - (a^2/36 - 1/2)/(2*sigma^(1/3)) - sigma^(1/3)/2- a/6 - (3^(1/2)*((a^2/36 - 1/2)/sigma^(1/3) - sigma^(1/3))*i)/2 - (a^2/36 - 1/2)/(2*sigma^(1/3)) - sigma^(1/3)/2

3、符号微积分

3.1 符号积分

3.1.1 limit函数

  • limit(f,x,a): 计算lim(x→a)f(x)
  • limit(f,a): 当默认变量趋向于 a 时的极限
  • limit(f): 计算 a=0 时的极限
  • limit(f,x,a,‘right’): 计算右极限
  • limit(f,x,a,‘left’): 计算左极限

​

>> syms x h n;
>> L=limit((log(x+h)-log(x))/h,h,0)L =1/x>> M=limit((1-x/n)^n,n,inf)M =exp(-x)>>

3.1.2 int函数

  • int(f,a,b)         %对f求定积分，下上限为a,b
  • int(f)               %对f求不定积分
  • int(f,'t')           %用t作为变量求f的不定积分

  • int(f,v,a,b): 计算定积分∫abf(v)dv

  • int(f,a,b): 计算关于默认变量的定积分

  • int(f,v): 计算不定积分∫f(v)dv

>> syms x         % 定义x为变量
y=x.^3+2*x.^2+2;  % 创建函数y=x^3+2*x^2+2
q=int(y,1,2)      % 用int()求积分上下限为2，1的积分q =125/12

举例：求下列不定积分

（1）​​

>> sym('x');
f=(3-x^2)^3;
int(f)%求不定积分ans =- x^7/7 + (9*x^5)/5 - 9*x^3 + 27*x

（2）​

>> sym('x');
f=sqrt(x^3+x^4);
g=int(f)%求不定积分
simplify(g)%化简结果g =(2*x*(x^4 + x^3)^(1/2)*hypergeom([-1/2, 5/2], 7/2, -x))/(5*(x + 1)^(1/2))ans =(2*x*(x^3*(x + 1))^(1/2)*hypergeom([-1/2, 5/2], 7/2, -x))/(5*(x + 1)^(1/2))>>

3.2 符号微分

  • g=diff(f,v)：求符号表达式 f 关于 v 的导数
  • g=diff(f)：求符号表达式 f 关于默认变量的导数
  • g=diff(f,v,n)：求 f 关于 v 的 n 阶导数

solve(f,v)：求方程关于指定自变量的解，f 是符号表达式(低版本还可以用字符串表示的方程、表达式)；solve 也可解方程组(包含非线性)；得不到解析解时，给出数值解。

举例1：分别求表达式f=x^x的导数和三次导数

>> syms x;f=x^x;diff(f)diff(f,3)ans =x*x^(x - 1) + x^x*log(x)ans =log(x)*(x*(x^(x - 2)*(x - 1) + x^(x - 1)*log(x)) + x^x/x + x^(x - 1) + log(x)*(x*x^(x - 1) + x^x*log(x))) - x^x/x^2 + (2*(x*x^(x - 1) + x^x*log(x)))/x + x*((x^(x - 3)*(x - 2) + x^(x - 2)*log(x))*(x - 1) + log(x)*(x^(x - 2)*(x - 1) + x^(x - 1)*log(x)) + x^(x - 2) + x^(x - 1)/x) + 2*x^(x - 2)*(x - 1) + 2*x^(x - 1)*log(x)>>

举例2：求下列函数的导数

​​求​​

>> syms x y;
f=sqrt(1+exp(x));
diff(f)ans =exp(x)/(2*(exp(x) + 1)^(1/2))>>

​​求​​​​

>> syms x y;
f=x*cos(x);
diff(f,x,2)
diff(f,x,3)ans =- 2*sin(x) - x*cos(x)ans =x*sin(x) - 3*cos(x)>>

举例3：求曲线y=x^3＋3x－2上与直线y=4x-1平行的切线的切

>>  x=sym('x');
y=x^3+3*x-2;    % 定义曲线
f=diff(y);      % 对曲线求导
g=f-4;
solve(g)        % 求方程的根ans =-3^(1/2)/33^(1/2)/3>>

举例4：假设水平初速度为2m/s的球，在高为10m的平台上抛出，重力加速度为9.8m/s^2，画出他的抛物线。

clc,clear
syms dt
%假设初速度为2m/s的球，在高为10m的平台上抛出，重力加速度为9.8m/s^2
v_x=2.*dt./dt
v_y=9.8.*dt
dx=v_x.*dt
dy=v_y.*dt.*dt
s=solve(int(dy,0,dt)-10,dt)
time=int(dy,0,1.4214) %可以解出掉地时间。为1.4214秒i=0
while i<1.4214i=i+0.05x=int(dx,0,i)y=10-int(dy,0,i)   hold ongrid onplot(x,y,'o')hold on
end

​

3.3 积分变换

通过积分运算，把一个函数f变成另外一个函数 F，变换过程是​

其中二元函数K(x,t)称为变换的核，变换的核决定了变换的不同名称。f为变换F的原函数，而F(t)成为f(x)的像函数，在一定条件下，两者是一一对应的，可以相互转化。

解微分方程思路：①不易直接求得解f，则对原方程进行变换； ②能从变换后的方程中求得解f，则对F进行逆变换， 当然，在选择变换的核时，应使得变换后的方程比 原方程容易求解。

常见的变换：傅立叶变换、拉氏变换、Z变换

3.3.1 Fourier变换

Fourier变换函数：

fourier(fx,x,t)  傅里叶变换ifourier(Fw,t,x) 傅里叶反变换

举例：求 y=|x|的傅立叶变换及其逆变换

>> syms x t;
y=abs(x);
Ft=fourier(y,x,t)%求y的傅立叶变换
fx=ifourier(Ft,t,x)%求Ft的傅氏逆变换Ft =-2/t^2fx =x*sign(x)>>

傅立叶变换是一种常见的分析方法，傅立叶变换将满足一定条件的函数表示为一些函数的加权和（或者积分）。可以分为四个类别：
1. 非周期连续性信号
对应于傅里叶变换，频域连续非周期
2. 周期性连续性信号
对应于傅立叶级数，频域离散非周期
3. 非周期离散信号
对应于DTFT（离散时间傅立叶变换），频域连续周期
4. 周期性离散信号
对应于DFT（离散时间傅立叶变换），频域离散周期

  • 傅里叶级数

傅立叶级数是将一个信号在一组正交基上进行分解的体现

​

  • 连续时间傅立叶变换

​

  • 离散时间傅立叶变换

离散时间傅立叶变换在频域上是连续的，但由于计算机无法表示无限长的时间片段，已经无法表示全部频率，一般取一定频域的分量。

​

  •  离散傅立叶变换

只有离散傅立叶变换在频域和时域都是离散的，即计算机可以处理的，因此DFT是可以实际进行编程并实用的。DFT的信号首先要进行截断，因为能处理的信号必须是有限的；然后对信号进行采样，对频谱进行离散化。

​

  •  二维傅立叶变换

​

  • 傅立叶变换实现

只有离散傅里叶变换才可以实现，在matlab中实现有fft，fft2进行傅里叶变换，同样可以手动进行变换。

一维傅立叶变换

基于FFT

%  xn是信号，n是坐标，N是点数
%  N =8;
%  n = [0:1:N-1];
%  xn = 0.5.^n;        % 指数信号
function [] = DFTusefft(xn,n,N)figure(1);Xk=fft(xn,N);      % 傅立叶变换subplot(211);stem(n,xn);title('原信号');subplot(212);stem(n,abs(Xk));title('FFT变换')
end

​

DFT公式

function [] = DFT(xn,n,N)Xk = zeros(1,N);    for k=1:Nsn =0.0;for i=1:Nsn = sn+xn(i)*exp(-j*2*pi*i*k/N);endXk(k) = sn;endfigure(2);subplot(211);stem(n,xn);title('原信号');subplot(212);stem(n,abs(Xk));title('DFT')
end

​

DTFT
由于DTFT的频域是连续的而且是无穷的，当我们选择的最高频域足够高时，可以基本代表信号特征，可以进行编程。

function [] = testDTFT(xn,n,N)figure(3);w=[-800:1:800]*4*pi/800;     %频域共-800----+800 的长度（本应是无穷，高频分量很少，故省去）    w = [-N/2:1:N/2]*4*pi*2/N;X=xn*exp(-j*(n'*w));         %求dtft变换，采用原始定义的方法，对复指数分量求和而得subplot(211)stem(n,xn);title('原始信号(指数信号)');subplot(212);plot(w/pi,abs(X));title('DTFT变换')
end

​

 二维傅立叶变换

clear
clc
img=imread('peppers.png');
subplot(2,2,1);imshow(img);title('原图');
f=rgb2gray(img);    %对于RGB图像必须做的一步，也可以用im2double函数
F=fft2(f);          %傅里叶变换
F1=log(abs(F)+1);   %取模并进行缩放
subplot(2,2,2);imshow(F1,[]);title('傅里叶变换频谱图');
Fs=fftshift(F);      %将频谱图中零频率成分移动至频谱图中心
S=log(abs(Fs)+1);    %取模并进行缩放
subplot(2,2,3);imshow(S,[]);title('频移后的频谱图');
fr=real(ifft2(ifftshift(Fs)));  %频率域反变换到空间域，并取实部
ret=im2uint8(mat2gray(fr));    %更改图像类型
subplot(2,2,4);imshow(ret),title('逆傅里叶变换');

​

3.3.2 Laplace变换

我们通常把f(t)的拉普拉斯变换写成F(s)，拉普拉斯变换的优点是能把微分方程写成代数方程，而求解代数方程则要简单的多。

Laplace变换函数：

laplace(fx,x,t) 拉普拉斯变换laplace(Ft,t,x) 反拉普拉斯变换

举例：计算y=x^2的拉普拉斯变换及其逆变换

>> x=sym('x');
y=x^2;
Ft=laplace(y,x,t)%对函数y进行拉氏变换
fx=ilaplace(Ft,t,x)%对Ft进行拉氏逆变换Ft =2/t^3fx =x^2>>
  • 幂函数

公式：​

二次方：

>> syms t;
laplace(t^2)ans =2/s^3>>

三次方：

>> syms t;
laplace(t^3)ans =6/s^4>>
  • 指数函数
>> syms t b;
laplace(exp(-b*t))ans =1/(b + s)>>
  • 正弦、余弦、双曲正弦、双曲余弦函数
>> syms w t;
laplace(sin(w*t))ans =w/(s^2 + w^2)>> laplace(cos(w*t))ans =s/(s^2 + w^2)>> laplace(sinh(w*t))ans =w/(s^2 - w^2)>> laplace(cosh(w*t))ans =s/(s^2 - w^2)>>

3.3.3 Z变换

Z变换函数：

ztrans(fn,n,z)  Z变换iztrans(Fz,z,n) 逆Z变换

举例：求数列fn=e^-n的Z变换及其逆变换

>> syms x z;
fn=exp(-n);
Fz=ztrans(fn,n,z)%对fn的Z变换
f=iztrans(Fz,z,n)%对Fz的逆Z变换Fz =z/(z - exp(-1))f =kroneckerDelta(n, 0) - exp(-1)*(exp(1)*kroneckerDelta(n, 0) - exp(-1)^n*exp(1))>>

3.4 符号求和

函数：symsum

功能：级数符号求和

语法：

symsum(S)

symsum(S,V)

symsum(S,a,b)

说明:函数symsum(S)中S为符号表达式,S相对于符号变量k的和,k取值从0到k-1.函数symsum(S,V)中指定S相对于变量V的和,V从0变到V-1。函数symsum(S,a,b)和symsum(s,v,a,b)指定符号表达式从v=a累加到v=b。

举例：求下列级数的和

​

​

​

​

>> n=sym('n');x=sym('x');
s1=symsum(1/n^2,n,1,inf)
s2=symsum((-1)^(n+1)/n,1,inf)
s3=symsum(n*x^n,n,1,inf)
s4=symsum(n^2,1,100)s1 =pi^2/6s2 =log(2)s3 =piecewise([abs(x) < 1, x/(x - 1)^2])s4 =338350>>

3.5 Taylor级数展开

函数：taylor

函数调用格式：

taylor(fcn,x,x0,'Order',6);%对函数fcn在点x0处，进行6阶泰勒展开；

举例：对y=exp(-x)进行4阶泰勒展开

>> syms x
y=exp(-x);
y1=taylor(y,x,2,'Order',4)y1 =exp(-2) - exp(-2)*(x - 2) + (exp(-2)*(x - 2)^2)/2 - (exp(-2)*(x - 2)^3)/6>>

4、符号方程的求解

4.1 线性方程组的符号求解

函数：linsolve(A,b)

举例：求线性方程组Ax=b的解

​

>> A=[34 8 4;3 34 3;3 6 8];
b=[4;6;2];
x=linsolve(A,b)
A\bx =0.06750.16140.1037ans =0.06750.16140.1037>>

4.2 非线性方程组的符号求解

函数：

y=solve(eq)

y=solve(eq,var)

y=solve(eq1,eq2,..var1,var2,..)

举例：解下列方程

​

​

​

​

x=solve('1/(x+2)+4x/(x^2-4)=1+2/(x-2)','x');
x=solve('x-(x^3-4*x-7)^(1/3)=1','x')
x=solve('2*sin(3*x-pi/4)=1')
x=solve('x+x*exp(x)-10','x')

4.3 常微分方程组的符号求解

4.3.1 dsolve 函数

dsolve函数用于求常微分方程组的精确解，也称为常微分方程的符号解。如果没有初始条件或边界条件，则求出通解；如果有，则求出特解。

D: 微分符号；D2表示二阶微分，D3表示三阶微分，以此类推。

函数格式：

Y = dsolve('eq1,eq2,…' , 'cond1,cond2,…' , 'Name')

其中，‘eq1,eq2,…’:表示微分方程或微分方程组;

’cond1,cond2,…’:表示初始条件或边界条件;

‘Name’:表示变量。没有指定变量时，matlab默认的变量为t；

举例1：求下列微分方程的通解

​

>> y1=dsolve('Dy-(x^2+y^2)/x^2/2','x')%解（1）y2=dsolve('Dy*x^2+2*x*y-exp(x)','x')%解（2）y1 =x-x*(1/(C6 + log(x)/2) - 1)y2 =-(C9 - exp(x))/x^2>>

举例2：求微分方程的特解

​

>> y1=dsolve('Dy=2*x*y^2','y(0)=1','x')%解（1）
y2=dsolve('Dy-x^2/(1+y^2)','y(2)=1','x')%解（2）y1 =-1/(x^2 - 1)y2 =(((x^3/2 - 2)^2 + 1)^(1/2) + x^3/2 - 2)^(1/3) - 1/(((x^3/2 - 2)^2 + 1)^(1/2) + x^3/2 - 2)^(1/3)>>

举例3：用微分方程的数值解法和符号解法解下列方程，并对结果进行比较。

​

>> dsolve('Dy+2*y/x-4*x','y(1)=2','x') %用符号方法得到方程的解析解ans =1/x^2 + x^2>>

为了求方程的数值解，需按要求建立一个函数文件fxyy.m

function f=fxyy(x,y)f=(4*x^2-2*y)/x %只能是y’=f(x,y)的形式return
[t,w]=ode45('fxyy',[1 2],2);

为了对两种结果进行比较，在同一个坐标系中 作出两种结果的图形。

x=linspace(1,2,100);y=x.^2+1./x.^2;%为作图把符号解的结果离散化plot(x,y,'b',t,w,'r-');%符号解图形用兰色圆点表示，数值解图形用红色实线表示。

举例4：求下列微分方程

​

>> [x,y]=dsolve('Dx=4*x-2*y','Dy=2*x-y','t')
[x,y]=dsolve('D2x-y','D2y+x','t')x =C22/2 + 2*C21*exp(3*t)y =C22 + C21*exp(3*t)x =2^(1/2)*C24*exp((2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 + sin((2^(1/2)*t)/2)/2) - 2^(1/2)*C23*exp((2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 - sin((2^(1/2)*t)/2)/2) + 2^(1/2)*C25*exp(-(2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 + sin((2^(1/2)*t)/2)/2) + 2^(1/2)*C26*exp(-(2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 - sin((2^(1/2)*t)/2)/2)y =2^(1/2)*C23*exp((2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 + sin((2^(1/2)*t)/2)/2) + 2^(1/2)*C24*exp((2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 - sin((2^(1/2)*t)/2)/2) - 2^(1/2)*C25*exp(-(2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 - sin((2^(1/2)*t)/2)/2) + 2^(1/2)*C26*exp(-(2^(1/2)*t)/2)*(cos((2^(1/2)*t)/2)/2 + sin((2^(1/2)*t)/2)/2)>>

4.3.2 ode函数

ode是Matlab专门用于解微分方程的功能函数。该求解器有变步长（variable-step）和定步长（fixed-step）两种类型。不同类型有着不同的求解器

ode45

一步算法：4、5阶龙格库塔方程：累计截断误差(Δx)^5；大部分尝试的首选算法

ode23

一步算法：2、3阶龙格库塔方程：累计截断误差(Δx)^3；适用于精度较低的情形

ode113

多步算法：Adams；计算时间比ode45短；刚性ode求解命令

ode23t

梯形算法；适度刚性情形

ode15s

多步法：Gear’s反向数值微分：精度中等；若ode45失效时，可以尝试使用

ode23s

一步法：2阶Rosebrock算法：精度低；当精度较低时，计算时间比ode15s短

ode23tb

梯形算法：精度低；当精度较低时，计算时间比ode15s短

下面以ode45为例具体介绍函数的使用方法：

函数格式：  

[T,Y] = ode45(‘odefun’,tspan,y0)

[T,Y] = ode45(‘odefun’,tspan,y0,options)

[T,Y,TE,YE,IE] = ode45(‘odefun’,tspan,y0,options)

sol = ode45(‘odefun’,[t0 tf],y0...)

  • odefun是函数句柄,可以是函数文件名,匿名函数句柄或内联函数名；
  • tspan 是求解区间 [t0 tf]，或者一系列散点[t0,t1,...,tf]；
  • y0 是初始值向量
  • T 返回列向量的时间点
  • Y 返回对应T的求解列向量
  • options 是求解参数设置,可以用odeset在计算前设定误差,输出参数,事件等
  • TE 事件发生时间
  • YE 事件发生时之答案
  • IE 事件函数消失时之指针i

5、简易绘图函数

5.1 二维绘图函数

matlab提供了一个ezplot函数绘制隐函数图形，有三种调用方式：

  1. 对于函数f=f(x)，ezplot函数的调用格式为：

    • ezplot(f)：在默认区间-2pi<=x<=2pi内绘制f=f(x)的图形
    • ezplot(f,[a,b])：在区间a<=x<=b内绘制f=f(x)的图形

  2. 对于隐函数f=f(x,y)，ezplot函数的调用格式为：

    • ezplot(f)：在默认区间-2pi<=x<=2pi和-2pi<=y<=2pi内绘制f(x,y)=0的图形
    • ezplot(f,[x_min,x_max,y_min,y_max])：在区间x_min<=x<=x_max和y_min<=y<=y_max内绘制f(x,y)=0的图形
    • ezplot(f,[a,b])：在区间a<=x<=b和a<=y<=b内绘制f(x,y)=0的图形

  3. 对于参数方程x=x(t)和y=y(t)，ezplot函数的调用格式为：

    • ezplot(x,y)：在默认区间0<=t<=2pi内绘制x=x(t)和y=y(t)的图形
    • ezplot(x,y,[t_min,t_max])：在区间t_min<=t<=t_max绘制x=x(t)和y=y(t)的图形

举例： 在极坐标下绘制函数表达式1+cost的二维图形

syms t
ezpolar(1+cos(t))

​

5.2 三维绘图函数

函数格式：

  • ezplot3(x,y,z) ：在默认区间0
  • ezplot3(x,y,z,[tmin,tmax])：在tmin < t < t max上绘制曲线 x = x(t), y = y(t), and z = z(t)
  • ezplot3(...,'animate')：产生一个空间曲线的动态轨迹

举例：绘制参变量函数：x=sin(t) ，y=cos(t)，z=t  区间：[0,6pi]

ezplot3('sin(t)','cos(t)','t',[0,6*pi])

​

5.3 等高线绘图函数

5.3.1 函数ezcontour

函数格式：

ezcontour(fun) 使用contour函数绘制 fun(x,y) 的等高线。fun 在默认域 -2π < x < 2π 和 -2π < y < 2π 上绘制。

ezcontour(fun,domain) 在指定的 domain 中绘制 fun(x,y)。domain 可以是 4×1 向量 [xmin, xmax, ymin, ymax] 或 2×1 向量 [min, max]（其中 min < x < max，min < y < max）。

ezcontour(...,n) 使用 n×n 网格在默认域中绘制 fun。n 的默认值为 60。

ezcontour(axes_handle,...) 将图形绘制到带有句柄 axes_handle 的坐标区中，而不是当前坐标区 中。

h = ezcontour(...) 将句柄返回给 h 中的等高线对象。

举例：绘制下列函数的等高线

​

将函数句柄 f 传递给 ezcontour。在 x 方向和 y 方向上同时指定一个从 -3 到 3 的域，并使用 49×49 计算网格。

f = @(x,y) 3*(1-x).^2.*exp(-(x.^2) - (y+1).^2) ...- 10*(x/5 - x.^3 - y.^5).*exp(-x.^2-y.^2) ...- 1/3*exp(-(x+1).^2 - y.^2);ezcontour(f,[-3,3],49)

​

5.3.2 函数contour

函数格式：

contour(Z) 绘制矩阵 Z 的等高线图，其中 Z 解释为有关 x-y 平面的高度。Z 必须至少是 2×2 矩阵，该矩阵包含至少两个不同值。x 值对应于 Z 的列索引，y 值对应于 Z 的行索引。自动选择等高线层级。

contour(Z,n) 以 n 个等高线层级绘制矩阵 Z 的等高线图，其中 n 为标量。自动选择等高线层级。

contour(Z,v) 绘制矩阵 Z 的等高线图，其中等高线位于单调递增向量 v 中指定的数据值。要在特定值位置显示单个等高线，请将 v 定义为一个二元素向量，并且两个元素都等于所需的等高线层级。例如，要在 k 层级绘制等高线，请使用 contour(Z,[k k])。指定向量 v 可将LevelListMode属性设置为 manual。

contour(X,Y,Z)、contour(X,Y,Z,n) 和 contour(X,Y,Z,v) 使用 X 和 Y 绘制 Z 的等高线图，以确定 x 和 y 值。

  • 如果 X 和 Y 为向量，则 length(X) 必须等于 size(Z,2) 且 length(Y) 必须等于 size(Z,1)。这些向量必须是严格递增或严格递减的，并且不能包含任何重复值。

  • 如果 X 和 Y 为矩阵，则其大小必须等于 Z 的大小。通常，应设置 X 和 Y 以使列严格递增或严格递减并且行是均匀的（或者使行严格递增或严格递减并且列是均匀的）。

contour(...,LineSpec) 使用LineSpec指定的线型和颜色绘制等高线。contour 忽略标记符号。

contour(...,Name,Value) 使用一个或多个属性名-属性值对组指定等高线属性。Name 是属性名且必须置于单引号 ('') 内。Value 是对应的值。例如，'LineWidth',2 将等高线宽度设置为 2。

contour(ax,...) 将图形绘制到 ax 指定的坐标区中，而不是当前坐标区 (gca) 中。

[C,h] = contour(...) 返回等高线矩阵 C（包含定义等高线的数据）和 Contour 对象 h。

举例1：使用 meshgrid 函数生成矩阵 X 和 Y。创建第三个矩阵 Z 并绘制其等高线。

x = linspace(-2*pi,2*pi);
y = linspace(0,4*pi);
[X,Y] = meshgrid(x,y);
Z = sin(X)+cos(Y);figure
contour(X,Y,Z)

​

举例2：将来自 peaks 函数的数据存储于矩阵 X、Y 和 Z 中。对 Z 中的数据绘制 20 条等高线。

[X,Y,Z] = peaks;
figure
contour(X,Y,Z,20)

​

举例3：设置矩阵 X、Y 和 Z。通过将 ShowText 属性设置为 on 来创建一个等高线图并显示等高线标签。

x = -2:0.2:2;
y = -2:0.2:3;
[X,Y] = meshgrid(x,y);
Z = X.*exp(-X.^2-Y.^2);figure
contour(X,Y,Z,'ShowText','on')

​

参考资料：

傅里叶变换及其实现（MATLAB）_新安浅滩的博客-CSDN博客_matlab傅里叶变换

（不推荐）易用的等高线绘图函数 - MATLAB ezcontour- MathWorks 中国matlab从无到有系列（三）：数值计算基础

matlab从无到有系列（四）：符号数学基础相关推荐

  1. matlab从无到有系列（七）：GUI程序设计

  2. Matlab系列之符号运算（上）

    Matlab系列之符号运算 前言 创建符号对象 基本操作 符号变量的基本操作 符号表达式的基本操作 四则运算 多项式的操作 符号表达式化简 符号表达式的替换 反函数求解 复合函数 更多精彩等你发现~ ...

  3. matlab 积分进阶教程,最适合小白的matlab教程系列_进阶系列二之微积分

    微积分篇 更多文章参考 符号变量 极限问题 求导数 求导数----偏导数 积分--符号解 积分数值解 多重积分 级数求和 泰勒级数展开 常微分方程 更多文章参考 最适合小白的matlab教程系列_基础 ...

  4. 最适合小白的matlab教程系列_进阶系列二之微积分

    微积分篇 更多文章参考 符号变量 极限问题 求导数 求导数----偏导数 积分--符号解 积分数值解 多重积分 级数求和 泰勒级数展开 常微分方程 更多文章参考 最适合小白的matlab教程系列_基础 ...

  5. sed修炼系列(四)：sed中的疑难杂症

    sed系列文章: sed修炼系列(一):花拳绣腿之入门篇 sed修炼系列(二):武功心法(info sed翻译+注解) sed修炼系列(三):sed高级应用之实现窗口滑动技术 sed修炼系列(四):s ...

  6. 机器学习入门系列四（关键词：BP神经网络）

    机器学习入门系列四(关键词:BP神经网络) 标签: 机器学习神经网络 2016-01-12 15:28 80人阅读 评论(0) 收藏 举报 本文章已收录于: 分类: 机器学习(3) 作者同类文章X 版 ...

  7. java面试常考系列四

    转载自 java面试常考系列四 题目一 大O符号(big-O notation)的作用是什么?有哪些使用方法? 大O符号描述了当数据结构里面的元素增加的时候,算法的规模或者是性能在最坏的场景下有多么好 ...

  8. [矩阵的QR分解系列四] QR(正交三角)分解

    QR分解 简介 QR分解 定义 存在和唯一性 存在性证明 唯一性证明 分解方法 施密特(Schmidt)方法 吉文斯(Givens)方法 豪斯霍尔德(Householder)方法 例子 施密特(Sch ...

  9. 效率系列(四) VS常用快捷键

    写在前面的话 :最近博主整理了一些关于 Visual Studio 2017 的常用快捷键,希望可以帮助到大家更高效更愉快的打码 1.打开 快捷键 描述 Ctrl + Shift + N 新建项目 C ...

最新文章

  1. ActiveMQ高级特性
  2. 使用Quartus建立第一个FPGA工程
  3. 1.html5究竟是什么
  4. 模拟银行取钱 (多线程操作的好处)
  5. oracle迁移父子数据
  6. How to create a hyperlink in SQL Server Reporting Services
  7. C语言变量unsafe,Java中Unsafe使用详解
  8. 网易云自动化爬虫方案，24小时爬取10万首音乐信息，30万条用户信息，100万个热评
  9. 华为薪资等级结构表_华为内部考核体系
  10. 如何学习硬件设计——实践篇
  11. kettle-java代码执行hive相关ktr时报错： database type with plugin id [HIVE2] couldn‘t be found!
  12. TIPTOP报表提示：foreach: 字元轉換至數字程序失敗.
  13. require.context()的用法
  14. 系统滴答定时器的应用
  15. pma连接，报错10061
  16. 修改Echarts源码实现柱状图的炫彩闪烁效果
  17. 今天，Java 12 正式发布了！ 你在用哪个版本？
  18. 通信算法之九十七：4G LTE通信系统 PHICH信道处理流程
  19. 防止小孩长时间沉迷电脑游戏-卓越电脑定时关机软件
  20. 离散数学及其应用_第1章_基础：逻辑和证明

热门文章

  1. 将本地仓库推送到远程SSHHTTPS
  2. 大学Java基础课程设计——网络聊天室
  3. 子域名收集原理与子域名爆破工具
  4. Ubuntu 10.04 更新源(ubuntu yuan)
  5. ZYNQ 图像处理之千兆网传（一）【寄存器级操作】
  6. win7桌面便笺：启动时提示windows7部分便笺的元数据已被损坏
  7. 西行漫记（6）：Diversity
  8. ITE平台开发 chapter5之jperf测试网络吞吐量
  9. sql modify的使用方法
  10. java个十百千万位余数_第一学期四年级数学期末试题