hi Roman,
' However, it may be possible to solve using a different method and obtain more accurate results. You should try to combine steps, so that less error accumulates from one step to another. This may make the problem harder. If you have exact data instead of floating point numbers that might help as well. When I have a minute I'll take a look at this and try to offer you some better suggestions.'
Thanks for this and look forward to your suggestions.
Also, there was another question in my first post. I need to switch between tension and compression cases. And the switching criteria is that 'the equilibirium does not produce any solutions'. Can I program this in Maple so that it switches automatically between compression and tension cases?
Many thanks
Kenny
hi Roman,
' However, it may be possible to solve using a different method and obtain more accurate results. You should try to combine steps, so that less error accumulates from one step to another. This may make the problem harder. If you have exact data instead of floating point numbers that might help as well. When I have a minute I'll take a look at this and try to offer you some better suggestions.'
Thanks for this and look forward to your suggestions.
Also, there was another question in my first post. I need to switch between tension and compression cases. And the switching criteria is that 'the equilibirium does not produce any solutions'. Can I program this in Maple so that it switches automatically between compression and tension cases?
Many thanks
Kenny
Hi Roman and Acer
Thanks for such quick replies.
The problem is here:
> solns := simplify(solve({eq1,eq2,eq3,eq4},{A1c,A2c,A3c,A4c}));
You are trying to symbolically solve a non-linear system with coefficients order 10^15 and one free parameter. Maple doesn't have a method for this (assuming there even is one) and as you can see the result is garbage:
> simplify(subs(solns, {eq1,eq2,eq3,eq4}));
I didnt realise this before. I just simply expected Maple to give me the correct answers. So what you are suggesting is that Maple is incapable of meeting my needs?
However, Using the solution of the axial force obtained from the equilibrium equation, i checked {eq1,eq2,eq3,eq4} and they were all satisfied.
> restart;
> Digits:=15:
Input
> T1:=50:T2:=100:b:=200:d:=400:l:=2000:E0:=200000:q:=-10:
> A:=b*d;k:=E0*A/2/l;r:=E0*(b*(d^3)/12)/2/l;rho:=(T1-T2)/d;alpha:=1.2*10^(-5):delta:=alpha*rho:
> rL:=r;rR:=r;kL:=k;kR:=k;
> Tc:=(T1+T2)/2;
> eta:=1+T/(2000*ln(T/1100));
Find y_bar, EA, EI,r,lamda
> y_bar:=evalf(d/3*(subs(T=T1,eta)+2*subs(T=T2,eta))/(subs(T=T1,eta)+subs(T=T2,eta)));
> EB:=evalf(int(y*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));
> EA:=evalf(E0*b*int((1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));
> EI:=evalf(E0*b*int(y^2*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));;
> M:=b*d^3/12+b*d*(d/2-y_bar)^2;
> epsilon_c:=evalf(alpha*subs(y=0,rho*y+T1-y_bar*rho));
> v:=A1c*cos(sqrt(N/EI)*z)+A2c*sin(sqrt(N/EI)*z)+A3c+A4c*z-(q*(z^2)/(2*N));
> eq1:=subs(z=l,v)=0;
> eq2:=subs(z=-l,v)=0;
> eq3:=EI*(subs(z=l,diff(v,z$2))+delta)+rR*subs(z=l,diff(v,z))=0;
> eq4:=-EI*(subs(z=-l,diff(v,z$2))+delta)+rL*subs(z=-l,diff(v,z))=0;
> solns:=simplify(solve({eq1,eq2,eq3,eq4},{A1c,A2c,A3c,A4c}));
> v_t:=subs(solns,A1c)*cos(sqrt(N/EI)*z)+subs(solns,A2c)*sin(sqrt(N/EI)*z)+subs(solns,A3c)+subs(solns,A4c)*z-(q*(z^2)/(2*N)):
> f:=1/4/l*int((diff(v_t,z))^2,z=-l..l);
> eq:=(-N/EA)+epsilon_c-N/(2*l)*(1/kL+1/kR)-f;
> N_sol:=fsolve(eq,N);
> v_sols:=evalf(subs(N=N_sol,v_t));
Check kinetic and static boundary equations eq1, eq2, eq3, eq4
> simplify(subs(z=-l,v_sols));simplify(subs(z=l,v_sols));
> simplify(EI*(subs(z=l,diff(v_sols,z$2))+delta)+rR*subs(z=l,diff(v_sols,z)));
> simplify(-EI*(subs(z=-l,diff(v_sols,z$2))+delta)+rL*subs(z=-l,diff(v_sols,z)));
So I have a deflected shape which satisfies equilibrium and all boundary conditions. It has to be correct, hasnt it?
I can approach the whole thing numerically but it would make little sense since I have an *exact* solution.
Cheers,
Kenny
Hi Roman and Acer
Thanks for such quick replies.
The problem is here:
> solns := simplify(solve({eq1,eq2,eq3,eq4},{A1c,A2c,A3c,A4c}));
You are trying to symbolically solve a non-linear system with coefficients order 10^15 and one free parameter. Maple doesn't have a method for this (assuming there even is one) and as you can see the result is garbage:
> simplify(subs(solns, {eq1,eq2,eq3,eq4}));
I didnt realise this before. I just simply expected Maple to give me the correct answers. So what you are suggesting is that Maple is incapable of meeting my needs?
However, Using the solution of the axial force obtained from the equilibrium equation, i checked {eq1,eq2,eq3,eq4} and they were all satisfied.
> restart;
> Digits:=15:
Input
> T1:=50:T2:=100:b:=200:d:=400:l:=2000:E0:=200000:q:=-10:
> A:=b*d;k:=E0*A/2/l;r:=E0*(b*(d^3)/12)/2/l;rho:=(T1-T2)/d;alpha:=1.2*10^(-5):delta:=alpha*rho:
> rL:=r;rR:=r;kL:=k;kR:=k;
> Tc:=(T1+T2)/2;
> eta:=1+T/(2000*ln(T/1100));
Find y_bar, EA, EI,r,lamda
> y_bar:=evalf(d/3*(subs(T=T1,eta)+2*subs(T=T2,eta))/(subs(T=T1,eta)+subs(T=T2,eta)));
> EB:=evalf(int(y*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));
> EA:=evalf(E0*b*int((1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));
> EI:=evalf(E0*b*int(y^2*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));;
> M:=b*d^3/12+b*d*(d/2-y_bar)^2;
> epsilon_c:=evalf(alpha*subs(y=0,rho*y+T1-y_bar*rho));
> v:=A1c*cos(sqrt(N/EI)*z)+A2c*sin(sqrt(N/EI)*z)+A3c+A4c*z-(q*(z^2)/(2*N));
> eq1:=subs(z=l,v)=0;
> eq2:=subs(z=-l,v)=0;
> eq3:=EI*(subs(z=l,diff(v,z$2))+delta)+rR*subs(z=l,diff(v,z))=0;
> eq4:=-EI*(subs(z=-l,diff(v,z$2))+delta)+rL*subs(z=-l,diff(v,z))=0;
> solns:=simplify(solve({eq1,eq2,eq3,eq4},{A1c,A2c,A3c,A4c}));
> v_t:=subs(solns,A1c)*cos(sqrt(N/EI)*z)+subs(solns,A2c)*sin(sqrt(N/EI)*z)+subs(solns,A3c)+subs(solns,A4c)*z-(q*(z^2)/(2*N)):
> f:=1/4/l*int((diff(v_t,z))^2,z=-l..l);
> eq:=(-N/EA)+epsilon_c-N/(2*l)*(1/kL+1/kR)-f;
> N_sol:=fsolve(eq,N);
> v_sols:=evalf(subs(N=N_sol,v_t));
Check kinetic and static boundary equations eq1, eq2, eq3, eq4
> simplify(subs(z=-l,v_sols));simplify(subs(z=l,v_sols));
> simplify(EI*(subs(z=l,diff(v_sols,z$2))+delta)+rR*subs(z=l,diff(v_sols,z)));
> simplify(-EI*(subs(z=-l,diff(v_sols,z$2))+delta)+rL*subs(z=-l,diff(v_sols,z)));
So I have a deflected shape which satisfies equilibrium and all boundary conditions. It has to be correct, hasnt it?
I can approach the whole thing numerically but it would make little sense since I have an *exact* solution.
Cheers,
Kenny
Posted this by mistakes please ignore
Posted this by mistakes please ignore
Hi roman, Thanks for your replies. The equilibirum equation is highly non-linear and cumbersome at best. The equation also changes from cases to case s. I attached here the maple working, which includes all steps that lead to the final equilibrium equation. > restart; > Digits:=15: Input > T1:=50:T2:=100:b:=200:d:=400:l:=2000:E0:=200000:q:=-10: > A:=b*d;k:=E0*A/2/l;r:=E0*(b*d^3/12)/2/l;rho:=(T1-T2)/d;alpha:=1.2*10^(-5):delta:=alpha*rho; > rL:=r;rR:=r;kL:=k;kR:=k; > Tc:=(T1+T2)/2; > eta:=1+T/(2000*ln(T/1100)); > Section properties > y_bar:=evalf(d/3*(subs(T=T1,eta)+2*subs(T=T2,eta))/(subs(T=T1,eta)+subs(T=T2,eta))); > EB:=evalf(int(y*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar)); > EA:=evalf(E0*b*int((1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar)); > EI:=evalf(E0*b*int(y^2*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));; > M:=b*d^3/12+b*d*(d/2-y_bar)^2; > epsilon_c:=evalf(alpha*subs(y=0,rho*y+T1-y_bar*rho)); > Solve for the coefficients A1,A2,A3,A4 > > v:=A1c*(EI/N)*exp(sqrt(N/EI)*z)+A2c*(EI/N)*exp(-sqrt(N/EI)*z)+A3c+A4c*z-(q*(z^2)/(2*N)); > eq1:=subs(z=l,v)=0; > eq2:=subs(z=-l,v)=0; > eq3:=EI*(subs(z=l,diff(v,z$2))+delta)+rR*subs(z=l,diff(v,z))=0; > eq4:=-EI*(subs(z=-l,diff(v,z$2))+delta)+rL*subs(z=-l,diff(v,z))=0; > solns:=simplify(solve({eq1,eq2,eq3,eq4},{A1c,A2c,A3c,A4c})); > v_t:=subs(solns,A1c)*(EI/N)*exp(sqrt(N/EI)*z)+subs(solns,A2c)*(EI/N)*exp(-sqrt(N/EI)*z)+subs(solns,A3c)+subs(solns,A4c)*z-(q*(z^2)/(2*N)): > > > Equilibrium equation > f:=1/4/l*int((diff(v_t,z))^2,z=-l..l); > eq:=(N/EA)+epsilon_c+N/(2*l)*(1/kL+1/kR)-f=0; > N_sol:=fsolve(eq,N); The final equation (eq=0) is highly non-linear and the function itself appears to be discontinuous. Many thanks, Kenny
Hi roman, Thanks for your replies. The equilibirum equation is highly non-linear and cumbersome at best. The equation also changes from cases to case s. I attached here the maple working, which includes all steps that lead to the final equilibrium equation. > restart; > Digits:=15: Input > T1:=50:T2:=100:b:=200:d:=400:l:=2000:E0:=200000:q:=-10: > A:=b*d;k:=E0*A/2/l;r:=E0*(b*d^3/12)/2/l;rho:=(T1-T2)/d;alpha:=1.2*10^(-5):delta:=alpha*rho; > rL:=r;rR:=r;kL:=k;kR:=k; > Tc:=(T1+T2)/2; > eta:=1+T/(2000*ln(T/1100)); > Section properties > y_bar:=evalf(d/3*(subs(T=T1,eta)+2*subs(T=T2,eta))/(subs(T=T1,eta)+subs(T=T2,eta))); > EB:=evalf(int(y*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar)); > EA:=evalf(E0*b*int((1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar)); > EI:=evalf(E0*b*int(y^2*(1+(rho*y+T1-rho*y_bar)/(2000*ln((rho*y+T1-rho*y_bar)/1100))),y=(y_bar-d)..y_bar));; > M:=b*d^3/12+b*d*(d/2-y_bar)^2; > epsilon_c:=evalf(alpha*subs(y=0,rho*y+T1-y_bar*rho)); > Solve for the coefficients A1,A2,A3,A4 > > v:=A1c*(EI/N)*exp(sqrt(N/EI)*z)+A2c*(EI/N)*exp(-sqrt(N/EI)*z)+A3c+A4c*z-(q*(z^2)/(2*N)); > eq1:=subs(z=l,v)=0; > eq2:=subs(z=-l,v)=0; > eq3:=EI*(subs(z=l,diff(v,z$2))+delta)+rR*subs(z=l,diff(v,z))=0; > eq4:=-EI*(subs(z=-l,diff(v,z$2))+delta)+rL*subs(z=-l,diff(v,z))=0; > solns:=simplify(solve({eq1,eq2,eq3,eq4},{A1c,A2c,A3c,A4c})); > v_t:=subs(solns,A1c)*(EI/N)*exp(sqrt(N/EI)*z)+subs(solns,A2c)*(EI/N)*exp(-sqrt(N/EI)*z)+subs(solns,A3c)+subs(solns,A4c)*z-(q*(z^2)/(2*N)): > > > Equilibrium equation > f:=1/4/l*int((diff(v_t,z))^2,z=-l..l); > eq:=(N/EA)+epsilon_c+N/(2*l)*(1/kL+1/kR)-f=0; > N_sol:=fsolve(eq,N); The final equation (eq=0) is highly non-linear and the function itself appears to be discontinuous. Many thanks, Kenny