# Manipulate an integral with variable in limit of the solution of DDEs

I have the following problem: I solved DDEs with some parameters which have two solutions: r and s. Now, I need to calculate and plot (in Manipulate) the following integral: And there are problems. Which function to use to calculate the integral? NIntegrate doesn't work because "x in limit is unexpected". If I try to use Integrate I obtain the following error:

Integral and error estimates are 0 on all integration subregions. Try increasing the value of the MinRecursion option. If value of integral may be 0, specify a finite value for the AccuracyGoal option.

But Integrate does not have MinRecursion and AccuracyGoal options. Also, I'm not sure if I put functions r and s in the integral correctly (see code below, these functions in the final Manipulate should be of h as in prob function?). How do I make Manipulate work?

CODE: Auxiliary functions:

hill[x_, n_] := x^n/(1 + x^n)
input = 0.4
Sort[{r, s} /.
NSolve[{input - r + s * hillweight * hill[r * s, n] == 0,
input - s + r * hillweight * hill[r * s, n] == 0, r >= 0,
s >= 0}, {r, s}, Reals]]


My DDEs solver:

modelsolver[σ_, τ_, hillweight_, end_, sigend_] :=
NDSolve[{r'[t] ==
input  + Piecewise[{{σ, t < sigend}, {0, t >= 0.5}}] -
r [t - τ] + s[t] * hillweight * hill[r[t] * s[t], 2],
s'[t] ==
input - s[t - τ] + r[t] * hillweight * hill[r[t] * s[t], 2],
r[t /; t <= 0] == First[steadypar[2, hillweight]][],
s[t /; t <= 0] == First[steadypar[2, hillweight]][]}, {r,
s}, {t, 0, end}]


Function that should calculate integral:

probtest[x_, β_, r_, s_] := 1/(1 + E^(-β * Integrate[r - s, {h, 0, x}]))


Alternative version:

prob[x_?NumericQ, β_, r_, s_] :=
1/(1 + E^(-β *
NIntegrate[r - s, {h, 0, x}, MinRecursion -> 100,
AccuracyGoal -> 100]))


Final Manipulate:

Simulation =
Manipulate[
Plot[Evaluate[({probtest[x, β, r[h], s[h]],
1 - probtest[x, β, r[h], s[h]]} /.
modelsolver[σ, τ, hillweight, end, 0.5][])], {x,
0, end}, PlotRange -> {0, 1}, MaxRecursion -> 15 ], {τ, 0,
2, 0.1}, {σ, 0 , 4, 0.01}, {hillweight, 0.85, 1,
0.01}, {end, 10, 200, 15}, {β, 1, 1.3, 0.05}]


Include the integral into your DDE. Here is a version of your NDSolve that does this:

modelsolver[σ_, τ_, hillweight_, end_, sigend_] := NDSolve[
{
r'[t] == input  + Piecewise[{{σ, t < sigend}, {0, t >= 0.5}}] - r[t - τ] + s[t] * hillweight * hill[r[t] * s[t], 2],
s'[t] == input - s[t - τ] + r[t] * hillweight * hill[r[t] * s[t], 2],
r[t /; t <= 0] == First[steadypar[2, hillweight]][],
s[t /; t <= 0] == First[steadypar[2, hillweight]][],
int'[t] == r[t] - s[t],
int[t /; t<=0] == 0
},
{r, s, int},
{t, 0, end}
]


For example:

sol = modelsolver[2, 1, .85, 10, .5];
Plot[int[t] /. sol, {t, 0, 10}, PlotRange->All] 