? ? ? ?此教程說明如何使用 MATLAB 構造幾種不同類型的微分方程并求解。MATLAB 提供了多種數值算法來求解各種微分方程：

• 初始值問題

• 邊界值問題

• 時滯微分方程

• 偏微分方程

初始值問題

vanderpoldemo?是用于定義 van der Pol 方程的函數

type vanderpoldemofunction dydt = vanderpoldemo(t,y,Mu)
dydt = [y(2); Mu*(1-y(1)^2)*y(2)-y(1)];

? ? ? ?該方程寫作包含兩個一階常微分方程 (ODE) 的方程組。將針對參數?μ?的不同值計算這些方程。為了實現更快的積分，您應該根據?μ?的值選擇合適的求解器。

? ? ? ?當?μ=1?時，任何 MATLAB ODE 求解器都能有效地求解 van der Pol 方程。ode45?求解器就是其中之一。該方程在域?[0,20]?中求解，初始條件為?y(0)=2?和?dydt?t=0=0。

tspan = [0 20];
y0 = [2; 0];
Mu = 1;
ode = @(t,y) vanderpoldemo(t,y,Mu);
[t,y] = ode45(ode, tspan, y0);% Plot solution
plot(t,y(:,1))
xlabel('t')
ylabel('solution y')
title('van der Pol Equation, \mu = 1')

? ? ? ?對于較大的?μ，問題將變為剛性。此標簽表示拒絕使用普通方法計算的問題。這種情況下，要實現快速積分，需要使用特殊的數值方法。ode15s、ode23s、ode23t?和?ode23tb?函數可有效地求解剛性問題。

? ? ? 當?μ=1000?時，van der Pol 方程的求解使用?ode15s，初始條件相同。您需要將時間范圍大幅度延長到?[0,3000]?才能看到解的周期性變化。

tspan = [0, 3000];
y0 = [2; 0];
Mu = 1000;
ode = @(t,y) vanderpoldemo(t,y,Mu);
[t,y] = ode15s(ode, tspan, y0);
plot(t,y(:,1))
title('van der Pol Equation, \mu = 1000')
axis([0 3000 -3 3])
xlabel('t')
ylabel('solution y')

邊界值問題

bvp4c?和?bvp5c?可以求解常微分方程的邊界值問題。

? ? ? 示例函數?twoode?將一個微分方程寫作包含兩個一階 ODE 的方程組。此微分方程為

type twoodefunction dydx = twoode(x,y)
%TWOODE Evaluate the differential equations for TWOBVP.
dydx = [ y(2); -abs(y(1)) ];

函數?twobc?求解該問題的邊界條件為：y(0)=0?和?y(4)=?2。

type twobcfunction res = twobc(ya,yb)
res = [ ya(1); yb(1) + 2 ];

? ? ? ?在調用?bvp4c?之前，您必須為要在網格中表示的解提供一個猜想值。然后，求解器就像對解進行平滑處理一樣修改網格。

bvpinit?函數以您可以傳遞給求解器?bvp4c?的形式設定初始猜想值。對于?[0 1 2 3 4]?的網格以及?y(x)=1?和?y'(x)=0?的常量猜想值，對?bvpinit?的調用為：

solinit = bvpinit([0 1 2 3 4],[1; 0]);

? ? ? ?利用這個初始猜想值，您可以使用?bvp4c?對該問題求解。使用?deval?計算?bvp4c?在某些點返回的解，然后繪制結果值。

sol = bvp4c(@twoode, @twobc, solinit);
xint = linspace(0, 4, 50);
yint = deval(sol, xint);
plot(xint, yint(1,:));
xlabel('x')
ylabel('solution y')
hold on

? ? ? 此特定的邊界值問題實際上有兩種解。通過將初始猜想值更改為?y(x)=?1?和?y'(x)=0，可以求出另一個解。

solinit = bvpinit([0 1 2 3 4],[-1; 0]);
sol = bvp4c(@twoode,@twobc,solinit);
xint = linspace(0,4,50);
yint = deval(sol,xint);
plot(xint,yint(1,:));
legend('Solution 1','Solution 2')
hold off

時滯微分方程

? ? ?dde23、ddesd?和?ddensd?可以求解具有各種時滯的時滯微分方程。示例?ddex1、ddex2、ddex3、ddex4?和?ddex5?構成了這些求解器的迷你使用教程。

ddex1?示例說明如何求解微分方程組

y′1(t)=y1(t?1)y′2(t)=y1(t?1)+y2(t?0.2)y′3(t)=y2(t).

您可以使用匿名函數表示這些方程

ddex1fun = @(t,y,Z) [Z(1,1); Z(1,1)+Z(2,2); y(2)];

問題的歷史解(t≤0?時)固定不變：

y1(t)=1y2(t)=1y3(t)=1.

您可以將歷史解表示為由 1 組成的向量。

ddex1hist = ones(3,1);

采用二元素向量表示方程組中的時滯。

lags = [1 0.2];

? ? ? 將函數、時滯、歷史解和積分區間?[0,5]?作為輸入傳遞給求解器。求解器在整個積分區間生成適合繪圖的連續解。

sol = dde23(ddex1fun, lags, ddex1hist, [0 5]);
plot(sol.x,sol.y);
title({'An example of Wille and Baker', 'DDE with Constant Delays'});
xlabel('time t');
ylabel('solution y');
legend('y_1','y_2','y_3','Location','NorthWest');

偏微分方程

? ? ?pdepe?使用一個空間變量和時間對偏微分方程求解。示例?pdex1、pdex2、pdex3、pdex4?和?pdex5?構成了?pdepe?的迷你使用教程。

此示例問題使用函數?pdex1pde、pdex1ic?和?pdex1bc。

pdex1pde?定義微分方程

type pdex1pdefunction [c,f,s] = pdex1pde(x,t,u,DuDx)
%PDEX1PDE Evaluate the differential equations components for the PDEX1 problem.
c = pi^2;
f = DuDx;
s = 0;

pdex1ic?設置初始條件

type pdex1icfunction u0 = pdex1ic(x)
%PDEX1IC Evaluate the initial conditions for the problem coded in PDEX1.
u0 = sin(pi*x);

pdex1bc?設置邊界條件

u(0,t)=0,

πe?t+??xu(1,t)=0.

type pdex1bcfunction [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
%PDEX1BC Evaluate the boundary conditions for the problem coded in PDEX1.
pl = ul;
ql = 0;
pr = pi * exp(-t);
qr = 1;

pdepe?需要提供空間離散?x?和時間向量?t(您要獲取解快照的時間點)。使用包含 20 個節點的網格求解此問題，并請求五個?t?值的解。提取解的第一個分量并繪圖。

x = linspace(0,1,20);
t = [0 0.5 1 1.5 2];
sol = pdepe(0,@pdex1pde,@pdex1ic,@pdex1bc,x,t);
u1 = sol(:,:,1);
surf(x,t,u1);
xlabel('x');
ylabel('t');
zlabel('u');

關注公眾號：MATLAB基于模型的設計 (ID：xaxymaker) ，每天推送MATLAB學習最常見的問題，每天進步一點點，業精于勤荒于嬉。

可保存后掃碼關注哦！

與50位技術專家面對面20年技術見證，附贈技術全景圖

總結

以上是生活随笔為你收集整理的matlab 解方程组_一文读懂MATLAB微分方程的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。