用拉普拉斯变换求解微分方程

打开实时脚本

通过使用符号数学工具箱™中的拉普拉斯变换来解决微分方程。有关拉普拉斯变换的简单例子，请参见拉普拉斯而且ilaplace．

定义:拉普拉斯变换

一个函数的拉普拉斯变换 $f （ t ）$ 是

$F （年代）＝ \int_{0}^{\infty} f （ t ） e^{- ts} dt ．$

概念:使用符号工作流

符号工作流保持自然符号形式的计算，而不是数字形式。这种方法可以帮助您理解解决方案的属性，并使用精确的符号值。只有在需要数值结果或不能以符号形式继续时，才可以用数字代替符号变量。详细信息请参见选择数字或符号算术．通常，这些步骤是:

申报的方程。
解决方程。
替代的价值观。
阴谋的结果。
分析的结果。

工作流程:用拉普拉斯变换求解RLC电路

申报的方程

你可以用拉普拉斯变换来解有初始条件的微分方程。例如，你可以解决电阻电感电容(RLC)电路，比如这个电路。

欧姆电阻: $R_{1} ， R_{2} ， R_{3.}$
电流(安培): $我_{1} ，我_{2} ，我_{3.}$
亨利电感: $l$
电容以法拉为单位: $C$
交流电压源(单位:伏): $E （ t ）$
电容器电量(库仑): $问（ t ）$

应用基尔霍夫电压和电流定律得到下列方程。

$\begin{array}{l} 我_{1} ＝我_{2} + 我_{3.} \\ l \frac{迪_{1}}{dt} + 我_{1} R_{1} + 我_{2} R_{2} ＝ 0 \\ E （ t ） + 我_{2} R_{2} - \frac{问}{C} - 我_{3.} R_{3.} ＝ 0 \end{array}$

代入关系 $我_{3.} ＝ d 问 / d t$ (即电容被充电的速率)代入上述方程，得到RLC电路的微分方程。

$\begin{array}{l} \frac{迪_{1}}{dt} - \frac{R_{2}}{l} \frac{dQ}{dt} ＝ - \frac{R_{1} + R_{2}}{l} 我_{1} \\ \frac{dQ}{dt} ＝ \frac{1}{R_{2} + R_{3.}} （ E （ t ） - \frac{问}{C} ） + \frac{R_{2}}{R_{2} + R_{3.}} 我_{1} \end{array}$

声明变量。因为物理量为正值，所以对变量设置相应的假设。让 $E （ t ）$ 交流电压为1v。

信谊lCI1 (t)问(t)年代R = sym(‘R % d’3 [1]);假设([t L C R] > 0) E(t) = 1*sin(t);%交流电压= 1v

声明微分方程。

dI1 = diff(I1,t);dQ = diff(Q,t);eqn1 = dI1 -(R(2)/L)*dQ == -(R(1)+R(2))/L*I1

eqn1 (t) =
                  
                    $\frac{\partial}{\partial t} 我_{1} （ t ） - \frac{R_{2} \frac{\partial}{\partial t} 问 （ t ）}{l} ＝ - \frac{我_{1} （ t ） (R_{1} + R_{2})}{l}$

eqn2 = dQ = = (1 / (R (2) + (3)) * (q / C)) + R (2) / (R (2) + (3)) * I1

eqn2 (t) =
                  
                    $\frac{\partial}{\partial t} 问 （ t ） ＝ \frac{罪 （ t ） - \frac{问 （ t ）}{C}}{R_{2} + R_{3.}} + \frac{R_{2} 我_{1} （ t ）}{R_{2} + R_{3.}}$

解决方程

计算的拉普拉斯变换eqn1而且eqn2．

eq1lt =拉普拉斯(eqn1,t,s)

eqn1LT =
                  
                    $年代 拉普拉斯 （ 我_{1} （ t ） ， t ， 年代 ） - 我_{1} （ 0 ） + \frac{R_{2} (问 （ 0 ） - 年代 拉普拉斯 （ 问 （ t ） ， t ， 年代 ）)}{l} ＝ - \frac{(R_{1} + R_{2}) 拉普拉斯 （ 我_{1} （ t ） ， t ， 年代 ）}{l}$

eqn2LT =拉普拉斯(eqn2,t,s)

eqn2LT =
                  
                    $年代 拉普拉斯 （ 问 （ t ） ， t ， 年代 ） - 问 （ 0 ） ＝ \frac{R_{2} 拉普拉斯 （ 我_{1} （ t ） ， t ， 年代 ）}{R_{2} + R_{3.}} + \frac{\frac{C}{{年代}^{2} + 1} - 拉普拉斯 （ 问 （ t ） ， t ， 年代 ）}{C (R_{2} + R_{3.})}$

