Nonlinear problems

We want to solve a nonlinear PDE.

A simple scalar PDE

We consider the simple PDE $$ - \Delta u + 5 u^2 = 1 \text{ in } \Omega $$ on the unit square $\Omega = (0,1)^2$. We note that this PDE can also be formulated as a nonlinear minimization problem (cf. unit 3.6 ).

In [1]:
import netgen.gui
%gui tk
import tkinter
from netgen.geom2d import unit_square
from ngsolve import *

mesh = Mesh (unit_square.GenerateMesh(maxh=0.3))

In NGSolve we can solve the PDE conveniently using the linearization feature of SymbolicBFI.

The "BilinearForm" (which is not bilinear!) needed in the weak formulation is $$ A(u,v) = \int_{\Omega} \nabla u \nabla v + 5 u^2 v - 1 v ~ dx \quad ( = 0 ~ \forall~v \in H^1_0) $$

In [2]:
V = H1(mesh, order=3, dirichlet=[1,2,3,4])
u = V.TrialFunction()
v = V.TestFunction()

a = BilinearForm(V)
a += SymbolicBFI( grad(u) * grad(v) + 5*u*u*v- 1 * v)

We use Newton's method and make the loop:

  • Given an initial guess $u^0$
  • loop over $i=0,..$ until convergence:
    • Compute linearization: $A(u^i) + \delta A(u^i) \Delta u^{i} = 0 $:
      • $f^i = A(u^i)$
      • $B^i = \delta A(u^i)$
      • Solve $B^i \Delta u^i = -f^i$
    • Update $u^{i+1} = u^i + \Delta u^{i}$
    • Evaluate stopping criteria
  • Evaluate $E(u^{i+1})$

As a stopping criteria we take $ \langle A (u^i),\Delta u^i \rangle = \langle A (u^i), A (u^i) \rangle_{(B^i)^{-1}}< \varepsilon$.

In [3]:
u = GridFunction(V)
Draw(u,mesh,"u")

res = u.vec.CreateVector()
du = u.vec.CreateVector()

for it in range(25):
    print ("Iteration",it)
    a.Apply(u.vec, res)
    a.AssembleLinearization(u.vec)
    #a.Assemble()
    
    du.data = a.mat.Inverse(V.FreeDofs()) * res
    u.vec.data -= du
    
    #stopping criteria
    stopcritval = sqrt(abs(InnerProduct(du,res)))
    print ("<A u",it,", A u",it,">_{-1}^0.5 = ", 
           stopcritval)
    if stopcritval < 1e-13:
        break
    
Iteration 0
<A u 0 , A u 0 >_{-1}^0.5 =  0.18745205801776055
Iteration 1
<A u 1 , A u 1 >_{-1}^0.5 =  0.0025712160474074893
Iteration 2
<A u 2 , A u 2 >_{-1}^0.5 =  5.324479829649632e-07
Iteration 3
<A u 3 , A u 3 >_{-1}^0.5 =  2.3040755730987054e-14