Message Boards Message Boards

Suggest NDSolve Method options when argument contains InverseFunction?

Mervin Hanson
Posted 9 years ago

I have a Beer's Law problem where the absorbing species is given as an InverseFunction. More details are given in the attached notebook. If you run the code or look at the attached notebook, you will see NDSolve produces and error. I have tried a few Method options to no avail.

Any suggestions?

The problem code is

ClearAll[ a, x, t, tmin, tmax, \[Sigma],  nN, \[Sigma]N, lgtI, lgtIo,
  lgtIx0t0, soln, tmp, time0]
tmin = 0.0; d = 1.0 10^-3; tmax = 2400.0;
\[Sigma] =  1.78 10^-17; nN =  
 1.0 10^21; \[Sigma]N = \[Sigma] nN; lgtIx0t0 = 3.1*10^15; 
lgtIo[x_] := 
 lgtIx0t0 Exp[-\[Sigma]N x] (*we are sure of values only at time 0*)(* \
but for non zero times *)
a[x_, t_] := InverseFunction[
   ( -(1/lgtI[x, t])
        6.610388104920668`*^14 (-1391.97544007514` Log[
          9.732394366197182` - 8.732394366197182` #1] + 
        1476.962545808129` Log[#1])) &] [t]      
eqns = {D[lgtI[x, t], x ] == - \[Sigma]N a[x, t] lgtI[x, t],
                    lgtI[0, t] == lgtIx0t0, lgtI[x, 0] == lgtIo[x]
                    };
time0 = AbsoluteTime[];
soln =
 NDSolve[eqns, lgtI, {x, 0.0, d}, {t, 0.00, tmax}]
Print[(AbsoluteTime[] - time0)/60 minutes]
Speak["Done"]
Plot3D[Evaluate[lgtI[x, t] /. soln], {x, 0, 0.000352}, {t, 0, tmax}, 
 PlotRange -> All,
 AxesLabel -> {"  x  ", "  t  ",  "     lgtI  "}, ImageSize -> Large]
Plot3D[Evaluate[a[x, t] /. soln], {x, 0, 0.000352}, {t, 0, tmax}, 
 PlotRange -> All,
 AxesLabel -> {"  x  ", "  t  ",  "  a  "}, ImageSize -> Large]
Attachments:
POSTED BY: Mervin Hanson
6 Replies
Sort By:

If don't misunderstand your code, the function you need to invert is of the form f[x, t]*g[#], so you only need to inverte the one-variable function

g[r_] = -1392 Log[691/71 - (620 r)/71] + 1477 Log[r]

If you are content with numerical methods, you can build an approximate inverse by interpolation in a suitable interval.
POSTED BY: Gianluca Gorni
Mervin Hanson
Mervin Hanson
Posted 9 years ago

I do not see how this helps. I need the inverse of what you write as f[x, t] * g[#]. Consider f == 2 and g= 1/(1 + #)& with the code

In[23]:= ClearAll[a, f, g, h]
g =  1/(1 + #) &
InverseFunction[g]
h = InverseFunction[2 g]
h[x]
h[0.5]


Out[24]= 1/(1 + #1) &

Out[25]= (1 - #1)/#1 &

Out[26]= InverseFunction[2 (1/(1 + #1) &)]

Out[27]= InverseFunction[2 (1/(1 + #1) &)][x]

InverseFunction[2 (1/(1 + #1) &)][0.5]

I need a number.
POSTED BY: Mervin Hanson

It seems to me that Mathematica cannot find a fomula for the inverse:

InverseFunction[-1392 Log[691/71 - (620 #)/71] + 1477 Log[#] &]
Reduce[-1392 Log[691/71 - (620 r)/71] + 1477 Log[r] == t, r, Reals]

However, you can plot the inverse

g[r_] = -1392 Log[691/71 - (620 r)/71] + 
  1477 Log[r];
ParametricPlot[{g[r], r}, {r, 0, 691/620}, 
     AspectRatio -> 1/GoldenRatio]

and build a numerical inverse using interpolation:

rmin = 1/620; rmax = 691/620 - rmin;
data = Table[
   {{g[r]}, r, 1/g'[r], -g''[r]/g'[r]^3}, {r, rmin, rmax, 
    rmin}]//N;
myInverse = Interpolation[data];
myInverse[100]

Plot[myInverse[t], {t, g[rmin], g[rmax]}]

I don't know if this gives enough precision for your purposes.
POSTED BY: Gianluca Gorni
Mervin Hanson
Mervin Hanson
Posted 9 years ago

It seems to me you are solving a problem I do not have. You give a method of finding the inverse of g. What I need is the inverse of the product ( 1/lgtI[x, t] ) g[t]. I think what my first response showed was that knowing the inverse of g was not enough.
POSTED BY: Mervin Hanson
Mervin Hanson
Mervin Hanson
Posted 9 years ago

This discussion is going in the wrong direction. NDSolve is asked to integrate a first order PDE for a function lgtI[x, t] involving the derivative wrt x and a function of x, t, and a[x, t]. The function a[x, t] is given as an inverse function.The domain of integration is {0 < x < d} and {0, t < tmax}. Mathematica has no problem here with the InverseFunction. NEVER HAVING A PROBLEM WITH THE INVERSEFUNCTION until x is about 0.00035 on its way to 0.001. Then an error occurs. My question is simply (?) How do I stop this error from appearing?
POSTED BY: Mervin Hanson
Mervin Hanson
Mervin Hanson
Posted 9 years ago

It is time to retire this discussion
POSTED BY: Mervin Hanson
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract