Problem using the bvp4c solver. First derivate seems to be constant on the graphic
1 view (last 30 days)
Show older comments
Daniel Clemente
on 21 Apr 2014
Commented: Daniel Clemente
on 22 Apr 2014
Hello there, I have implemented a code to simulate the deflection of a beam. It works fine but physically there seems to be a problem with my non-linear solution. I wonder if the bvp4c solver is considering my first derivate correctly. It only appears on the non-linear solution so the graphics from the linear(in blue) and non linear(in red) should be somehow different, but the red looks as if they were multiplied by a constant.
It may be a syntax problem but I have tried to find it for quite a while and can´t. My question is very simple: is the "Non linear solver" correctly considering that y(2) is my dy/dx? The code is below:
function 2modificado= 2modificado(qo,L)
global qmax Lg EI;
qmax = qo;
Lg = L;
Em=200*10^9
%This part below is not relevant for the problem
for i=1:3
if i==1
b=0.055;
h=0.1;
Im=b*h^3/12
end
if i==2
Im=2.5175*10^(-6)
end
if i==3
re=0.065;
ri=0.060;
Im=pi/4*((re^4)-(ri^4))
end
EI=(Em*Im)*(1000^2)
x0 = [0 0];
%Längenintervall
xspan = [0 L];
% This part below is:
%Call Solver -> Linear
solinit = bvpinit(linspace(0,L,10),[0 0]);
solL=bvp4c(@bvp4odeL, @bvp4bcL, solinit);
xint=linspace(0,L)
Sxint=deval(solL,xint)
%Plot result
figure
plot(xint,Sxint(1,:))
if i==1 title('y(x) in a Retangular Profile')
end
if i==2 title('y(x) in a Iprofile(IPE 100)')
end
if i==3 title('y(x) in a Round hollow section')
end
set(gca,'yDir','reverse')
xlabel('[mm]');
ylabel('[mm]');
hold on;
%Plot result
plot(xint,Sxint(1,:));
%Call Solver Non-linear
solinitNL = bvpinit(linspace(0,L,10),[0 0]);
solNL=bvp4c(@bvp4odeNL, @bvp4bcNL, solinitNL);
xintNL=linspace(0,L)
SxintNL=deval(solNL,xint)
plot(xintNL,SxintNL(1,:),'-r')
set(gca,'yDir','reverse')
end
return
%---------------------------------------------------------------
function dydx = bvp4odeL(x,y)
global qmax Lg EI;
dydx = [y(2) (qmax*x*(x-Lg)/(2*EI))]
function res=bvp4bcL(ya,yb)
res=[ya(1) yb(1)]
function dydx = bvp4odeNL(x,y)
global qmax Lg EI;
dydx = [y(2) (qmax*x*(x-Lg)/(2*EI))*(1+y(2)^2)^1.5]
function res=bvp4bcNL(ya,yb)
res=[ya(1) yb(1)]
0 Comments
Accepted Answer
Sean de Wolski
on 21 Apr 2014
One thing sticks out: dydx should be a column vector
dydx = [y(2); (qmax*x*(x-Lg)/(2*EI))*(1+y(2)^2)^1.5]
Otherwise, try putting a function in there that you know for sure should make it highly non-linear.
More Answers (0)
See Also
Categories
Find more on Boundary Value Problems in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!