Three differential equations

Question

I have the next figure with inflow and outflow rates:

And I need to calculate the time for which the pollution will be 50% from the initial value. Assuming that all lakes have the same pollution concentration p initially. The pollution flows from a lake to another in the chain and only inflows not from a lake are clean water.

I wrote the equations for first 3 lakes:

s'[t] + 15/2900 * s[t] == 0

m'[t] + 38/1180 * m[t] == 0

h'[t] + 68/850 * h[t] - 15/2900 * s[t] - 38/1180 * m[t] == 0

And because all the lakes initially have the same pollution concentration p => s(0)==p, m(0)==p, h(0)==p, right?

First equation:

eqS = s'[t] == -15/2900*s[t]
solS = DSolve[{eqS, s[0] == p}, s, t]
s[t_] = s[t] /. First@solS

Second equation:

eqM = m'[t] == -38/1180*m[t]
solM = DSolve[{eqM, m[0] == p}, m, t]
m[t_] = m[t] /. First@solM

Third equation:

eqH = h'[t] == 15/2900*s[t] + 38/1180*m[t] - 68/850*h[t]
solH = DSolve[{eqH, h[0] == p}, h, t]
h[t_] = h[t] /. First@solH

And to find the time for the first and second equation I used:

Solve[s[t] == 0.5*p]
Solve[m[t] == 0.5*p]

And it worked ok, but trying to do the same thing for h[t] takes too much time to compute, and I wonder if my solution is ok.

Answer 1

The problem is that the function h is too complicated for solving it symbolically with Solve. You can use FindRoot to obtain a numerical solution as follows:

eq = Simplify[h[t]/p == 1/2]
FindRoot[eq, {t, 0}]

{t -> 18.0728}

By the way: You can also solve the system of differential equations at once; this might get handy for more complicated task, e.g. if there were cyclic dependencies between the seas:

ClearAll[s, m, h]
des = {
   s'[t] + 15/2900*s[t] == 0,
   m'[t] + 38/1180*m[t] == 0,
   h'[t] + 68/850*h[t] - 15/2900*s[t] - 38/1180*m[t] == 0
   };
initials = {s[0] == p, m[0] == p, h[0] == p};
{s, m, h} = DSolveValue[Join[des, initials], {s, m, h}, t];
h

Function[{t}, ( E^(-2 t/25) (7867 + 20615 E^(141 t/2950) + 2115 E^(217 t/2900)) p)/30597]