MATLAB: Where are the bugs for this ODE finite difference problem that solve using Newton Raphson method

finite differnce methodode

Could someone provide me help to solve this Euler-Bernoulli beam equation by using finite difference method and Newton Raphson please.

With boundary value of y(0) = 0 and dy/ds(L) = 0

I am continually getting answers that are nowhere near the results from the bvp4c command.

Best Answer

Here is my attempt

clc,clear
N = 10;   % 10 Set parameters
L = 16*0.3048; %meter
b = 19.625 * 0.0254; d = 1.625 * 0.0254;
p = 26.5851; I = (b*(d^3))/12; h = L/N; F = 0.7436; E = 68.9e6;
alpha = (h^2)/(E*I);
w = p*b*d;
fode = @(s,f) [f(2); -1/E/I*(w*(L-s)+F)*cos(f(1))];
fbound= @(ya,yb) [ya(1)-0; yb(2)-0];
ss = linspace(0,5);
finit = [0 0];
solinit = bvpinit(ss,finit);
sol = bvp4c(fode,fbound,solinit);
plot(sol.x,sol.y(1,:))

Related Solutions

MATLAB: Help regarding use of BVP4C in solving an third order ODE.

Hi,

this works:

init=bvpinit(linspace(0,4,8),[0 0 0]);
sol=bvp4c(@rhs_bvp,@bc_bvp,init);
plot(sol.x,sol.y(1,:),sol.x,sol.y(2,:),sol.x,sol.y(3,:))
function bc=bc_bvp(yl,yr) %ya(1), ya(2), yb(2)
bc=[yl(1); yl(2); yr(2)-1];
end
function rhs=rhs_bvp(x,y)
rhs=[y(2); y(3); y(2)-y(1)*y(3)-1];
end

Best regards

Stephan

MATLAB: If function error, here i want MatLab to find the values of V5L that are between 1 and 4.5 and name it V, and use it in another upcoming equation

This works for me:

Q = 3.3;
D5L = [12.75*0.0254 14*0.0254 16*0.0254 18*0.0254 20*0.0254 22*0.0254 24*0.0254 26*0.0254 28*0.0254 30*0.0254 32*0.0254 34*0.0254 36*0.0254 38*0.0254 40*0.0254 42*0.0254 44*0.0254 46*0.0254 48*0.0254 52*0.0254 56*0.0254 60*0.0254 64*0.0254 68*0.0254 72*0.0254 76*0.0254 80*0.0254];
V5L = ((4.*Q)./(pi.*D5L.^2));
V = V5L((V5L > 1) & (V5L < 4.5));
Do = sqrt((4.*Q)./(pi.*V))

producing:

Do =
    1.0160    1.0668    1.1176    1.1684    1.2192    1.3208    1.4224    1.5240    1.6256    1.7272    1.8288    1.9304    2.0320

EDIT — (27 Jan 2019 at 15:09)

This runs for me without error (in R2018b):

Q=3.3;
%gamma=8.43 kN/m3 = 859.62076503 Kg/m3 (*1000/9.81)
rouh=859.62076503;
g=9.81;
Gamma=(rouh.*g);
Zs=0;
Ze=0;
L=1285000;
Ks=(0.000045);
Er=(1E-8);
viscosity=(3.83.*1E-3);
% material grade X70
strength= 482000000;
%API 5L
% % D5L=[12.75*0.0254 14*0.0254 16*0.0254 18*0.0254 20*0.0254 22*0.0254 24*0.0254 26*0.0254 28*0.0254 30*0.0254 32*0.0254 34*0.0254 36*0.0254 38*0.0254 40*0.0254 42*0.0254 44*0.0254 46*0.0254 48*0.0254 52*0.0254 56*0.0254 60*0.0254 64*0.0254 68*0.0254 72*0.0254 76*0.0254 80*0.0254]
% % V5L=((4.*Q)./(pi.*D5L.^2));
% % if V=V5L(V5L > 1 & abs(V5L - 4.5<1e-5));
% % end
% % Do=((4.*Q)./(pi.*V)).^(0.5)
D5L = [12.75*0.0254 14*0.0254 16*0.0254 18*0.0254 20*0.0254 22*0.0254 24*0.0254 26*0.0254 28*0.0254 30*0.0254 32*0.0254 34*0.0254 36*0.0254 38*0.0254 40*0.0254 42*0.0254 44*0.0254 46*0.0254 48*0.0254 52*0.0254 56*0.0254 60*0.0254 64*0.0254 68*0.0254 72*0.0254 76*0.0254 80*0.0254];
V5L = ((4.*Q)./(pi.*D5L.^2));
V = V5L((V5L > 1) & (V5L < 4.5));
Do = sqrt((4.*Q)./(pi.*V));
Re=(V.*Do./viscosity);
Fo=0.01;
for I=0:1:1E6;
    Fn=(1./(-4.*log10((Ks./(3.71.*Do))+((1.26)./(Re.*sqrt(Fo)))))).^2;
    E=abs((Fn-Fo)/Fn);
    if  E<=Er;
        Fn=Fn(end);
        display(E)
        break
    end
    if E>Er;
        Fo=Fn;
    end
end
Hloss=((Fn.*L.*Q.^2)./(12.*Do.^5));
Hpump=Hloss+Ze-Zs;
%calculating Hloss every 128500M
LC=1:128500:1285000;
for i=1:1:length(LC);
    Hpump=Hloss+Ze-Zs;
    PRESSURE=rouh.*g.*Hpump.*1E-5;
    numberofpumps=ceil(PRESSURE./90);
    pressureperpump=(PRESSURE./numberofpumps);
    HlossN=((Fn.*LC(i).*Q.^2)./(12.*sqrt(Do)));
    TEL(i,:)=(Hpump-HlossN+Zs);
    HGL(i,:)=(TEL(i,:)-((V.^2)./(2.*g)));
    PE=Gamma.*(HGL-Ze);
    T=((PE(i).*Do)./(2.*strength));
end
for ii=1:1:length(Do)
    subplot(4,4,ii)
    plot(LC,TEL(:,ii))
    hold on
    plot(LC,HGL(:,ii))
    xlabel('LC')
    ylabel('HGL & TEL')
end

Best Answer

Related Solutions

MATLAB: Help regarding use of BVP4C in solving an third order ODE.

MATLAB: If function error, here i want MatLab to find the values of V5L that are between 1 and 4.5 and name it V, and use it in another upcoming equation

Related Question