We have a thin crystal of thickness d illuminated uniformly from the left at an intensity lgti[0, t] of unity. The crystal is composed of a photo reactive species a which absorbs light. Inside the crystal at location x and time t the photo chemical reaction leads to a local concentration of a given as a[x, t] in the code below.

Code to generate the concentration, a[x, t], and light intensity, lgti[x, t], within the crystal is straightforward and is shown.

The cause of the errors generated using the code is that in some cases a[x, t] does not return a number.These exceptional cases are near a known limit of the InverseFunction and a[x, t] in these cases could be given the values unity if found. My problem is that using NumberQ inside a Module definition of a[x, t] always gives false because ligti[x, t] is unevaluated.

The starting code is shown but see the attached notebook.

ClearAll[a, x, t, eqns, \[Sigma]N, lgti, soln, t0, d, tmax]

\[Sigma]N = 1.78*10^4; d = 0.001; tmax = 3000.0;

a[x_, t_] := InverseFunction[(1.8996253051848473`/lgti[x, t] *
      (162.99559471365637` Log[1 + 8.98898678414097` (1 - #1)] - 
        172.98458149779734` Log[#1]) ) &
                                                            ] [t]

eqns = {D[lgti[x, t], x ] == - \[Sigma]N a[x, t ] lgti[x, t],(*    Beer's Law *)
                   lgti[0, t] == 1, 
   lgti[x, 0] == Exp[-\[Sigma]N x]
                    };

t0 = AbsoluteTime[];

soln = NDSolve[eqns, lgti, {x, 0, d}, {t, 0, tmax}, 
              MaxStepFraction -> 0.01]   [[1]];

Print[ToString[(AbsoluteTime[] - t0)/60] <> "   minutes"]

Any advice on how to code a[x, t] so that lgti[x, t] appears as a number within the body of the code would be welcome.

An alternate approach would also be well received.

Attachments: