solving simultaneous equations of fourth order PDE in MATLAB

5 次查看（过去 30 天）

Devansh Gupta 2021-3-13

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/771978-solving-simultaneous-equations-of-fourth-order-pde-in-matlab

编辑： David Goodmanson 2021-3-18

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 个评论
显示 8更早的评论隐藏 8更早的评论

darova 2021-3-15

Yes, the answer is yes. You need to express

and use ode45

Devansh Gupta 2021-3-15

编辑：Devansh Gupta 2021-3-15

在 MATLAB Online 中打开

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

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/771978-solving-simultaneous-equations-of-fourth-order-pde-in-matlab#answer_649962

编辑：David Goodmanson 2021-3-18

在 MATLAB Online 中打开

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.