0
$\begingroup$

Bug introduced in 12.3 or earlier and persisting through 12.3.1 or later


I've checked the question about NDSolve::ndinid

This is the equation(Sorry for the complexity of this equation) and it's already written by Piecewise and on the right side as Michael E2 mentioned. However, the problem still exist.

μm=0.67;
ks=0.28;
x2max=2039;
x3max=939.5;
x4max=1026;
x5max=360.9;
x7max=80;
kgm=0.53;
kpm=0.14;

NR[expr_]:=Piecewise[{{1,expr>0},{0,expr<= 0}}];
μ = (μm*x2[t])/(
   x2[t] + ks) (1 - x7[t]/x7max) (1 - x2[t]/x2max) (1 - 
     x3[t - τ1]/x4max) (1 - x4[t]/x4max) (1 - x5[t]/x5max);
q20 = m2 + μ/y2 + dq2*x2[t]/(x2[t] + k2max);
q2 = (k1*x2[t])/(x2[t] + k2) + 
   k3 (x2[t] - x6[t - τ2])*NR[x2[t] - x6[t - τ2]];
q3 = l1*(k4*x8[t - τ4])/(x8[t - τ4] + k5) + 
   l2*k6*(x8[t - τ4] - x3[t - τ1])*
    NR[x8[t - τ4] - x3[t - τ1]];
q4 = m4 + μ*y4;
q5 = m5 + μ*y5;
l1=1;
l2=1;

cond = {x1[t /; t <= 0] == 0.102, x2[t /; t <= 0] == 418.2609, 
   x3[t /; t <= 0] == 0, x4[t /; t <= 0] == 0, x5[t /; t <= 0] == 0, 
   x6[t /; t <= 0] == 0, x7[t /; t <= 0] == 0, x8[t /; t <= 0] == 0};
pars={k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,m2,m4,m5,y2,y4,y5,dq2,k2max};
minPara = 
  Thread[pars -> {
     v1 -> 74.15520144717327, v2 -> 20.601402223057743`, 
      v3 -> 202.89918221914687`, v4 -> 124.79211581419787`, 
      v5 -> 2.4441045037299847`, v6 -> 47.960088480721915`, 
      v7 -> 9.272966935510649, v8 -> 55.015434647393015`, 
      v9 -> 4.369360890981919, v10 -> 2103.346395453384, 
      v11 -> 8.387912952673362, v12 -> 268.3929837388525, 
      v13 -> 33.53996289971859, v14 -> -3.224365183029838, 
      v15 -> 11.193999240900595`, v16 -> 33.05890696229195, 
      v17 -> 88.20748789196996, v18 -> -3.3751062489269783`, 
      v19 -> -3.673156600729279, v20 -> 3.8923944544765714`, 
      v21 -> -11.123011817112472`, v22 -> 31.8854734957005, 
      v23 -> 11.36136104963341, v24 -> 55.43326765393106, 
      v25 -> 9.250141726334826}[[;; , 2]]];

eqn = {
x1'[t] == μ*x1[t],
x2'[t] == -q2*x1[t],
x3'[t] == q3*x1[t],
x4'[t] == q4*x1[t],
x5'[t] == q5*x1[t],
x6'[t] == 
 1/k7*((k8*x2[t])/(x2[t] + k9) + 
     k10 (x2[t] - x6[t - τ2])*NR[x2[t] - x6[t - τ2]] - 
     q20) - μ*x6[t - τ2],
x7'[t] == (k11*x6[t - τ2])/(
  kgm (1 + x7[t - τ3]/k12) + x6[t - τ2]) - (
  k13*x7[t - τ3])/(
  kpm + x7[t - τ3] (1 + x7[t - τ3]/k14)) - μ*
   x7[t - τ3],
x8'[t] == (k13*x7[t - τ3])/(
  kpm + x7[t - τ3] (1 + x7[t - τ3]/k14)) - 
  l1*k15*x8[t - τ4]/(x8[t - τ4] + k16) - 
  l2*k17 (x8[t - τ4] - x3[t - τ1]) NR[
    x8[t - τ4] - x3[t - τ1]] - μ*x8[t - τ4]
   }/. minPara;

tauNDSolve[v1_, v2_, v3_, v4_] := 
  NDSolveValue[{eqn /. {τ1 -> v1, τ2 -> v2, τ3 -> 
       v3, τ4 -> v4}, cond}, {x1, x2, x3, x4, x5, x6, x7, 
    x8}, {t, 0, 7}];

The error is

NDSolveValue::ndinid: Initial condition Sign[x8\$...] is not on the range specified by the discrete variable NDSolve`s\$..

Is there any modification to make it work? Grateful for any suggestion!

$\endgroup$
3
  • 1
    $\begingroup$ There is probably some constant missing in the first DE, x1'[t] == *x1[t]. $\endgroup$
    – Domen
    Jul 30 '21 at 16:47
  • $\begingroup$ Sorry, it's μ * x1[t]. $\endgroup$ Jul 30 '21 at 17:49
  • $\begingroup$ I think you should report it to WRI. I’m not sure “delayed” events (such as represented by your Piecewise) are supported. There’s nothing about it in the docs that I’ve found that indicate whether it is or isn’t. I was unable to get even a simple system with a delayed event to work. If unsupported, better parsing might give an error message to that effect. $\endgroup$
    – Michael E2
    Jul 30 '21 at 23:36
1
$\begingroup$

The problem lies in your definition of NR. Apparently, Mathematica has some problems with these piecewise-like functions (see here and here).

Try changing it to:

NR[expr_ /; expr > 0] := 1;
NR[expr_ /; expr <= 0] := 0;

This works for me in Mathematica 12.3 - I can evaluate it at say $(1,1,1,1$).

tauNDSolve[1, 1, 1, 1]

Mathematica graphics

$\endgroup$
2
  • $\begingroup$ Oh, right, I didn't notice it stopped at 2.4 ... $\endgroup$
    – Domen
    Jul 30 '21 at 18:27
  • 1
    $\begingroup$ That the computation fails at about 2.4 has nothing to do with the original issue raised in this question. Instead, it fails because 0.14 + (1 - 0.3101385678220016` x7[t - v3]) x7[t - v3]` vanishes there, causing x7'[t] and x8[t]` to become singular. This is an inherent problem with the ODEs for the parameters chosen. Therefore, the answer above correctly identifies problem and provides an effective workaround. $\endgroup$
    – bbgodfrey
    Jul 31 '21 at 15:26
0
$\begingroup$

As @bbgodfrey pointed out, the DDE in tauNDSolve[1, 1, 1, 1] has a singularity at t = ~2.4. This is coincidentally near a discontinuity event due to NR[], so I mistakenly thought that was the problem. But with that cleared up, my original idea of turning off discontinuity processing works:

tauNDSolve[v1_, v2_, v3_, v4_] := 
 NDSolveValue[{eqn /. {\[Tau]1 -> v1, \[Tau]2 -> v2, \[Tau]3 -> 
      v3, \[Tau]4 -> v4}, cond}, {x1, x2, x3, x4, x5, x6, x7, x8}, {t,
    0, 7}, Method -> "DiscontinuityProcessing" -> False]

sol = tauNDSolve[1, 1, 1, 1]

Mathematica graphics

I still think the error message, which contains internal local symbols, should be reported to WRI. That shouldn't happen.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.