这个函数解决只求解符号变量。因此，使用解决，第一个代换拉普拉斯(I1 (t), t, s)而且拉普拉斯(Q (t), t, s)有了变量I1_LT而且Q_LT．

信谊I1_LTQ_LTeqn1LT = subs(eqn1LT，[拉普拉斯(I1,t,s)拉普拉斯(Q,t,s)]，[I1_LT Q_LT])

eqn1LT =
                  
                    $我_{1 ， LT} 年代 - 我_{1} （ 0 ） + \frac{R_{2} (问 （ 0 ） - 问_{LT} 年代)}{l} ＝ - \frac{我_{1 ， LT} (R_{1} + R_{2})}{l}$

eqn2LT = subs(eqn2LT，[拉普拉斯(I1,t,s)拉普拉斯(Q,t,s)]，[I1_LT Q_LT])

eqn2LT =
                  
                    $问_{LT} 年代 - 问 （ 0 ） ＝ \frac{我_{1 ， LT} R_{2}}{R_{2} + R_{3.}} - \frac{问_{LT} - \frac{C}{{年代}^{2} + 1}}{C (R_{2} + R_{3.})}$

解的方程I1_LT而且Q_LT．

eqns = [eqn1LT eqn2LT];vars = [I1_LT Q_LT];[I1_LT, Q_LT] = solve(eqns,vars)

I1_LT =
                  
                    $\frac{l 我_{1} （ 0 ） - R_{2} 问 （ 0 ） + C R_{2} 年代 + l {年代}^{2} 我_{1} （ 0 ） - R_{2} {年代}^{2} 问 （ 0 ） + C l R_{2} {年代}^{3.} 我_{1} （ 0 ） + C l R_{3.} {年代}^{3.} 我_{1} （ 0 ） + C l R_{2} 年代 我_{1} （ 0 ） + C l R_{3.} 年代 我_{1} （ 0 ）}{({年代}^{2} + 1) (R_{1} + R_{2} + l 年代 + C l R_{2} {年代}^{2} + C l R_{3.} {年代}^{2} + C R_{1} R_{2} 年代 + C R_{1} R_{3.} 年代 + C R_{2} R_{3.} 年代)}$

Q_LT =
                  
                    $\frac{C (R_{1} + R_{2} + l 年代 + l R_{2} 我_{1} （ 0 ） + R_{1} R_{2} 问 （ 0 ） + R_{1} R_{3.} 问 （ 0 ） + R_{2} R_{3.} 问 （ 0 ） + l R_{2} {年代}^{2} 我_{1} （ 0 ） + l R_{2} {年代}^{3.} 问 （ 0 ） + l R_{3.} {年代}^{3.} 问 （ 0 ） + R_{1} R_{2} {年代}^{2} 问 （ 0 ） + R_{1} R_{3.} {年代}^{2} 问 （ 0 ） + R_{2} R_{3.} {年代}^{2} 问 （ 0 ） + l R_{2} 年代 问 （ 0 ） + l R_{3.} 年代 问 （ 0 ）)}{({年代}^{2} + 1) (R_{1} + R_{2} + l 年代 + C l R_{2} {年代}^{2} + C l R_{3.} {年代}^{2} + C R_{1} R_{2} 年代 + C R_{1} R_{3.} 年代 + C R_{2} R_{3.} 年代)}$

计算 $我_{1}$ 而且 $问$ 通过计算的拉普拉斯逆变换I1_LT而且Q_LT．简化结果。抑制输出，因为它很长。

I1sol = ilaplace(I1_LT,s,t);Qsol = ilaplace(Q_LT,s,t);I1sol = simplify(I1sol);Qsol = simplify(Qsol);

替代值

在绘制结果之前，用电路元件的数值代替符号变量。让 $R_{1} ＝ 4 Ω$ ， $R_{2} ＝ 2 Ω$ ， $R_{3.} ＝ 3. Ω$ ， $C ＝ 1 / 4 F$ ， $l ＝ 1 ． 6 H$ ．假设初始电流为 $我_{1} （ 0 ）＝ 2 一个$ 初始电荷是 $问（ 0 ）＝ 2 C$ ．

vars = [R L C I1(0) Q(0)];数值= [4 2 3 1.6 1/4 2 2];I1sol = subs(I1sol,vars,values)

