I am dealing with a 3D nonlinear PDE which looks like, with u - scalar solution, f - vector solution:
div(f) = 0;
where norm2(f) is the L2 norm. and a b c are predefined constant.
As you can see dividing f by its norm can cause dividing-by-zero, as a can be nonzero. In simple numpy scripts I can just do a if statement to make such case output a 0 vector. i.e. I tried
sele = FiniteElement(“CG”, mesh.ufl_cell(), 1)
vele = VectorElement(“CG”, mesh.ufl_cell(), 1)
mixed_space = FunctionSpace(mesh, MixedElement([sele, vele]))
U = Function(mixed_space)
u,f = split(U)
V = TestFunction(mixed_space)
v,tau = split(V)
However this woudn’t compile with error
UFLException: Shape mismatch between conditional branches.
How should I get around with this or does fenics have its all function to deal with it?
Thanks a lot in advance!