Numerical issues with logarithm

I need to solve a Poisson-equation that’s coupled with nonlinear equations containing logarithms.
That code is unfortunately really unstable because I keep getting negative values during iterations that make the logarithm diverge.

I got it working by cutting the logarithmic function off using a conditional statement and Taylor expansion at that point – however, this solution seems not to work when i extend this problem to pure Neumann conditions. Is there any way to restrict the value range of a function to positive values? Or are there other possibilities to avoid this problem?

Thanks!

After several workarounds, I finally found a solution that works nicely!
The general issue was that, in my variational form, there are terms of the form

\ln\left(\frac{u}{1-u}\right).

The Newton solver sometimes produced negative values of u as a numerical artifact, which led to the failure of the solver.

I got rid of diverging solutions by forcing the argument of the logarithm to remain positive (DOLFIN_EPS) with a conditional statement for u>0:

ln ( conditional( gt(u, 0), u/(1-u), DOLFIN_EPS ) )

I think that this is rather straightforward, but I didn’t think of it. I had the idea when I saw this post:

So thank you, @kamensky!

1 Like