I1sol =
                  
                    $\frac{200 因为 （ t ）}{8161} + \frac{405 罪 （ t ）}{8161} + \frac{16122 e^{- \frac{81 t}{40}} (\cosh （ \frac{\sqrt{1761} t}{40} ） - \frac{742529 \sqrt{1761} \sinh （ \frac{\sqrt{1761} t}{40} ）}{14195421})}{8161}$

Qsol = subs(Qsol,vars,values)

Qsol =
                  
                    $\frac{924 罪 （ t ）}{8161} - \frac{1055 因为 （ t ）}{8161} + \frac{17377 e^{- \frac{81 t}{40}} (\cosh （ \frac{\sqrt{1761} t}{40} ） + \frac{1109425 \sqrt{1761} \sinh （ \frac{\sqrt{1761} t}{40} ）}{30600897})}{8161}$

阴谋的结果

绘制电流I1sol和电荷Qsol．通过使用两个不同的时间间隔显示瞬态和稳态行为: $0 \leq t \leq 15$ 而且 $2 \leq t \leq 25$ ．

subplot(2,2,1) fplot(I1sol，[0 15])“当前”) ylabel (“I1 (t)”)包含(“t”) subplot(2,2,2) fplot(Qsol，[0 15])“充电”) ylabel (“问(t)”)包含(“t”fplot(I1sol，[2 25])“当前”) ylabel (“I1 (t)”)包含(“t”)文本(3,-0.1,瞬态的-0.07)文本(15日,“稳定状态”fplot(Qsol，[2 25]) title(“充电”) ylabel (“问(t)”)包含(“t”)文本(3,0.35,瞬态的0.22)文本(15日,“稳定状态”）

图中包含4个轴对象。标题为Current的坐标轴对象1包含一个functionline类型的对象。标题为Charge的Axes对象2包含一个functionline类型的对象。标题为Current的坐标轴对象3包含3个类型为functionline, text的对象。标题为Charge的Axes对象4包含3个functionline, text类型的对象。

分析结果

最初，电流和电荷呈指数下降。然而，从长期来看，它们是震荡的。这些行为分别被称为“瞬态”和“稳态”。使用符号结果，您可以分析结果的属性，这在数值结果中是不可能的。

视觉检查I1sol而且Qsol．它们是项的和。通过使用找到这些项孩子们．然后，通过把它们画出来，找出这些项的贡献[0 15]．图中显示了暂态项和稳态项。

I1terms = children(I1sol);I1terms = [I1terms{:}];Qterms = children(Qsol);Qterms = [Qterms{:}];图;subplot(1,2,1) fplot(I1terms，[0 15]) ylim([-0.5 2.5]) title(“目前的条件”) subplot(1,2,2) fplot(Qterms，[0 15]) ylim([-0.5 2.5]) title(“收费条款”）

图中包含2个轴对象。标题为Current terms的Axes对象1包含3个functionline类型的对象。标题为Charge terms的Axes对象2包含3个functionline类型的对象。

图表显示I1sol有瞬态项和稳态项，而Qsol有一个暂态项和两个稳态项。从目测，注意I1sol而且Qsol有一个术语包含经验值函数。假设这一项导致了瞬态指数衰减。把暂态项和稳态项分开I1sol而且Qsol通过检查经验值使用有．

I1transient = I1terms(有(I1terms，“经验”）)

I1transient =
                  
                    $\frac{16122 e^{- \frac{81 t}{40}} (\cosh （ \frac{\sqrt{1761} t}{40} ） - \frac{742529 \sqrt{1761} \sinh （ \frac{\sqrt{1761} t}{40} ）}{14195421})}{8161}$

I1steadystate = I1terms(~有(I1terms，“经验”）)

I1steadystate =
                  
                    $（ \begin{array}{c} \frac{200 因为 （ t ）}{8161} & \frac{405 罪 （ t ）}{8161} \end{array} ）$

同样的,单独的Qsol变成瞬态和稳态。这个结果演示了符号计算如何帮助您分析问题。

Qtransient = Qterms(有(Qterms，“经验”）)

Qtransient =
                  
                    $\frac{17377 e^{- \frac{81 t}{40}} (\cosh （ \frac{\sqrt{1761} t}{40} ） + \frac{1109425 \sqrt{1761} \sinh （ \frac{\sqrt{1761} t}{40} ）}{30600897})}{8161}$

Qsteadystate = Qterms(~有(Qterms，“经验”）)

Qsteadystate =
                  
                    $（ \begin{array}{c} - \frac{1055 因为 （ t ）}{8161} & \frac{924 罪 （ t ）}{8161} \end{array} ）$