Ode45 solves an equation that containing a definite integral term
조회 수: 18 (최근 30일)
이전 댓글 표시
Hi, i have a problem on the following equation when solved by Ode45, which contains a definite integral term. I dont know how to transform it so that it can be solved by ode45.
can someone help me ?
Thank you in advance
댓글 수: 2
답변 (2개)
Torsten
2019년 9월 6일
편집: Torsten
2019년 9월 6일
function main
a = ...;
b = ...;
c = ...;
d = ...;
u0 = 1;
usol = fzero(@(u)fun(u,a,b,c,d),u0);
fun_ode = @(t,y)[y(2);y(3);y(4);-usol*y(3)-y(1)^2];
y0 = [a;b;c;d];
tspan = [0 1];
[T,Y] = ode45(fun_ode,tspan,y0);
plot(T,Y)
end
function res = fun(u,a,b,c,d)
fun_ode = @(t,y)[y(2);y(3);y(4);-u*y(3)-y(1)^2;y(2)^2];
y0 = [a;b;c;d;0];
tspan = [0 1];
[T,Y] = ode45(fun_ode,tspan,y0);
res = Y(end,5)-u;
end
댓글 수: 4
Arpan Laha
2021년 5월 20일
Dear Torsten,
Please correct me if I am wrong.
fun_ode = @(t,y)[y(2);y(3);y(4);-u*y(3)-y(1)^2;y(2)^2]; solves the function provided by Xuan replacing the definite integral term by indefinite integral. The y we got as solution from res will not be the same if solved by taking definite integral. Lets term this as Yindf. lets term the actual solution as Ydef. For obvious reason Yindf~=Ydef. Then how can we substitute the value of u into the actual equation ? Thanks
Torsten
2021년 5월 20일
편집: Torsten
2021년 5월 20일
For a given value of u, you determine the solution y of the differential equation y''''+u*y''+y^2=0 in the interval [0;1] (components 1-4 of fun_ode).
Simultaneously, you integrate y'^2 in the interval [0;1] (component 5 of fun_ode).
Usually, u will not be equal to component 5 of fun_ode, evaluated at x=1.
So, fzero must be used to adjust u such that the two numbers become equal.
yuan bin
2023년 1월 11일
편집: yuan bin
2023년 1월 11일
refer idsolver
IDSOLVER: A general purpose solver for nth-order integro-differential equations
My code is below:
[x ,nominales] = ode45(@model0,[0,1],[a,b,c,d]);
nominales = [x, nominales];
Tol = 1e-8;
st.TolQuad = 1e-8;
% Iterative solution
error = 1e3;
iteration = 1;
fprintf(' Error convergence\n ');
fprintf(' ================= \n');
fprintf(' Iteration Error \n');
while error > Tol
st.nominales = nominales;
S = ode45(@(x,y)model(x,y,st),[0,1],[a,b,c,d]);
x = nominales(:,1);%time points
R = deval(S, x);
y = R.';
error = sum((y(:,1)-nominales(:,2)).^2);
fprintf(' %4i %8.2e\n',[iteration error]);
alpha = 0.5;
nominales(:,2) = (1-alpha)*nominales(:,2)+alpha*y(:,1);
nominales(:,1) = x;
iteration = iteration+1;
end
warning('on')
plot(x,y);
% legend('y','dy','d2y','d3y');
disp('done');
function dy = model0(x,y)
% Initial guess generator
dy = zeros(4,1);
% y-> [y,dy d2y d3y];
dy(1) = y(2);
dy(2) = y(3);
dy(3) = y(4);
dy(4) = -(y(1)^2+ y(3)* 0); % 0 for initiation value
end
function dy = model(x,y,st)
nominales = st.nominales;
TolQuad = st.TolQuad;
% Interpolation step
ys = @(s) interp1(nominales(:,1),nominales(:,3),s);
% Integro-differential equation
dy = zeros(4,1);
dy(1) = y(2);
dy(2) = y(3);
dy(3) = y(4);
dy(4) = -(y(1)^2 + y(3) *quadl(@(s) ys(s).*ys(s) ,0,1,TolQuad));
end
댓글 수: 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!