# NIntegrate InterpolatingFunction Kernel Crash

I would like to Investigate if the edges of this 2D heat equation are insulated by the Mathematica assumed Neumann conditions.

aap = NDSolve[{\!$$\*SubscriptBox[\(\[PartialD]$$, $$t$$]$$f[x, y, t]$$\) == 0.05 \!$$\*SubsuperscriptBox[\(\[Del]$$, $${x, y}$$, $$2$$]$$f[x, y, t]$$\),
f[x, y, 0] == E^(-40 (x^2 + y^2))},
f, {t, 0, 1}, {x, y} \[Element] Disk[]][[1]]


The problem is that Mathematica seems unable to Integrate the solution to the problem.

    Export["test.gif",
Table[Plot3D[aap[x, y, i], {x, y} \[Element] Disk[],
PlotRange -> {Automatic, Automatic, {0, 2}}], {i, .1, 5, 0.2}],
"AnimationRepetitions" -> \[Infinity], "ImageSize" -> 20]


I try

NIntegrate[f[x, y, 0] /. aap, {x, y} \[Element] Disk[]]


But my kernel just quits, same thing on Mathematica Online.

Do you also see this?/Is someone able to help?

edit:

I now did get the result one time:

In[5]:= NIntegrate[f[x, y, 0] /. aap, {x, y} \[Element] Disk[]]

During evaluation of In[5]:= Outer::heads: Heads Alternatives and List at positions 3 and 2 are expected to be the same.

During evaluation of In[5]:= NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

During evaluation of In[5]:= NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 0.07854482887754174 and 2.6180493354023447*^-7 for the integral and error estimates.

Out[5]= 0.0785448


which is the result I need. After that I managed to crash the kernel again and I am back at the same point...

edit (2) the workaround solution as mention below works! The following now shows there is indeed no outflow of "heat".

aap = NDSolveValue[{D[f[x, y, t], t] ==
0.05 Laplacian[f[x, y, t], {x, y}],
f[x, y, 0] == 40/\[Pi] E^(-40*(x^2 + y^2))},
f, {t, 0, 10}, {x, y} \[Element] Disk[]];


and plot

Plot[NIntegrate[
aap[x, y, t], {x, y} \[Element] aap["ElementMesh"]], {t, 0, 10},
PlotRange -> {0, 2}]


to get:

aap = NDSolveValue[{D[f[x, y, t], t] ==
0.07854480753557473
`