# Having trouble with Nonlinear Poisson

**URL:** <https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657>\
**Category:** Uncategorized\
**Created:** [March 16, 2020, 6:18pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657 "2020-03-16T18:18:41Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![Sterling](https://avatars.discourse-cdn.com/v4/letter/s/f475e1/32.png) [@Sterling](https://fenicsproject.discourse.group/u/Sterling)\
**Post date:** [March 16, 2020, 6:18pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657/1 "2020-03-16T18:18:41Z")

</div>

Hi, I am having trouble solving the following equation:  
 -\nabla^2u + A = \frac{B}{u^2}  
where A and B are constants.  
This is calculated in the space between two spheres, one located within a larger one.  
I set boundary conditions (zero) at surfaces, and I confirmed that the code works for regular electrostatic equation (right hand side of equation being constant instead of 1/u^2).  
I also set g=0 as another boundary condition.

```python
u = fn.Function(V)
v = fn.TestFunction(V)
f = B/u**2
g = fn.Constant(0)

a = (fn.inner(fn.grad(u),fn.grad(v))+ A*v)*fn.dx
L = f*v*fn.dx + g*v*fn.ds
F = a-L

fn.solve(F==0, u, bcs)

```

However, this returns an error from Newtonian solver that suggests the solution didn’t converge.  
The output from Newtonian solver only shows r (abs) =±nan, so something is wrong.  
I know the solution exists, since this is a well known physical model. Can someone point out what I’m doing wrong with this?

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [March 16, 2020, 6:25pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657/2 "2020-03-16T18:25:40Z")

</div>

Well, this problem is not well defined, as it is only determined up to a constant when solving Poisson with pure Neumann conditions, see: [Poisson with pure Neumann conditions](https://bitbucket.org/fenics-project/dolfin/src/ac79caa94a72906b0e56b14e38e250c86620f553/python/demo/documented/neumann-poisson/demo_neumann-poisson.py.rst?at=master)

---

<div class="post-metadata">

**Author:** ![Sterling](https://avatars.discourse-cdn.com/v4/letter/s/f475e1/32.png) [@Sterling](https://fenicsproject.discourse.group/u/Sterling)\
**Post date:** [March 16, 2020, 6:50pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657/3 "2020-03-16T18:50:08Z")

</div>

Hi dokken, Thanks for your comment.  
I have both Dirichlet and Newmann boundary.

Following is an example of working code.

```python
from dolfin import *

# Create mesh and define function space
mesh = UnitSquareMesh(32, 32)
V = FunctionSpace(mesh, "Lagrange", 1)

# Define Dirichlet boundary (x = 0 or x = 1)
def boundary(x):
    return x[0] < DOLFIN_EPS or x[0] > 1.0 - DOLFIN_EPS

# Define boundary condition
u0 = Constant(0.0)
bc = DirichletBC(V, u0, boundary)

# Define variational problem
u = Function(V)
v = TestFunction(V)
f = 1 # <-----------------want to use 1/u**2
g = Constant(0)
a = inner(grad(u), grad(v))*dx
L = f*v*dx + g*v*ds
F = a-L

# Compute solution
solve(F==0, u, bc)

# Plot solution
plot(u, interactive=True)

```

Does what I’m saying make sense or am I doing something completely wrong?

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [March 16, 2020, 7:08pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657/4 "2020-03-16T19:08:18Z")

</div>

Your zero initial guess makes the problem singular. Using for instant `u=project(Constant(1), V)` as an initial guess solves the problem

---

<div class="post-metadata">

**Author:** ![Sterling](https://avatars.discourse-cdn.com/v4/letter/s/f475e1/32.png) [@Sterling](https://fenicsproject.discourse.group/u/Sterling)\
**Post date:** [March 16, 2020, 7:51pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657/5 "2020-03-16T19:51:25Z")

</div>

Now it worked! Thank you so much for your help, @dokken !

---

<div class="post-metadata">

**Author:** ![Sterling](https://avatars.discourse-cdn.com/v4/letter/s/f475e1/32.png) [@Sterling](https://fenicsproject.discourse.group/u/Sterling)\
**Post date:** [March 16, 2020, 7:57pm UTC](https://fenicsproject.discourse.group/t/having-trouble-with-nonlinear-poisson/2657/6 "2020-03-16T19:57:41Z")

</div>

In case anyone else is interested, here is my final working example code:

```python
from dolfin import *

# Create mesh and define function space
mesh = UnitSquareMesh(32, 32)
V = FunctionSpace(mesh, "Lagrange", 1)

# Define Dirichlet boundary (x = 0 or x = 1)
def boundary(x):
    return x[0] < DOLFIN_EPS or x[0] > 1.0 - DOLFIN_EPS

# Define boundary condition
u0 = Constant(0.0)
bc = DirichletBC(V, u0, boundary)

# Define variational problem
u = project(Constant(1),V)
v = TestFunction(V)
f = 1/u**2
g = Constant(0)
a = inner(grad(u), grad(v))*dx
L = f*v*dx + g*v*ds
F = a-L

# Compute solution
solve(F==0, u, bc)

# Plot solution
plot(u, interactive=True)

```

Also just a quick note: As it always does, initial guess requires to be within the similar order of magnitude to the final solution (the solution did not converge when the initial guess was off by factor of 10 for more complicated case when I tested). This tolerance is probably model-dependent, but I haven’t tested with other methods.