In solving ordinary differential equations of elasticity, I got the following hint“Unable to find symbolic solution”.How should I modify the solution

Question

Guoyao.Li 2021-8-23

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1438644-in-solving-ordinary-differential-equations-of-elasticity-i-got-the-following-hint-unable-to-find-sy

评论： Wan Ji 2021-8-23

syms sigma_theta(r) sigma_r(r) p p0 E mu a b
u=r/E*(-mu*sigma_r+(1-mu)*sigma_theta);
eqn=[2*(sigma_theta-sigma_r)/r==diff(sigma_r,r,1),...
    diff(u,r,1)==1/E*(sigma_r-2*mu*sigma_theta)];
cond=[sigma_theta(a)==-p,sigma_theta(b)==-p0];
S=dsolve(eqn,cond);
Warning: Unable to find symbolic solution.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Wan Ji 2021-8-23

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1438644-in-solving-ordinary-differential-equations-of-elasticity-i-got-the-following-hint-unable-to-find-sy#answer_772504

That means matlab can not find symbolic solution for the ode with the given boundary condintions (Given values at two ends).

Take a look at the general ode, which is always put forward with a given value at only one end called initial conditions.

2 个评论
显示无隐藏无

Guoyao.Li 2021-8-23

在 MATLAB Online 中打开

Thanks!

I try to make two equations into one and get the solution as follows.

syms sigma_theta(r) sigma_r(r) p p0 E mu a b

sigma_theta=1/2*r*diff(sigma_r)+sigma_r;

u=r/E*(-mu*sigma_r+(1-mu)*sigma_theta);

eqn=[diff(u,r,1)==1/E*(sigma_r-2*mu*sigma_theta)];

cond=[sigma_r(a)==-p,sigma_r(b)==-p0];

sigma_r=dsolve(eqn,cond);

sigma_theta=1/2*r*diff(sigma_r);

simplify(sigma_theta)

ans =

simplify(sigma_r)

ans =

Wan Ji 2021-8-23

在 MATLAB Online 中打开

Yes, that worked! This is only because there is a relation between sigma_r and sigma_theta.

I have also found the relation:

syms sigma_theta(r) sigma_r(r) p p0 E mu d1 d2 f(r) c1 c2 s1 s2
syms a b real positive 
% d1 = diff(sigma_r,r,1);
% d2 = diff(sigma_theta,r,1);
u=r/E*(-mu*sigma_r+(1-mu)*sigma_theta);
eqn=[2*(sigma_theta-sigma_r)/r==diff(sigma_r,r,1);...
    diff(u,r,1)==1/E*(sigma_r-2*mu*sigma_theta)]
eq1 = subs(eqn, [diff(sigma_r), diff(sigma_theta)],[d1,d2]);
[d1, d2] = solve(eq1,d1,d2)
d1/d2 % 由于你的sigma_r 和 sigma_theta 对r 对r求导后比例是-2
% 所以令f = sigma_r - sigma_theta
% 有 diff(f,r,1)== -3*f/r 求解就行了
q = dsolve(diff(f,r,1)==(-3*f/r))
% p = int(q/r), get p = c1/(3*r^3)
% p = int(-2*q/r), get p = (2*C1)/(3*r^3)
eq = -s1/(3*r^3) + s2; % 指定两个待定系数；这里的eq就是sigma_theta的表达式
[c1, c2] = solve([subs(eq,r,a)+p; subs(eq, r,b)+p0],[s1,s2]) % 使用边界条件求解两个系数
sigma_theta = subs(eq,[s1,s2],[c1,c2]) % 求得的c1,c2代入eq中
sigma_theta = subs(sigma_theta,[a,b,p0,p],[1,2,5,10]); % 赋予a=1; b=2; p0=5；p=10
rr = 1:0.01:2;
sigma_theta_val = eval(subs(sigma_theta,r,rr));
plot(rr, sigma_theta_val)
xlabel('r'); ylabel('\sigma_\theta')