Passing RootSearch through NIntegrate?

Question

I am trying to integrate a function that is itself a function of RootSearch (from the Wolfram Library Archive). Specifically, for a given value of parameter t0, I am using RootSearch to calculate a particular value (FinalTime1). This value is then used as a parameter for Function1. I am then trying to integrate Function1 over all values of t0.

The code I have so far is as follows.

Needs["RootSearch`"]
a = 0.03; b = d = 1; d2 = 0; c = 0.5;
System1 = 
   ParametricNDSolveValue[{X'[t] == 
         Piecewise[{{X[t]bH[t]/HMax - X[t]d(X[t] + Y[t])/K1 - X[t]d2,
             t >= t0}, {X[t]bH[t]/HMax - X[t]d(X[t])/K1 - X[t]d2, 
            t < t0}}], 
        H'[t] == 
         Piecewise[{{Y[t]b(1 - cg)(1 - H[t]/HMax)g - aH[t], t >= t0},
          {-aH[t], t < t0}}], 
        Y'[t] == 
         Piecewise[{{Y[t]b(1 - cg)H[t]/HMax + 
             Y[t]b(1 - cg)g(1 - H[t]/HMax) - Y[t]d(X[t] + Y[t])/K1,
             t >= t0}}], X[0] == K1, H[0] == HMax, Y[0] == 1}, {X[t], 
        H[t], Y[t]}, {t, 0,5000}, {g, K1, HMax, t0}];
FinalTime1[g_, t0_, K1_, HMax_] := Module[{YTrajectory, Times1},
   YTrajectory = (System1[g, K1, HMax, t0][[3]]) /. t -> t + t0;
   Times1 = RootSearch[YTrajectory == 100, {t, 0.00000001, 500}];
   Max[Flatten[Times1][[All, 2]]]
]
Function1[g_, t0_, K1_, HMax_] := 
  NIntegrate[
    System1[g, K1, HMax, t0][[3]], {t, t0, t0 + FinalTime1[g, t0, K1, HMax]}
  ];
NIntegrate[Function1[1, t0, 2000, 2000], {t0, 0, 100}]

When I try this, I get a RootSearch error that says "False is not a well-formed equation." Is there a way around this problem?

Answer 1

RootSearch is a Mathematica package I wrote a long time ago and it can be downloaded here. To install the definitions with a recent version of Wolfram Language, I had to open (RootSearch.m) with WolframLanguage and click the button in the top right to (Run all code). To get this working, you should modify your code to use this:

RootSearch[YTrajectory[t] == 100, {t, 0.00000001, 500}];

where YTrajectory[t] returns a number when t is a real number. You may need to define YTrajectory[t] as follows

YTrajectory[t_?NumericQ]:= (* your code *)

Your code was too confusing for me to give more specific advice.