# Taking a derivative, then squaring a two dimensional numerical function

Context

I am interested in first finding an interpolating function of the solution to the linearly damped wave equation. Here is the solution to the LDWE with smooth square inital function :

ub[epsilon_] :=
Interpolation[{{{-5.944 - epsilon}, 0, 0}, {{-5.944 + epsilon}, 1,0},
{{5.944 - epsilon},     1, 0}, {{5.944 + epsilon}, 0, 0}},
"ExtrapolationHandler" -> {0 &, "WarningMessage" -> False}]

square = ub[0.0001]; initialWavePulseSquare[x_] := square[x];

exactSolSquare[x_, t_, c_, τ_] := {pEvaluate = Sqrt[1 - (x - s)^2/(c t)^2];
limitInt1 = x - c t; limitInt2 = x + c t;
1/2 (initialWavePulseSquare[x - c t] + initialWavePulseSquare[x + c t])*Exp[-t/(2 τ)]
+ 1/(4 c τ) Exp[-t/(2 τ)]*
Quiet[NIntegrate[
initialWavePulseSquare[s]*(BesselI[0, pEvaluate*t/(2 τ)] +
(1/pEvaluate)* BesselI[1, pEvaluate*t/(2 τ)] ) , {s, limitInt1,
limitInt2}]]}

solutionRewritten[x_, t_] := exactSolSquare[x, t, 1, 25]


Here, I used

square = ub[0.0001];
initialWavePulseSquare[x_] := square[x]


to define a "smooth" rectangular initial wave pulse. Then, I defined the solution to the LDWE as exactSolSquare[x, t, c, τ] where c and tau are speed and damping constant, respectively. Then, I was able to create a numerical two dimensional solution to the LDWE defined as solutionRewritten[x_, t_]

Objective

I want to ultimately take the integral of (Derivative[0,1][solutionRewritten][x,t])^2 from t=0 to t=300. However, I have hard time even getting the InterpolatingFunction.

Question

How can I create an interpolating function of two dimensional numerical data and evaluate the integral of (Derivative[0,1][solutionRewritten][x,t])^2 from t=0 to t=300.

Attempt 1

I tried something like this and MM threw a handful of errors:

Interpolation[Table[Evaluate[{{x, t}, solutionRewritten[x, t]}], {x, -200, 200, .5},
{t, 0.01, 300, 0.5}] ]


I previously asked similar question before. But this time, I know the solution of NDSolve explicitly. Here is the link: Integrating Squared of Interpolating Function with respect to one variable

Attempt 2

I tried using a smaller boundary to see if the code even runs. And as expected, it does run and gives me an InterpolatingFunction.

solnInterpolation = ListInterpolation[Table[First[solutionRewritten[x, t]], {x, -20, 20}, {t, 0, 20}]]


However, it outputs the following error: Power::infy: "Infinite expression 1/0 encountered". I tried to suppress the error messages and re-evaluated the code with larger boundary as follows:

Off[Power::infy]; solnInterpolation = ListInterpolation[Table[First[solutionRewritten[x, t]], {x, -200, 200}, {t, 0, 300}]]


However, MM evaluates this line of code runs indefinitely and doesn't output an InterpolatingFunction. Is there a way to speed this process with pretty good precision?

Thank You very much for your help!