Solve differential equation to get variable
显示 更早的评论
% I am trying to proceed to get the solution for a here by equating differentiation of u_a withrespect to a equal to Zero, but got error like!
%Warning: Solutions are only valid under certain conditions. To include parameters and %conditions in the solution, specify the 'ReturnConditions' value as
%'true'.
%> In sym/solve>warnIfParams (line 478)
%In sym/solve (line 357)
%In HW2_Pro2_VirtualForce (line 55)
%>>
clc;
close;
clear;
syms d_F ;
syms d_F_c;
syms dM_B;
%syms Theta_B;
syms x real;
syms u_mid;
syms u_c;
syms u_max;
syms a;
%syms w_x;
syms M_z(x);
syms dM_z(x);
syms dP;
syms u_a;
%syms n;
%Q = 2.5; %kip/ft
L = 30.; % ft
%a = 18; %assume maximum deflection at 15', 20'
E = 29000.*12^2; %ksi to kips/ft2
I = 1890./(12)^4; %ft^4
E_I = E.*I; %kip*ft^2;
%% deflection at mid of AB
w_x = @(x) 1*x/(L/3)*heaviside(x-L/3)*(1-heaviside(x-2*L/3));
w_x_1 = w_x(x);
eq1 = diff(M_z(x),2) == -w_x_1;
sol1(x) = dsolve(eq1, M_z(0)==-150, M_z(L) == 250, 'IgnoreAnalyticConstraints',true);
set(gca, 'XAxisLocation', 'origin', 'YAxisLocation', 'origin')
fplot(x, sol1(x), [0,L]);
figure();
%moment_x = M_z(x)/sol1;
fplot(x, w_x(x), [0,L]);
origin = [0 0];
figure();
%%
%apply virtual load at distance a from support A
%a = L/2;
%dP = -1.;
dM_z= @(x,a) (dP*x*(L-a)/L)-dP.*(x-a).*heaviside(x-a);
%fplot(x, dM_z(x,a), [0,L])
dW = @(x,a) M_z(x)*dM_z(x,a);
dW_int = @(a) 1/E_I.*(int(dW, x, 0 ,L));
dW_ext = dP*u_a;
eqn = dW_int == dW_ext;
sol = solve(eqn, u_a);
u_a_diff = diff(u_a, a);
eqn2 = u_a_diff == 0;
sol_a = solve(eqn2, a);
a_solve = double(sol_a);
采纳的回答
Hi Milan,
I modified as follows. Not sure if this is what you're looking for:
syms d_F ;
syms d_F_c;
syms dM_B;
%syms Theta_B;
syms x real;
syms u_mid;
syms u_c;
syms u_max;
syms a;
%syms w_x;
syms M_z(x);
syms dM_z(x);
syms dP;
syms u_a;
%syms n;
%Q = 2.5; %kip/ft
L = 30.; % ft
%a = 18; %assume maximum deflection at 15', 20'
E = 29000.*12^2; %ksi to kips/ft2
I = 1890./(12)^4; %ft^4
E_I = E.*I; %kip*ft^2;
%% deflection at mid of AB
w_x = @(x) 1*x/(L/3)*heaviside(x-L/3)*(1-heaviside(x-2*L/3));
w_x_1 = w_x(x);
eq1 = diff(M_z(x),2) == -w_x_1;
sol1(x) = dsolve(eq1, M_z(0)==-150, M_z(L) == 250, 'IgnoreAnalyticConstraints',true)
sol1(x) is the solution for M_z(x), so make that assignment because M_z(x) is used below.
M_z(x) = sol1(x);
% set(gca, 'XAxisLocation', 'origin', 'YAxisLocation', 'origin')
% fplot(x, sol1(x), [0,L]);
% figure();
%moment_x = M_z(x)/sol1;
% fplot(x, w_x(x), [0,L]);
% origin = [0 0];
% figure();
%%
%apply virtual load at distance a from support A
%a = L/2;
%dP = -1.;
Use symfun objects in instead of anonymous function
%dM_z= @(x,a) (dP*x*(L-a)/L)-dP.*(x-a).*heaviside(x-a);
dM_z(x,a) = (dP*x*(L-a)/L)-dP.*(x-a).*heaviside(x-a)
%fplot(x, dM_z(x,a), [0,L])
%dW = @(x,a) M_z(x)*dM_z(x,a);
dW(x,a) = M_z(x)*dM_z(x,a);
dW_int(a) = 1/E_I.*(int(dW(x,a), x, 0 ,L));
dW_ext = dP*u_a;
eqn = dW_int == dW_ext;
Solve eqn for u_a
sol = solve(eqn, u_a, 'ReturnConditions',true)
We see six possible solutions depending on the value of a, all with dP ~= 0.
sol.conditions
Assume that a is bracketed by 10 < a < 20 because there was a comment upstream that a = 18. Assume that dP ~= 0
assumeAlso(10 < a < 20);
assumeAlso(dP ~= 0);
sol = solve(eqn, u_a, 'ReturnConditions',true);
Under those assumptions, we have one solution for u_a
sol.u_a
Differentiate u_a wrt a
u_a_diff = diff(sol.u_a, a)
And solve for a
eqn2 = u_a_diff == 0;
sol_a = solve(eqn2, a);
a_solve = double(sol_a)
3 个评论
Milan
2022-10-30
Hard to know what's going on w/o seeing the actual code that generated the error message. See below for what I think the goal is. Not that u_a_1 is not defined in the original code, so I had to guess what you're trying to do.
syms d_F ;
syms d_F_c;
syms dM_B;
%syms Theta_B;
syms x real;
syms u_mid;
syms u_c;
syms u_max;
syms a;
%syms w_x;
syms M_z(x);
syms dM_z(x);
syms dP;
syms u_a;
%syms n;
%Q = 2.5; %kip/ft
L = 30.; % ft
%a = 18; %assume maximum deflection at 15', 20'
E = 29000.*12^2; %ksi to kips/ft2
I = 1890./(12)^4; %ft^4
E_I = E.*I; %kip*ft^2;
%% deflection at mid of AB
w_x = @(x) 1*x/(L/3)*heaviside(x-L/3)*(1-heaviside(x-2*L/3));
w_x_1 = w_x(x);
eq1 = diff(M_z(x),2) == -w_x_1;
sol1(x) = dsolve(eq1, M_z(0)==-150, M_z(L) == 250, 'IgnoreAnalyticConstraints',true);
M_z(x) = sol1(x);
dM_z(x,a) = (dP*x*(L-a)/L)-dP.*(x-a).*heaviside(x-a);
dW(x,a) = M_z(x)*dM_z(x,a);
dW_int(a) = 1/E_I.*(int(dW(x,a), x, 0 ,L));
dW_ext = dP*u_a;
eqn = dW_int == dW_ext;
assumeAlso(10 < a < 20);
assumeAlso(dP ~= 0);
sol = solve(eqn, u_a, 'ReturnConditions',true);
Define u_a
u_a(a) = sol.u_a;
u_a_diff = diff(u_a, a);
eqn2 = u_a_diff == 0;
a_solve = solve(eqn2, a);
a_solve = double(a_solve)
M_max = double(M_z(a_solve))
%u_max = u_a_1(a_solve) % u_a_1 is not defined!
u_max = double(u_a(a_solve))
Milan
2022-10-30
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Programming 的更多信息
Translated by 
