2
$\begingroup$

I am trying to solve an energy functional equation using Mathematica's NDSolve function, but the kernel is temperamental. It almost always crashes after a second of trying to execute the code, but sometimes it will complete and return an interpolationfunction, which is what I want.

The minimal working example is:

op = Det[Inactive[Grad][Inactive[Grad][u[x, y], {x, y}], {x, y}]]

Ω = Disk[{0, 0}, 4];

usol = 
  NDSolve[{op == 0, DirichletCondition[u[x, y] == y, True]}, u, {x, y} ∈ Ω, 
  Method -> {"TimeIntegration" -> "ExplicitEuler"}]

The differential operator is the determinant of the hessian for U(x,y):

$\big(\partial^2_x u(x,y) \big)\big( \partial^2_y u(x,y)\big) - \big(\partial_x \partial_y u(x,y)\big)^2 = f(x,y)$

Its a purely spatial DiffEq with no affine parameterization. I want the code to be able to solve the equation for almost any (smooth) source function $f(x,y)$ in the NDSolve expression argument.

When I change $f(x,y)=0$ to any other number, say $f(x,y)=1$, the kernel crashes.

Is this due to the method argument in NDSolve? I tried using "ExplicitEuler" as well as "Extrapolation", but to no avail. Perhaps the PDE is too complex for Mathematica to handle?

I would expect at least some errors if it wouldn't be able to solve it. I appreciate any help!

Improve this question
$\endgroup$
3
  • 4
    $\begingroup$ What equation are you trying to solve? ExplicitEuler is a time integration method, but the setup looks like a stationary (2D) problem. $\endgroup$
    – user21
    Nov 13, 2020 at 15:53
  • $\begingroup$ Do not use the tag 'bugs'. It is reserved for community use. It will be applied to your question when and only when it has be verified that that the issue you raise is a Mahthematica bug and not your bug. $\endgroup$
    – m_goldberg
    Nov 13, 2020 at 15:57
  • 1
    $\begingroup$ @user21, I didn't realize 'time integration' literally meant solving a time-parameterized equation. I thought it was some backend terminology. The exact equation is (d_x^2 f(x,y))*(d_y^2 f(x,y)) - (d_x d_y f(x,y))^2 = g(x,y) subject to some dirichlet boundary conditions. In the sample code, g(x,y)=0. I checked the docs and I tried using PDEDiscretization but still suffering from the same described errors. Is there a better NDsolve method? $\endgroup$
    – shanedrum
    Nov 14, 2020 at 12:24

1 Answer 1

Reset to default
3
$\begingroup$

The crash is of course a bug and is fixed in the upcoming version 12.2. There you will get an error message:

enter image description here

Unfortunately, currently this can not be solved with easily with the finite element method.

Det[Grad[Grad[u[x, y], {x, y}], {x, y}]] // InputForm

-Derivative[1, 1][u][x, y]^2 + 
 Derivative[0, 2][u][x, y]*Derivative[2, 0][u][x, y]

This has a second order spatial derivative as a coefficient. Currently the FEM can have a maximum of first order derivative of the dependent variable as a coefficient. You may be able to work around this by using something like

-Derivative[1, 1][u][x, y]^2 + 
     gg[x, y]*Derivative[2, 0][u][x, y]

And put that in a look and replace gg with the second derivative of the previous solution. But that has it's own set of problems.

Improve this answer
$\endgroup$
1
  • $\begingroup$ Yes, I see. I read more into the documentation for the FEM, and I realized that was what Mathematica was angry about. I realized after I wrote out the code using functions within NDSolve FEM package, specifically InitializePDECoefficients. I guess I will have to try and construct an analytic solution myself. Good thing the PDE has intuitive meaning. I appreciate your help! $\endgroup$
    – shanedrum
    Nov 16, 2020 at 12:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.