Bisection implementation in optimizing non linear ode solution
조회 수: 5 (최근 30일)
이전 댓글 표시
Hey guys. I have this ODE with the following bc's:
My code includes the shooting method and an intellectual guess of the derivative of rho at x=0 (transforming one BVP problems to one IVP problem), using the Implicit Euler method to write down the derivatives, solve the algebraic system with the Newton-Raphson method and make sure that I converge with infinity norm condition. I try to optimize my results by conducting bisection on the difference between a current value of rho to the its value at the bc that I chose (minimizing f1). This is the implicit Euler code:
function[x,y,f1,BISEC]=implicit_euler_method_DFT(y10,y20,x0,xf,dx,m,n)%i.g=-0.4394
x=x0:dx:xf;
y=zeros(2,length(x));
y(1,1)=y10;
y(2,1)=y20;
B=1;k=1;bc=0.25;
for i=1:1:length(x)-1
y(:,i+1)=y(:,i);
J=Jacovian(y(:,i));
g=[y(1,i+1)-dx*y(2,i+1)-y(1,i);
y(2,i+1)-dx*(2*B/k*y(1,i+1)+1/k*log(abs(y(1,i+1)))-B/k)-y(2,i)];
d=gausselemination(J,g);
f1=y(1,i+1)-bc;
BISEC=bisectionDFTbc(f1);
while (infinity_norm(d)>10^(-m))||(infinity_norm(g)>10^(-n))
y(:,i+1)=y(:,i+1)+d;
J=Jacovian(y(:,i+1));
g=[y(1,i+1)-dx*y(2,i+1)-y(1,i);
y(2,i+1)-dx*(2*B/k*y(1,i+1)+1/k*log(abs(y(1,i+1)))-B/k)-y(2,i)];
d=gausselemination(J,g);
f1=y(1,i+1)-bc;
BISEC=bisectionDFTbc(f1);
end
end
end
where x is the x tilde basically, y is rho tilde, y10 is the initial value for rho, y20 is the initial guess (-0.4394) for its derivative at x=0,x0=0, xf can be 1 or more, dx is to our desire, m and n are tolerances. The 2nd code for the bisection is:
function [x,iter_num]= bisectionDFTbc(a,b,n,m)
iter_num=0;%iteration number
x=(a+b)/2;
while sqrt((a-b)^2)>10^(-m)||sqrt(f1(x)^2)>10^(-n)
iter_num=iter_num+1;
x=(a+b)/2;
if f1(a)*f1(x)<0
b=x;
else
a=x;
end
end
x=(b+a)/2;
end
I have two main problems:
1.The y values go to infinity very fast - I don't get a stable solution.
2.Obviously I don't connect the two functions well.
I would appreciate your help.
Thanks,
Roi
댓글 수: 0
답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!