精华内容
下载资源
问答
  • 如何利用matlab求解方程

    万次阅读 多人点赞 2018-04-22 12:46:54
    如何利用matlab求解方程1. 前言作为三大数学软件之一,matlab在数值计算方法的能力首屈一指。求解方程是工科学习和工程计算中最基础、最常见的问题。掌握利用现代化工具求解方程的方法对于提升我们的工科素养至关...

    如何利用matlab求解方程

    1.    前言

    作为三大数学软件之一,matlab在数值计算方法的能力首屈一指。求解方程是工科学习和工程计算中最基础、最常见的问题。掌握利用现代化工具求解方程的方法对于提升我们的工科素养至关重要。为此,本篇将对matlab中求解方程的方法进行介绍。

    2.    用法

    求解过程

    2.1    指明变量

    告诉电脑方程中所含有的变量,包括参数和未知变量。比如:所求解的方程为:,很显然该方程中有a,b,c,x符号变量,因此该步骤的写法为:

    syms a b c x  

     

    2.2    指明方程,未知数和限制条件(非必需)

    eqns

    方程,如果超过一个,则放在[ ]中,并用逗号隔开。如:

    Vars

    待求的未知数

    Names-value(非必需)

    Names:‘returnConditions’

    是否返回出含有参数的通解。’true‘为返回,’false‘为否,即给出一个特解;

    Name: 'IgnoreAnalyticConstraints'

    是否给出解的最简形式。

    ‘true‘为是,‘false’为否

    Name:'PrincipalValue'

    是否仅给出一个解。False为返回所有的解,true为仅返回一个解;

    Name:’Real’

    是否仅返回实数解

    2.3    获得所求方程的解

    如果为多个函数,该解存储形式为结构体。

     

    3.    具体实例

    3.1   求解sin(x)=1的通解

     

    具体代码:

    syms x 

    [x,params,conds]=solve(sin(x)==1,'ReturnConditions', true) 

     

    结果

    solx =pi/2+2*pi*

    params =

    conds =in(k,'integer')

     

    可以看出,该方程的通解为:

     

    3.2   求解以下方程:

    代码:

    syms a b c y x

    [x,y]=solve([a*x^2+b*y+c==0,a*x+2*y==4],[x,y])

     

    结果:

    x =

     ((a*b)/4-(-(a*(- a*b^2+32*b +16*c))/16)^(1/2))/a

     ((a*b)/4+(-(a*(- a*b^2+32*b +16*c))/16)^(1/2))/a

     

    y =

     (-(a*(- a*b^2+32*b +16*c))/16)^(1/2)/2-(a*b)/8+2

     2-(-(a*(- a*b^2+32*b +16*c))/16)^(1/2)/2-(a*b)/8

     

    即:




    展开全文
  • 本文是科学计算与MATLAB语言专题六第2小节的学习笔记,如果大家有时间的话,可以去听听课,没有的话,可以看看我的笔记,还是很不错的。 1 直接法 1.线性方程组的直接解法 高斯(Gauss)消去法 列主元消去法 矩阵的...

    前言

    今天我们要说的就是数值微积分,赶紧看看他和高等数学中的微积分有什么区别吧。本文是科学计算与MATLAB语言专题六第2小节的学习笔记,如果大家有时间的话,可以去听听课,没有的话,可以看看我的笔记,还是很不错的。

    1 直接法

    1.线性方程组的直接解法
    高斯(Gauss)消去法
    列主元消去法
    矩阵的三角分解法
    高斯(Gauss)消去法是一个经典的直接法,由它改进得到的列主元消去法,是目前计算机上求解线性方程组的标准算法,其特点就是通过消元将一般线性方程组的求解问题转化为三角方程组的求解问题。此外,还有矩阵的三角分解法等许多直接求解算法。
    (1)利用左除运算符的直接解法
    MATLAB提供了一个左除运算符“\”用于求解线性方程组,它使用列主元消去法,使用起来十分方便。对于线性方程组Ax=bAx=b,可以利用左除运算符反斜杠求解,b左除以A可获得线性方程组的数值解x。
    Ax=bAx=bx=A\bx=A\backslash b
    注意,如果矩阵A是奇异的或接近奇异的,则MATLAB会给出警告信息。
    例1 用左除运算符求解下列线性方程组。
    {2x1+x25x3+x4=13x15x2+7x4=92x2+x3x4=6x1+6x2x34x4=0\left\{ \begin{aligned} 2x_1+x_2-5x_3+x_4=13\\ x_1-5x_2+7x_4=-9\\ 2x_2+x_3-x_4=6\\ x_1+6x_2-x_3-4x_4=0 \end{aligned} \right.

    A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];
    b=[13,-9,6,0]';
    x=A\b
    

    (2)利用矩阵分解求解线性方程组矩阵分解是设计算法的重要技巧,是指将一个给定的矩阵分解成若干个特殊类型矩阵的乘积,从而将一个一般的矩阵计算问题转化为几个易求的特殊矩阵的计算问题。通过矩阵分解方法求解线性方程组的优点是运算速度快,可以节省存储空间。
    LU分解
    QR分解
    Cholesky分解
    ①LU分解的基本思想
    矩阵的LU分解就是将一个n阶矩阵表示为一个下三角矩阵和一个上三角矩阵的乘积。线性代数中已经证明,只要方阵是非奇异的,LU分解总是可以进行的。
    A=LUAx=bLUx=bA=LU →Ax=b→LUx=b{Ly=bUx=y\left\{ \begin{aligned}Ly=b\\Ux=y \end{aligned} \right.
    对于三角方程很容易求解,于是可以首解向量y使Ly=b,再求解Ux=y,从而达到求解线性方程组Ax=b的目的。
    ②MATLAB的LU分解函数
    LU分解函数是根据列主元LU分解算法定义的,具有较好的数据稳定性。lu函数有两种调用格式:
    [L,U]=lu(A):
    产生一个上三角阵U和一个变换形式的下三角阵L,使之满足A=LU。注意,这里的矩阵A必须是方阵。
    [L,U,P]=lu(A):
    产生一个上三角阵U和一个下三角阵L以及一个置换矩阵P,使之满足PA=LU。同样,矩阵A必须是方阵。
    当使用第一种格式时,矩阵L往往不是一个下三角阵,但可以通过行交换成为一个下三角阵。
    ③用LU分解求解线性方程组
    Ax=bAx=bLUx=bLUx=bx=U\(L\bx=U \backslash (L\backslash b)

    Ax=bAx=bPAx=PbPAx=PbLUx=PbLUx=Pbx=U\Lbx=U \backslash(L*b)
    通过LU分解后可以大大提高运算速度。
    例2 用LU分解求解例1中的线性方程组。

    A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];
    b=[13,-9,6,0]';
    [L,U]=lu(A);
    x=U\(L\b)
    
    x =
      -66.5556
       25.6667
      -18.7778
       26.5556
    

    2 迭代法

    2.线性方程组的迭代解法
    迭代法是一种不断用变量的原值推出它的新值的过程,是用计算机解决问题的一种基本方法。
    {10x1x2=9x1+10x22x3=72x2+10x3=6\left\{ \begin{aligned} 10x_1-x_2=9\\ -x_1+10x_2-2x_3=7\\ -2x_2+10x_3=6 \end{aligned} \right.    \implies {x2=10x19x1=10x22x37x3=(6+2x2)10\left\{ \begin{aligned} x_2&=10x_1-9\\ x_1&=10x_2-2x_3-7\\ x_3&=(6+2x_2)10 \end{aligned} \right.    \implies{x1k+1=10x2k2x3k7x2k+1=10x1k9x3k+1=0.6+0.2x2k\left\{ \begin{aligned} x_1^{k+1}&=10x_2^{k}-2x_3^{k}-7\\ x_2^{k+1}&=10x_1^{k}-9\\ x_3^{k+1}&=0.6+0.2x_2^{k} \end{aligned} \right.
    (1)雅可比(Jacobi)迭代法
    Ax=b,A=DLU(DLU)x=bAx=b, A=D-L-U, (D-L-U)x=b
    D=[a11a22ann]D=\left [ \begin{array}{cccc} a_{11} & & & \\ & a_{22} & & \\ & & \ddots & \\ & & & a_{nn} \end{array} \right]
    L=[0a210a31a320an1an20]L=-\left [ \begin{array}{cccc} 0 & && & \\ a_{21} & 0 && & \\ a_{31} & a_{32} &0& & \\ \vdots &\vdots&\vdots &\ddots&\\ a_{n1} & a_{n2}&\cdots& \cdots&0 \end{array} \right]
    U=[0a12a13a1n0a23a2nan1n0]U=-\left [ \begin{array}{cccc} 0 & a_{12} &a_{13}& \cdots & a_{1n} \\ & 0 &a_{23}& \cdots &a_{2n} \\ & &\ddots& \ddots & \vdots \\ & & &\ddots&a_{n-1n}\\ &&& &0 \end{array} \right]

    求解公式为:
    x=D1L+Ux+D1bx=D^{-1}(L+U)x+D^{-1}b
    与之对应的迭代公式为:
    xk+1=D1L+Uxk+D1bx^{k+1}=D^{-1}(L+U)x^{k}+D^{-1}b     (B=D1(L+U),f=D1b\implies(B=D^{-1}(L+U),f=D^{-1}b)     xk+1=Bxk+f\implies x^{k+1}=Bx^{k}+f
    雅可比迭代法的函数文件jacobi.m:

    function [y,n]=jacobi(A,b,x0,ep)
    D=diag(diag(A));
    L=-tril(A,-1);
    U=-triu(A,1);
    B=D\(L+U);
    f=D\b;
    y=B*x0+f;
    n=1;
    while norm(y-x0)>=ep
        x0=y;
        y=B*x0+f;
        n=n+1;
    end
    

    (2)高斯-赛德尔(Gauss-Serdel)迭代法
    Dxk+1=L+Uxk+bDxk+1=Lxk+1+Uxk+bDLxk+1=Uxk+bxk+1=DL1U×0+UL1bB=DL1U,f=DL1bxk+1=Bxk+fDx^{k+1}=(L+U)x^{k}+b →Dx^{k+1}=Lx^{k+1}+Ux^k+b→ (D-L)x^{k+1}=Ux^{k}+b \\ x^{k+1}=(D-L)^{-1}U×0+(U-L)^{-1}b →B=D-L^{-1}U,f=D-L^{-1}b→ x^{k+1}=Bx^{k}+f
    Gauss-Serdel迭代法的函数文件gauseidel.m

    function [y,n]=gauseidel(A,b,x0,ep)
    D=diag(diag(A));
    L=-tril(A,-1); 
    U=-triu(A,1);
    B=(D-L)\U;
    f=(D-L)\b;
    y=B*x0+f;
    n=1;
    while norm(y-x0)>=ep
        x0=y;
        y=B*x0+f;
        n=n+1;
    end
    
    
    

    例3 分别用雅可比迭代法和高斯-赛德尔迭代法求解线性方程组。设迭代初值为0,迭代精度为10^(-6)。
    [421243133]\left[\begin{matrix} 4&-2&-1\\ -2&4&3\\ -1&-3&3 \end{matrix}\right][x1x2x3]\left[\begin{matrix} x_1\\ x_2\\ x_3 \end{matrix}\right]=[150]\left[\begin{matrix} 1\\ 5\\ 0 \end{matrix}\right]

    A=[4,-2,-1;-2,4,3;-1,-3,3];
    b=[1,5,0]';
    [x,n]=jacobi(A,b,[0,0,0]',1.0e-6)
    [x,n]=gauseidel(A,b,[0,0,0]',1.0e-6)
    
    x =
        0.9706
        0.8529
        1.1765
    n =
        35
    x =
        0.9706
        0.8529
        1.1765
    n =
        16
    

    例4 分别用雅可比迭代法和高斯-赛德尔迭代法求解下列线性方程组,看是否收敛。
    [122111221]\left[\begin{matrix} 1&-2&2\\ 1&1&1\\ 2&2&1 \end{matrix}\right][x1x2x3]\left[\begin{matrix} x_1\\ x_2\\ x_3 \end{matrix}\right]=[976]\left[\begin{matrix} 9\\ 7\\ 6 \end{matrix}\right]

    A=[1,2,-2;1,1,1;2,2,1];
    b=[9;7;6];
    [x,n]=jacobi(A,b,[0;0;0],1.0e-6)
    [x,n]=gauseidel(A,b,[0;0;0],1.0e-6)
    
    x =
       -27
        26
         8
    n =
         4
    x =
       NaN
       NaN
       NaN
    n =
        1012
    

    小结

    直接法:以矩阵初等变换为基础,可以求得方程组的精确解;占用的内存空间大、程序实现较为复杂;一般适合求解低阶稠密线性方程组。
    迭代法:从给定初始值逐步逼近精确解的过程,求解过程占用存储空间小,程序设计简单;适用于求解大型稀疏矩阵线性方程组;要考虑算法的收敛性。
    大家学会了吗?另外分享给大家一个Markdown的公式神器- Katex\ Katex ,为了编辑出优美的公式,大家可以多看看官方支持文档.最后没有一键三连,但欢迎大家
    点赞👍,收藏⭐,转发🚀,
    如有问题、建议,请您在评论区留言💬哦。

    展开全文
  • MATLAB求解微分方程

    2021-05-23 20:48:31
    下面介绍如何Matlab 来计算微分方程(组)的数值解。 Euler折线法 微分方程的基本数值解法——Euler折线法。 步骤: 分割求解区间,差商代替微商,解代数方程 例子: 解: 等距剖分:a=x0<x1<x2<⋯&l

    求微分方程的解

    自牛顿发明微积分以来,微分方程在描述事物运动规律上已发挥了重要的作用。实际应用问题通过数学建模所得到的方程,绝大多数是微分方程。由于实际应用的需要,人们必须求解微分方程。然而能够求得解析解的微分方程十分有限,绝大多数微分方程需要利用数值方法来近似求解。下面介绍如何用 Matlab 来计算微分方程(组)的数值解。

    Euler折线法

    微分方程的基本数值解法——Euler折线法。

    步骤:
    分割求解区间,差商代替微商,解代数方程

    例子:
    在这里插入图片描述解:

    等距剖分:a=x0<x1<x2<<xn1<xn=ba=x_0<x_1<x_2< \cdots <x_{n-1}<x_n=b

    步长:h=xk+1xk=(ba)/n,k=0,1,2,...,n1h=x_{k+1}-x_k=(b-a)/n,k=0,1,2,...,n-1

    y(xk+1)y(xk)+hy(xk) y(x_{k+1})\approx y(x_k)+hy\prime(x_k)

    取步长 h = (2 - 0)/n = 2/n,得差分方程
    在这里插入图片描述
    当 h=0.4,即 n=5 时,代码如下:

    clear
    a=0;  b=2;
    h=0.4;
    n=(b-a)/h+1; 
    x(1)=0;  y(1)=1;
    for i=1:n-1        
        y(i+1)=y(i)+h*(y(i)+2*x(i)/y(i)^2);
        x(i+1)=x(i)+h;
    end
    y1=((5/3)*exp(3*x)-2*x-2/3).^(1/3);
    plot(x,y1,x,y);
    h = legend('解析解','Euler解');
    

    在这里插入图片描述

    解得:
    y=(53e3x2x23)1/3 y=(\frac{5}{3}e^{3x}-2x-\frac{2}{3})^{1/3}

    Runge-Kutta方法

    龙格-库塔法(Runge-Kutta methods)是用于非线性常微分方程的解的重要的一类隐式或显式迭代法。

    同样是上面的例子:
    在这里插入图片描述
    在这里插入图片描述

    clear;
    a=0; b=2; h=0.4;
    n=(b-a)/h+1; 
    x(1)=0; y(1)=1;  
    for i=1:n-1
        L1=y(i)+2*x(i)/y(i)^2;
        t=x(i)+h/2;t1=y(i)+h*L1/2;
        L2=t1+2*t/t1^2;t1=y(i)+h*L2/2;
        L3=t1+2*t/t1^2;t1=y(i)+h*L3/2;
        L4=t1+2*(x(i)+h)/t1^2;
        y(i+1)=y(i)+h*(L1+2*L2+2*L3+L4)/6; 
        x(i+1)=x(i)+h;
    end
    plot(x,y)
    

    在这里插入图片描述

    MATLAB自带的ODE求解器

    ODE (常微分方程)

    没有一种算法可以有效地解决所有的 ODE 问题,因此MATLAB 提供了多种ODE求解器,对于不同的ODE,可以调用不同的求解器。
    在这里插入图片描述

    基本语法
    [T,Y] = solver(odefun,tspan,y0)

    其中 y0 为初值条件,tspan为求解区间;Matlab在数值求解时自动对求解区间进行分割,T (向量) 中返回的是分割点的值(自变量),Y (向量) 中返回的是解函数在这些分割点上的函数值。
    solver 为Matlab的ODE求解器(可以是 ode45、ode23、ode113、ode15s、ode23s、ode23t、ode23tb)

    展开全文
  • 文章目录前言1 常微分方程1.1 常微分方程的概念1.2 常微分方程数值求解的一般概念2 常微分方程数值求解函数3 刚性问题结语 前言 1 常微分方程 1.1 常微分方程的概念 1.2 常微分方程数值求解的一般概念 求解常微分...

    前言

    今天我们要说的就是利用MATLAB对常微分方程进行数值求解。本文是科学计算与MATLAB语言专题六第5小节的学习笔记,如果大家有时间的话,可以去听听课,没有的话,可以看看我的笔记,还是很不错的。

    1 常微分方程

    1.1 常微分方程的概念

    在初等数学中就有各种各样的方程,比如线性方程、二次方程、高次方程、指数方程、对数方程、三角方程和方程组等等。这些方程都是要把研究的问题中的已知数和未知数之间的关系找出来,列出包含一个未知数或几个未知数的一个或者多个方程式,然后取求方程的解。
    但是在实际工作中,常常出现一些特点和以上方程完全不同的问题。比如:物质在一定条件下的运动变化,要寻求它的运动、变化的规律;某个物体在重力作用下自由下落,要寻求下落距离随时间变化的规律;火箭在发动机推动下在空间飞行,要寻求它飞行的轨道,等等,要以现有数据求得出形式上的函数解析式,而不是以已知函数来计算特定的未知数。
    物质运动和它的变化规律在数学上是用函数关系来描述的,因此,这类问题就是要去寻求满足某些条件的一个或者几个未知函数。也就是说,凡是这类问题都不是简单地去求一个或者几个固定不变的数值,而是要求一个或者几个未知的函数。
    解这类问题的基本思想和初等数学解方程的基本思想很相似,也是要把研究的问题中已知函数和未知函数之间的关系找出来,从列出的包含未知函数的一个或几个方程中去求得未知函数的表达式。但是无论在方程的形式、求解的具体方法、求出解的性质等方面,都和初等数学中的解方程有许多不同的地方。
    在数学上,解这类方程,要用到微分和导数的知识。因此,凡是表示未知函数的导数以及自变量之间的关系的方程,就叫做微分方程。1

    1.2 常微分方程的定义

    定义1:凡含有参数,未知函数和未知函数导数 (或微分) 的方程,称为微分方程,有时简称为方程,未知函数是一元函数的微分方程称作常微分方程,未知函数是多元函数的微分方程称作偏微分方程。微分方程中出现的未知函数最高阶导数的阶数,称为微分方程的阶。定义式如下:
    F(x,y,y,y,,y(n))=0F(x,y,y',y'',…,y^{(n)})=0  
    定义2:任何代入微分方程后使其成为恒等式的函数,都叫做该方程的解
    若微分方程的解中含有任意常数的个数与方程的阶数相同,且任意常数之间不能合并,则称此解为该方程的通解(或一般解).
    当通解中的各任意常数都取特定值时所得到的解,称为方程的特解
    一般地说,nn阶微分方程的解含有 nn个任意常数。
    也就是说,微分方程的解中含有任意常数的个数和方程的阶数相同,这种解叫做微分方程的通解。通解构成一个函数族。
    如果根据实际问题要求出其中满足某种指定条件的解来,那么求这种解的问题叫做定解问题,对于一个常微分方程的满足定解条件的解叫做特解。
    对于高阶微分方程可以引入新的未知函数,把它化为多个一阶微分方程组。2

    1.3 常微分方程数值求解的一般概念

    求解常微分方程初值问题就是寻找函数y(t)y(t)使之满足如下方程:
    y=f(ty),t0tb,y(t0)=y0y'=f(t,y),t_0≤t≤b, y(t_0)=y_0
    所谓其数值解法,就是求y(t)y(t)在离散结点t处的函数近似值yy。的方法,ynyxny_n\approx y_{xn}
    这些近似值称为常微分方程初值问题的数值解。相邻两个结点之间的距离称为步长。
    单步法:在计算yn+1y_{n+1}时只用到前一步的yny_n,因此在有了初值之后就可以逐步往下计算,其代表是龙格-库塔(Runge-Kutta)法。$$
    多步法:在计算yn+1y_{n+1}时,除了用到前一步的值yny_n之外,还要用到ynpp=1,2,,k,k>0y_{n-p}(p=1,2,…,k,k>0)的值,即前面的k步。其代表就是亚当斯Adams(Adams)法。

    2 常微分方程数值求解函数

    MATLAB提供了多个求常微分方程初值问题数值解的函数,一般调用格式为:
    [t,y]=solver(filename,tspan,y0,option)
    其中,t和y分别给出时间向量和相应的数值解。
    solver为求常微分方程数值解的函数。
    filename是定义f(t,y)的函数名,该函数必须返回一个列向量。
    tspan形式为[to,tf],表示求解区间。
    y0是初始状态向量。
    option是可选参数,用于设置求解属性,常用的属性包括相对误差值Re1To1(默认值是10310^{-3})和绝对误差值AbsTol(默认值是10610^{-6})。
    常微分方程数值求解函数的统一命名格式:odennxx
    ode是0rdinary Differential Equation的缩写,是常微分方程的意思。
    nn是数字,代表所用方法的阶数。例如,ode23采用2阶龙格-库塔RungeKutta(Runge-Kutta)算法,用3阶公式做误差估计来调节步长,具有低等精度。
    ode45采用4阶龙格-库塔算法,用5阶公式做误差估计来调节步长,具有中等精度。
    xx是字母,用于标注方法的专门特征。
    ode15s、ode23s中的“s”代表(Stiff),表示函数适用于刚性方程。
    求常微分方程数值解的函数

    求解函数 采用方法 适用场合
    ode23 2阶或3阶Runge-Kutta算法,低精度 非刚性
    ode45 4阶或5阶Runge-Kutta算法,中精度 非刚性
    ode113 Adams算法,精度可到10-3~10-6 非刚性,计算时间比ode45短
    ode15s Gear’s反向数值微分算法,中精度 刚性
    ode23s 2阶Rosebrock算法,低精度 刚性,当精度较低时,计算时间比ode15s短
    ode23tb 梯形算法,低精度 刚性,当精度较低时,计算时间比ode15s短

    例1求解微分方程初值问题,并与精确解yt=t+1+1y(t)=\sqrt {t+1}+1进行比较。
    {y=y2t24(t+1),0t10y(0)=2 \left\{ \begin{aligned} y'&=\frac{y^2-t-2}{4(t+1)},0≤t≤10\\ y(0)&=2 \end{aligned} \right.

    f=@(t,y)(y^2-t-2)/4/(t+1);
    [t,y]=ode23(f,[0,10],2);
    y1=sqrt(t+1)+1;
    plot(t,y,'b:',t,y1,'r');
    

    6.51
    例2已知一个二阶线性系统的微分方程为:
    {d2xdt2+ax=0,a>0x(0)=0,x(0)=1 \left\{ \begin{aligned} \frac{d^2x}{dt^2} +ax&=0,a>0\\ x(0)&=0 ,x'(0)=1 \end{aligned} \right.
    取a=2,绘制系统的时间响应曲线和相平面图。令x2=x,x1=x’,则得到系统的状态方程:
    {x2=x1x1=ax2x2(0)=0,x1(0)=1 \left\{ \begin{aligned} x_2'&=x_1\\ x_1'&=-ax_2\\ x_2(0)&=0,x_1(0)=1\\ \end{aligned} \right.

    3 刚性问题

    有一类常微分方程,其解的分量有的变化很快,有的变化很慢,且相差悬殊,这就是所谓的刚性问题(Stiff)。
    对于刚性问题,数值解算法必须取很小步长才能获得满意的结果,导致计算量会大大增加。解决刚性问题需要有专门方法。非刚性算法可以求解刚性问题,只不过需要很长的计算时间。

    f=@(t,x)[-2,0;0,1]*[x(2);x(1)];
    [t,x]=ode45(f,[0,20],[1,0]); 
    subplot(1,2,1);
    plot(t,x(:,2)); 
    subplot(1,2,2);
    plot(x(:,2),x(:,1));
    

    6.52
    例3假如点燃一个火柴,火焰球迅速增大直至某个临界体积,然后,维持这一体积不变,原因是火焰球内部燃烧耗费的氧气和从球表面所获氧气达到平衡。其简化模型如下:
    {y=y2y3y(0)=λ,0t2/λ\left\{ \begin{aligned} y'&=y2-y3\\ y(0)&=\lambda ,0≤t≤2/\lambda \end{aligned} \right.
    其中,yty(t)代表火焰球半径,初始半径是λ\lambda,它很小。分析入的大小对方程求解过程的影响。

    lamda=0.01;
    f=@(t,y) y^2-y^3;
    tic;[t,y]=ode45(f,[0,2/lamda],lamda);toc 
    disp(['ode45计算的点数' num2str(length(t))]);
    
    时间已过0.005774秒。
    ode45计算的点数157
    

    tic和toc函数用来记录微分方程求解命令执行的时间,使用tic函数启动计时器,使用toc函数显示从计时器启动到当前所经历的时间。最后还输出计算的点数,运行结果表明这时常微分方程不算很刚性。
    注意:
    1 这里的代码易出现错误。例如1和l。
    2 λ\lambda的拼写应为lambda。

    lamda=1e-5;
    f=@(t,y)y^2-y^3;
    tic;[t,y]=ode45(f,[0,2/lamda],lamda);toc 
    disp(['ode45计算的点数' num2str(length(t))]);
    
    时间已过 1.729913 秒。
    ode45计算的点数120565
    

    这时计算时间明显加长,计算的点数剧增,表明这时常微分方程表现为刚性。

    lamda=1e-5;
    f=@(t,y)y^2-y^3;
    tic;[t,y]=ode15s(f,[0,2/lamda],lamda);toc 
    disp(['ode15s计算的点数' num2str(length(t))]);
    

    时间已过 0.321166 秒。
    ode15s计算的点数98
    对于刚性问题,我们需要改变求解算法。我们选择以“s”结尾的函数,他们专门用于求解刚性方程。计算时间大大缩短,计算的点数大大减少。

    结语

    欢迎大家点赞👍,收藏⭐,转发🚀
    如有问题、建议,请您在评论区留言💬哦。


    1. 百度百科-常微分方程。 ↩︎

    2. 百度百科-常微分方程。 ↩︎

    展开全文
  • ![图片说明](https://img-ask.csdn.net/upload/202010/30/1604047742_249951.jpg) 二阶微分方程如何求解,对x,不胜感激
  • 萌新希望通过fsolve求解一个方程组在不同参数条件下的解,但是没有找到很好的办法来控制这些参数。 如图,问题的具体情境是这样,根据空气动力学知识,我列出了一个方程组,需要求解...
  • 该程序利用矢量化松弛算法来求解线性偏微分方程。 请提交您对改进程序的任何建议(我最终希望包含一个收敛常数以实现 SOR 方法)。 例子: %这是一个如何将程序用于网格空间中具有随机点电荷的平行板电容器的示例: ...
  • 如何利用MATLAB建立Lotka-Volterra模型及其改进模型

    千次阅读 多人点赞 2020-08-14 19:13:25
    今天我们要说的就是利用MATLAB对常微分方程进行数值求解的实例应用-Lotka-Volterra模型及其改进模型。20世纪40年代,Lotka(1925)和Volterra(1926)奠定了种间竞争关系的理论基础,他们提出的种间竞争方程对现代...
  •   本文主要介绍如何利用MATLAB解特征方程,并将特征根的分布画在坐标轴上,便于分析系统的稳定性    我们知道,一旦求出系统的闭环特征根就很容易判定系统的稳定性,但是对于高阶系统,闭环特征根求起来是很...
  • 本教程展示了如何使用 MATLAB 求解器 DDE23 求解具有恒定延迟的延迟微分方程 (DDE)。 求解器在 MATLAB 6.5 及更高版本中可用。 本教程简要讨论了解决ODE之间的区别和 DDE,并描述了 DDE23 中使用的技术。 求解器的...
  • 一、前言 做CFD的气动声学仿真的时候,尤其是出于气动声学本身对高精度格式的需要,我们通常要将仿真的...这篇博文记录我如何利用matlab计算声波(acoustic wave)在某一时间下的压力脉动。 二、LEE(线性化欧拉方
  • Matlab方程式求根

    2021-03-24 08:16:27
    Matlab方程式求解定义符号与求解一元方程与函数1、定义符号`syms/sym()`:2、求解方程`solve()`: f(x)=xsinx-x=03、符号表示的方程: y=ax^2^ +b=04、利用符号直接算微分`diff()`:5、求原函数(积分)`int()`:...
  • 采用ODE23或者另外的MATLAB ODE求解器求解方程系统,定义起始和停止时间以及初识的状态向量。例如: t0 = 5 ; % 起始时间 tf = 20 ; % 停止时间 x0 = [1 –1 3] ; % 初识条件 [t , s] = ode23 ; x = ...
  • ode45求解微分方程

    2021-03-27 19:15:23
    如何利用matlab中的ode45来求解一阶和二阶微分方程 建议大家先看ode45的使用文档,再来看我的博客,拿我这个当个练习做吧,文档链接在文章底部 1.求解一阶微分方程 这是随机共振中的Langevin方程,现在我们来求解x ...
  • 实验三 求代数方程的近似根(解) 数学实验 问题背景和实验目的 实验三近似求解代数方程 解方程代数方程是最常见的...本实验主要介绍一些有效的求解方程的数值方法对分法迭代法 和 牛顿法同时要求大家学会如何利用Matlab
  • 求微分方程的解 自牛顿发明微积分以来微分方程在描述事物运动规律上已发挥了...本实验主要探讨如何Matlab 来计算微分方程组的数值解并重点介绍一个求解微分方程的基本数值解法Euler折线法 问题背景和实验目的 考虑
  • 如何利用电子计算机来快速、有效地求解线性方程组是数值线性代数研究的核心问题,而且也是目前人在继续研究的重大课题之一。 Cholesky分解法又叫做平方根法,是求解对称正定线性方程组最常用的方法之一。A是一个对称...
  • 如何利用电子计算机来快速、有效地求解线性方程组是数值线性代数研究的核心问题,而且也是目前人在继续研究的重大课题之一。 考虑非奇异的线性代数方程组Ax=b,令A=D-L-U,其中D=diag(diag(A)),L=-tril(A,-1),U=-...
  • 如何利用电子计算机来快速、有效地求解线性方程组是数值线性代数研究的核心问题,而且也是目前人在继续研究的重大课题之一。 考虑非奇异的线性代数方程组Ax=b,令A=D-L-U,其中D=diag(diag(A)),L=-tril(A,-1),U=-...
  • 8.8 符号方程求解 符号运算不仅能够求解方程还可以求解方程组下面介绍如何利用函数solve( )求解符号代数方程组函数dsolve( )求解微分方程组 8.8.1 符号代数方程组的求解 在MATLAB中利用函数solve( )求解一般符号代数...
  • 第11章至第15章为应用篇,介绍如何利用MATLAB求解实际的数学建模问题,给出了蚂蚁算法、模拟退火算法、神经元网络、图论算法和遗传算法等详细的算法原理、问题描述、数学模型建立与求解、模型验证和仿真代码的全部...
  • 此提交展示了如何利用 MATLAB 环境中存在的报告工具来比较所有这些方法的解决方案,并展示每种方法的优势: 1. 使用符号数学工具箱的符号完整解决方案。 2. 使用蒙特卡罗仿真对 MATLAB 函数进行参数扫描。 3. 为不同...
  • matlab中的矩阵

    千次阅读 2009-12-13 22:27:00
    我们知道,求解线性方程组是线性代数课程中的核心内容,而矩阵又在求解线性方程组的过程中扮演着举足轻重的角色。下面我们就利用科学计算软件MATLAB来演示如何使用矩
  • Matlab中的矩阵

    千次阅读 2012-12-01 22:45:06
    下面我们就利用科学计算软件MATLAB来演示如何使用矩阵,同时,也使学生对线性代数的认识更加理性。 一、矩阵的构造 在MatLab中,构造矩阵的方法有两种。一种是直接法,就是通过键盘输入的方式直接构造矩阵。另一种是...
  • MATLAB 数值计算

    2010-06-04 21:32:12
    但与一般数值计算教科书不同,本章的讨论重点是:如何利用现有的世界顶级数值计算资源MATLAB。至于数学描述,本章将遵循“最低限度自封闭”的原则处理,以最简明的方式阐述理论数学、数值数学和MATLAB计算指令之间的...
  • 新手上路,matlab基础(3)-matlab中的矩阵

    千次阅读 2007-08-19 09:42:00
    matlab中的矩阵我们知道,求解线性方程组是线性代数课程中的核心内容,而矩阵又在求解线性方程组的过程中扮演着举足轻重的角色。下面我们就利用科学计算软件MATLAB来演示如何使用矩阵,同时,也使学生对线性代数的认识...
  • 我学习了LMI理论,Matlab中LMI工具箱的使用,研究了常见的控制问题与LMI关系以及其表达式,并研究了基于LMI方法的鲁棒控制器设计问题,推导了如何将鲁棒控制器设计问题转化为LMI形式,给出了通过求解LMI方程构造控制...

空空如也

空空如也

1 2 3
收藏数 48
精华内容 19
关键字:

如何利用matlab求解方程

matlab 订阅