solving simultaneous equations of fourth order PDE in MATLAB

11 ビュー (過去 30 日間)
Devansh Gupta 2021 年 3 月 13 日

I'm trying to solve this system of ordinary differential equation with the boundary conditions : w(0)=w(L)=u(0)=u(L)=w'(0)=w'(L)=u'(0)=u'(L)=0.
I'm new to MATLAB and am unable to solve it using ODE toolbox.
Can someone with knowledge on the topic help me with this?
10 件のコメント表示非表示 9 件の古いコメント
Devansh Gupta 2021 年 3 月 15 日
I've expressed it like this
cons1 = ((A(1,1)*D(1,1))-(B(1,1)^2))/(A(1,1));
cons2 = (B(1,1)/A(1,1));
syms w(x)
Dq = diff((q*x),x);
Dw = diff(w,x);
D2w = diff(w,x,2);
D3w = diff(w,x,3);
D4w = diff(w,x,4);
ode = cons1 * D4w == -cons2 * Dq;
but I'm unable to apply the boundary conditions to this expression.

サインインしてコメントする。

回答 (1 件)

David Goodmanson 2021 年 3 月 17 日

Hello Devansh
Since you have boundary conditions at both x=0 and x=L you can use bvp4c to solve this. But I believe you have too many boundary conditions. Dropping the subscrips on A,B,D and denoting px by p and pz by q, then
Au'' - Bw''' + p(x) = 0 (1)
Bu''' - Dw'''' + q(x) = 0 (2)
and substituting as before gives
w'''' = (Aq(x) - Bp'(x))/(AD-B^2) (3)
u'' = (Bw'''-p(x))/A (4)
( u''' is not involved any more since the substitution process assures that if (3) and (4) work, so does (2) ).
There is a fourth order and a second order equation, so there will be six bc's and not eight. If you keep
w(0) = w'(0) = w(L) = w'(L) = 0,
then there are two bc's left for u. The code below assumes
u(0) = u(L) = 0
but see the note at the end. I assumed p, p' and q to be known algebraic functions.
L = 1.3;
% p,q,A,B,D are in function below
xvec = linspace(0,L,100);
solinit = bvpinit(xvec, @guess);
sol = bvp4c(@fun, @bcs, solinit);
x = sol.x';
% wu ~~ [w w' w'' w''' u u']
wu = sol.y';
w = wu(:,1:4);
u = wu(:,5:6);
figure(1); grid on
plot(x,wu)
legend('w','w''','w''''','w''''''','u','u''')
function dbydx = fun(x,wu)
% wu ~~ [w w' w'' w''' u u']
A = 1;
B = 2;
D = 1;
q = cos(7*x);
p = x^(3/2);
pprime = (3/2)*x^(1/2);
wiv = (A*q-B*pprime)/(A*D-B^2);
uii = (B*wu(4) -p)/A;
dbydx = [wu(2:4); wiv; wu(6); uii];
end
function bcvec = bcs(ya,yb)
bcvec = [ya([1 2 5]); yb([1 2 5])];
end
function g = guess(x)
g = [0 0 0 0 0 0]';
end
For bc's on u, you can specify values of u at both ends, or a value of u at one end and a value of u' at the other, but you can't specify values of u' at both ends. That's because integrating both sides of (4) from 0 to L gives
u'(L) - u'(0) = integral(right hand side)
and that condition can't be met with arbitrary choices of u' at each end.

サインインしてコメントする。

R2019a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by