Non-linear PDE in C++ - Writing the UFL file.

Question

Hello everyone,

This is my first time using Fenics and I've run into a bit of trouble.

I am trying to construct the UFL file to setup solving the following PDE:

$$ \nabla^2 \phi = -\frac{\Lambda^5}{\phi^2} + \frac{\rho}{M} $$

where $ \Lambda $ and $M$ are constants and $\rho$ is the source term.

I'm trying to do this by defining the semi-linear form F = (...) as for the dolfin implementation, but I'm running into the error:

"FFC: Applying nonlinear operator to expression depending on form argument v_1".

Here is my UFL:

Elements

element = FiniteElement("Lagrange", triangle, 1)

Trial and Test Functions

u = TrialFunction(element)
v = TestFunction(element)
f = Coefficient(element)

Parameters

Lambda = Constant(triangle)
M = Constant(triangle)

Semilinear Form

F = inner(grad(u), grad(v))dx - ((fv/M)dx) + (pow(Lambda, 5.0)v/pow(u, 2.0))*dx

I've tried my best to piece together what this should look like based on the documentation, but giving that this is my first project I definitely have a lot to learn.

Any advice would be hugely appreciated.

Thanks!

Answer 1

If you look at the nonlinear demos, e.g.
https://bitbucket.org/fenics-project/dolfin/src/26b682abe61b59ff952c3fe0ca6cde307e8ae16a/demo/documented/nonlinear-poisson/python/demo_nonlinear-poisson.py?at=master&fileviewer=file-view-default

you will notice that u is a Function, not a TrialFunction.

Comments

Thank you for your response. I should have mentioned that I've tried that too, upon switching TrialFunction to Function, I get

"NameError: Name 'Function' is not defined".

Any ideas?

---

Right, use Coefficient instead in C++.
Your forms must be linear in TrialFunction and TestFunction,
as explained in the documentation on UFL.
Function(Python) and Coefficient(C++) must be used
to express nonlinearities. Look at the C++ version of
nonlinear Poisson.

---

Hi,

Thanks, the use of Coefficient allows me to compile everything correctly. Sadly, I'm running into the following error upon execution:

" Coefficient number 1 ("Lambda") has not been set. "

Here is the segment of my main cpp file, where as far I'm concerned I've defined both of the constants within the necessary form:

      // Create mesh and function space

      chameleonSetup::FunctionSpace V(cMesh);

      // Define boundary condition
      Constant u0(0.0);
      DirichletBoundary boundary;
      DirichletBC bc(V, u0, boundary);

     // Create residual form defining the problem

      chameleonSetup::ResidualForm F(V);

      Source f;
      Function u(V);    
      Constant Lambda(10.0);
      Constant M(1.0);

      F.Lambda = Lambda;
      F.M = M;
      F.u = u;
      F.f = f;


      // Jacobian Form

      chameleonSetup::JacobianForm J(V, V);
      J.u = u;

      // Solver Parameters
      Parameters params("nonlinear_variational_solver");
      Parameters newton_params("newton_solver");
      newton_params.add("relative_tolerance", 1e-6);
      params.add(newton_params);

      // Solve Variational Problem

      solve(F == 0, u, bc, J, params);

Do I need to set the value of Lambda elsewhere?

Thanks

---

You need to assign all the Coefficients and Constnats to J also (ie, if they survive a derivative).
So I think you need:

J.Lambda = Lambda;