matlab 解方程组_一文读懂MATLAB微分方程
? ? ? ?此教程說明如何使用 MATLAB 構(gòu)造幾種不同類型的微分方程并求解。MATLAB 提供了多種數(shù)值算法來求解各種微分方程：
初始值問題
邊界值問題
時滯微分方程
偏微分方程
初始值問題
vanderpoldemo?是用于定義 van der Pol 方程的函數(shù)
type vanderpoldemofunction dydt = vanderpoldemo(t,y,Mu)dydt = [y(2); Mu*(1-y(1)^2)*y(2)-y(1)];
? ? ? ?該方程寫作包含兩個一階常微分方程 (ODE) 的方程組。將針對參數(shù)?μ?的不同值計算這些方程。為了實現(xiàn)更快的積分，您應該根據(jù)?μ?的值選擇合適的求解器。
? ? ? ?當?μ=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')
? ? ? ?對于較大的?μ，問題將變?yōu)閯傂浴４藰撕灡硎揪芙^使用普通方法計算的問題。這種情況下，要實現(xiàn)快速積分，需要使用特殊的數(shù)值方法。ode15s、ode23s、ode23t?和?ode23tb?函數(shù)可有效地求解剛性問題。
? ? ? 當?μ=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?可以求解常微分方程的邊界值問題。
? ? ? 示例函數(shù)?twoode?將一個微分方程寫作包含兩個一階 ODE 的方程組。此微分方程為
type twoodefunction dydx = twoode(x,y)%TWOODE Evaluate the differential equations for TWOBVP.
dydx = [ y(2); -abs(y(1)) ];
函數(shù)?twobc?求解該問題的邊界條件為：y(0)=0?和?y(4)=?2。
type twobcfunction res = twobc(ya,yb)res = [ ya(1); yb(1) + 2 ];
? ? ? ?在調(diào)用?bvp4c?之前，您必須為要在網(wǎng)格中表示的解提供一個猜想值。然后，求解器就像對解進行平滑處理一樣修改網(wǎng)格。
bvpinit?函數(shù)以您可以傳遞給求解器?bvp4c?的形式設(shè)定初始猜想值。對于?[0 1 2 3 4]?的網(wǎng)格以及?y(x)=1?和?y'(x)=0?的常量猜想值，對?bvpinit?的調(diào)用為：
solinit = bvpinit([0 1 2 3 4],[1; 0]);? ? ? ?利用這個初始猜想值，您可以使用?bvp4c?對該問題求解。使用?deval?計算?bvp4c?在某些點返回的解，然后繪制結(jié)果值。
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?構(gòu)成了這些求解器的迷你使用教程。
ddex1?示例說明如何求解微分方程組
y′1(t)=y1(t?1)y′2(t)=y1(t?1)+y2(t?0.2)y′3(t)=y2(t).
您可以使用匿名函數(shù)表示這些方程
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];? ? ? 將函數(shù)、時滯、歷史解和積分區(qū)間?[0,5]?作為輸入傳遞給求解器。求解器在整個積分區(qū)間生成適合繪圖的連續(xù)解。
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?構(gòu)成了?pdepe?的迷你使用教程。
此示例問題使用函數(shù)?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?設(shè)置初始條件
type pdex1icfunction u0 = pdex1ic(x)%PDEX1IC Evaluate the initial conditions for the problem coded in PDEX1.
u0 = sin(pi*x);
pdex1bc?設(shè)置邊界條件
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 個節(jié)點的網(wǎng)格求解此問題，并請求五個?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');
關(guān)注公眾號：MATLAB基于模型的設(shè)計 (ID：xaxymaker) ，每天推送MATLAB學習最常見的問題，每天進步一點點，業(yè)精于勤荒于嬉。
可保存后掃碼關(guān)注哦！
與50位技術(shù)專家面對面20年技術(shù)見證，附贈技術(shù)全景圖總結(jié)
以上是生活随笔為你收集整理的matlab 解方程组_一文读懂MATLAB微分方程的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： wavecn 2.0.0.5 正式版_微
- 下一篇： python通过url下载文件不可读_p