# How can I obtain the domain of a resulting Interpolating Function?

I have a problem with obtaining the domain of resulting Interpolating Function in my calculations. I want to solve the following differential equation

V[u_, v_] =ProductLog[E^(1/2 (-2 - u + v))]/(1 + ProductLog[E^(1/2 (-2 - u + v))])^4;
sol = First[NDSolve[{-4*D[S[u, v], u, v] == V[u, v]* S[u, v],
S[u, 0] == Exp[-(100/18)], S[0, v] == Exp[-((v - 10)^2/18)]},S, {u,100,200}, {v,100,200},
Method -> {"MethodOfLines","SpatialDiscretization" -> { "TensorProductGrid",
"MaxStepSize" -> 1}}, AccuracyGoal -> 1]]


As you can see, I determined the range of u and v from 100 to 200, but if you run the code, you will find that the resulting Interpolating Function range will be from 0 to 200 for v. How can I fix this problem? Thanks all.

• It seems to be the presence of the statement S[u, 0] == Exp[-(100/18)] in your NDSolve command that places the lower bound on v. If you change that to S[u, 10] == ..., then the interpolating function changes accordingly. I can't speak to whether this is the intended behavior - it seems to treat the bounds on u and v differently. Commented Aug 13, 2018 at 20:48
• @JasonB. Actually, S[u, v] == Exp[-(v - 10)^2/18] is the boundary condition of the problem which I used S[u, 0] and S[0, v] in calculations. Can I use S[u,100] == Exp[-(90)^2/18] as a new boundary condition to obtain the result that I am looking for? Commented Aug 13, 2018 at 21:16
• Can't you just try the new BCs and see what happens? (I'm pretty sure that the reason for the v domain is that in the method of lines, the spatial variable's derivatives are approximated on a spatial grid between the boundaries. S[u, 0] specifies one boundary to be v == 0. Now why can't NDSolve discard all the points outside the specified interval? I don't know. I suppose one cannot have an arbitrary boundary coincide with a grid point. Maybe that's why.) Commented Aug 13, 2018 at 23:01
• @Mehrab Because the boundary conditions are at u == 0 and v == 0, NDSolve must integrate to those surfaces, whether or not it includes them in the InterpolatingFunction. That being the case, simply ignore the unwanted part of the InterpolatingFunction. By the way, "TemporalVariable" -> v yields a domain of {{0., 200.}, {100., 200.}}, and "TemporalVariable" -> u yields a domain of {{100., 200.}, {0., 200.}}. But, why worry about this idiosyncrasy? It does not seem to cause any harm. Moreover, it is easy enough to trim the domain by postprocessing, if desired. Commented Aug 14, 2018 at 3:43
• @Mehrab I should add that AccuracyGoal -> 1 does not yield a very accurate answer. Just omit this option. Commented Aug 14, 2018 at 3:57

You can consider solving the problem in two steps, and then compare for different options AccuracyGoal. Below is the code and comparison for the four options AccuracyGoal. In Fig. it is clear that the oscillations of the solution are a consequence of the numerical error. These oscillations disappear as the accuracy increases.

V[u_, v_] =
ProductLog[
E^(1/2 (-2 - u + v))]/(1 + ProductLog[E^(1/2 (-2 - u + v))])^4;
sol = First[
NDSolve[{-4*D[S[u, v], u, v] == V[u, v]*S[u, v],
S[u, 0] == Exp[-(100/18)], S[0, v] == Exp[-((v - 10)^2/18)]},
S, {u, 0, 200}, {v, 0, 200},
Method -> {"MethodOfLines",
"SpatialDiscretization" -> {"TensorProductGrid",
"MaxStepSize" -> 1}}, AccuracyGoal -> 15]]
S0[u_, v_] := S[u, v] /. sol
sol1 = First[
NDSolve[{-4*D[S[u, v], u, v] == V[u, v]*S[u, v],
S[u, 100] == S0[u, 100], S[100, v] == S0[100, v]},
S, {u, 100, 200}, {v, 100, 200},
Method -> {"MethodOfLines",
"SpatialDiscretization" -> {"TensorProductGrid",
"MaxStepSize" -> 1}}, AccuracyGoal -> 15]]
{Plot3D[Evaluate[S[u, v] /. sol1], {u, 100, 200}, {v, 100, 200},
Mesh -> None, PlotRange -> All, PlotLabel -> ""],
Plot3D[Evaluate[S[u, v] /. sol], {u, 100, 200}, {v, 100, 200},
Mesh -> None, PlotRange -> All]}


• Thanks a lot. Although you answered the question, unfortunately, my problem is still there. For more details please have a look at my comments to Godfrey given above. Commented Aug 14, 2018 at 14:50
• Mehrab, you can not break the 2D region into intervals, you can break the region into subdomains. Thus, there is no benefit of partitioning into intervals, since it is required to specify a function on the boundaries of the domain, and not on the boundaries of the interval. Commented Aug 14, 2018 at 15:08
• Ok, thanks. I guess I should look for a supercomputer for my calculation and there is no other choice. Commented Aug 14, 2018 at 15:31
• A supercomputer is not enough to solve a problem. We need a good algorithm. With a good algorithm, we do not need a supercomputer. I'll see what else can be done in this task. Commented Aug 14, 2018 at 17:17
• Yes, you are right. Thanks for your kind effort. How can I become an expert in Mathematica? Do I need to pass some courses? Or I just need to read a book. Do you have any suggestion? Commented Aug 14, 2018 at 18:02