Problems with ode23s

2 次查看（过去 30 天）
Pablo R
Pablo R 2017-9-26
回答： Jan 2017-9-26
I have this warning on matlab
Warning: Failure at t=4.447832e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.580187e-14) at time t.
the important code is this and i dont know how to solve it:
t0 = 0.001;
tf = 144;
Opciones = odeset('RelTol', 1e-4, 'AbsTol', 1e-7);
[t,X] = ode23s('SedimenEstNN', [t0,tf], [X(1) X(2) X(3) X(4)...
X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12) X(13) X(14)...
X(15) X(16) X(17)], Opciones);
function dX=SedimenEstNN(t,X)
global A H K N P Va be
U = SedimenEntNN(t);
Qd = P*U(1);
Eda = Va*U(2)/(24*Qd*X(16));
dX(17) = Eda;
Edad =X(17)/t;
%more code...
function U=SedimenEntNN(t)
global A H K N P Va be
A = 200; H = 2; K = 7e-13;
N = 6; be = 2.2; Va = 600;
w = 2*pi/24*(t-10);
R = 0.75;
P = 0.0156;
if t>118
U(1) = 100; U(2) = 2500; U(3) = R*U(1);
else
U(1) = 100 + 50*sin(w);
U(2) = 2500 + 500*sin(w);
U(3) = R*U(1);
end
thank you before all and if i need to put more code tell me please
  1 个评论
Jan
Jan 2017-9-26
[X(1) X(2) X(3) X(4)...
X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12) X(13) X(14)...
X(15) X(16) X(17)]
?? What about:
X(1:17)

请先登录，再进行评论。

请先登录，再回答此问题。

回答（2 个）

Torsten
Torsten 2017-9-26
编辑：Torsten 2017-9-26
It's hard to tell why ode23s fails for your problem.
Some options how to proceed:
1. Make sure that you transposed the dX-array at the end of function "SedimenEstNN".
2. The vector of initial guesses [X(1) X(2) X(3) X(4) X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12) X(13) X(14) X(15) X(16) X(17)] should be a column vector, not a row vector. Did you give values to the X(i)'s before you called ode23s ?
3. Check the solution up to the time where the solver fails.
4. Check the dX values at the time where the solver fails.
5. Choose smaller tolerances.
6. Try a stiff solver (e.g. ode15s).
Good luck.
Best wishes
Torsten.
  1 个评论
Pablo R
Pablo R 2017-9-26
Hi Torsten,
i tried all of that but still having the same problem. i show you by the way,
1)
dX = [dX(1); dX(2); dX(3); dX(4); dX(5); dX(6); dX(7); dX(8);...
dX(9); dX(10); dX(11); dX(12); dX(13); dX(14); dX(15);...
dX(16); dX(17)];
2)
X(1) = 100; X(2) = 75; X(3) = 2; X(4) = 2500;
X(5) = 100; X(6) = 75; X(7) = 3; X(8) = 100;
X(9) = 75; X(10) = 3; X(11) = 2500; X(12) = 2500;
X(13) = 223; X(14) = 175; X(15) = 30; X(16) = 6000;
X(17) = 1e-3;
5) same problem
6) same problem
i will see the values at that time
but thank you so much for your answer!

请先登录，再进行评论。

Jan
Jan 2017-9-26
[X(1) X(2) X(3) X(4)...
X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12) X(13) X(14)...
X(15) X(16) X(17)]
?? What about:
X(1:17)
Do not provide the function to be integrated as string 'SedimenEstNN', but as function handle @SedimenEstNN. This is the preferred method since Matlab 6.5 (over 15 years now...).
The function to be integrated is not smooth at t=118. This is a really bad idea and Matlab's integrators are not designed to handle such jumps. See http://www.mathworks.com/matlabcentral/answers/59582#answer_72047
Use the reached time as endpoint and look at the trajectory:
t0 = 0.001;
tf = 4.447832e+00 - 10 * eps(4.447832);
Opciones = odeset('RelTol', 1e-4, 'AbsTol', 1e-7);
[t,X] = ode23s('SedimenEstNN', [t0,tf], [X(1) X(2) X(3) X(4)...
X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12) X(13) X(14)...
X(15) X(16) X(17)], Opciones);
plot(t, X); % Or X.'? I cannot remember this detail
Do you see a pole? If so, it is a mathematical problem and not a fault of the integrator.

类别

在 Help Center 和 File Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by