Get system of four equations solved

Is it possible to solve a system of four (some nonlinear) PDE for four independent functions?
Something like this

    from dolfin import *
    omega = UnitCubeMesh(8, 8, 8)
    omega_b = BoundaryMesh(omega, 'exterior')
    finite_element = FiniteElement("Lagrange", omega.ufl_cell(), 2)

    state_space = FunctionSpace(omega, finite_element * finite_element *  finite_element * finite_element)
    v, h, tau, nu = TestFunctions(state_space)
    state = Function(state_space)
    a, b, c, d = split(state)
    a_eq = inner(grad(a), grad(v)) * dx + inner(a ** 4 - b, v) * dx
    b_eq = inner(grad(b), grad(h)) * dx + inner(b - a **4, h) * dx
    c_eq = inner(grad(c), grad(tau)) * dx + inner(c, a ** 3 * tau) * dx 
    d_eq = inner(grad(d), grad(nu)) * dx + inner(d, b ** 3 * nu) * dx 
    a_src = 0.3 * v * ds
    b_src = c * h * ds
    solve(
        a_eq + b_eq - a_src - b_src + c_eq + d_eq == 0,
        state,
        form_compiler_parameters={"optimize": True, 'quadrature_degree': 3}
    )

This code generates exception on line
v, h, tau, nu = TestFunctions(state_space)

z1, z2, x1, x2 = TestFunctions(state_space)
ValueError: not enough values to unpack (expected 4, got 2)

Thanks in andvance.

I think you need something like this

V = TestFunctions(state_space)
v, h, tau, nu = split(V)

I hope this helps.

v, h, tau, nu = split(V)
  File "/usr/local/lib/python3.6/dist-packages/ufl/split_functions.py", line 64, in split
    element = v.ufl_element()
AttributeError: 'tuple' object has no attribute 'ufl_element'

Use a MixedElement([finite_element]*4). The tensor product notation you use for the elements probably isn’t what you intend to do.

from dolfin import *
omega = UnitCubeMesh(8, 8, 8)
omega_b = BoundaryMesh(omega, 'exterior')
finite_element = FiniteElement("Lagrange", omega.ufl_cell(), 2)

state_space = FunctionSpace(omega, MixedElement([finite_element]*4))
v, h, tau, nu = TestFunctions(state_space)
1 Like

Exactly what I wanted! But why the formulation for system of two equations is

finite_element = FiniteElement("Lagrange", omega.ufl_cell(), 2)
state_space = FunctionSpace(omega, finite_element * finite_element)

that presented in doc samples is working as intended?

Hi, see the Fenics tutorial book, page 76-7 ([https://fenicsproject.org/tutorial/]). I think it will provide you with the explanation you are after.

It boils down to that writing ‘element=P1 * P1 * P1’ is equivalent to writing ‘element = (P1 * P1) * P1’ so the result will be a mixed element consisting of two subsystems, the first of which in turn consists of two scalar subsystems.

1 Like