Weak formulation for square of gradient

variational formulation
2

Welcome robert. First: You can use LaTeX for formulas using $ signs.

Concerning your specific problem:

-\operatorname{div}\operatorname{grad}u+\operatorname{grad} u\cdot \operatorname{grad} u=f

The weak form is as usual obtained through multiplying by test function v and integration (assuming homogenous Neumann boundary conditions)

\int_\Omega \operatorname{grad}v\cdot\operatorname{grad}u\,\text dx +\int_\Omega v(\operatorname{grad} u)^2 \,\text dx =\int_\Omega fv\,\text dx\,.

You are on the right track in your approach, but the syntax solve(a==L, ...) is used for linear problems and the one you want to solve is nonlinear. But we can simply define a nonlinear variational form as

F(u,v) = \int_\Omega \operatorname{grad}v\cdot\operatorname{grad}u\,\text dx +\int_\Omega v(\operatorname{grad} u)^2 \,\text dx -\int_\Omega fv\,\text dx

such that your solution has to satisfy F(u,v)=0. FEniCS will automatically compute the Jacobian and iteratively solve your problem when solve is called in the right way (cf. my answer here: Function() v.s. TrialFunction())

Consider the following code:

from fenics import *

# Create mesh and define function space
mesh = UnitSquareMesh(8, 8) 
V = FunctionSpace(mesh, 'P', 1)

# Define boundary condition
u_D = Expression('1 + x[0]*x[0] + 2*x[1]*x[1]', degree=2) 
def boundary(x, on_boundary):
   return on_boundary 
bc = DirichletBC(V, u_D, boundary)

# Define variational problem

u = Function(V) 
v = TestFunction(V) 
f = Constant(-6.0) 
F = dot(grad(u), grad(v)) * dx + dot(grad(u), grad(u)) * v * dx - f*v*dx

solve(F==0, u, bcs = bc)

For introductory reading into nonlinear problems in FEniCS I recommend
http://home.simula.no/~hpl/homepage/fenics-tutorial/release-1.0-nonabla/webm/nonlinear.html
as a first read.

3 Likes