Code Translation from v 0.7.3 to 0.8.0 : Computation issues

dolfinx
1

Translate from v 0.7.3 to v 0.8.0 a mixed space code

I have a porous media code which works in dolfinx v 0.7.3 and computes in approximately 4 minutes.

Based on previous posts, I understood that the changes were in the element / mixed_element now from basix.ufl and FunctionSpace replaced by functionspace.

Therefore I changed:

# Updated mesh space
updated_mesh_space    = FunctionSpace(mesh, mesh.ufl_domain().ufl_coordinate_element())
# Parameter space
DG0_space = FunctionSpace(mesh, ("DG", 0))
# Mixed Space (R2,R) -> (u,p)
CG1_v     = VectorElement("CG", mesh.ufl_cell(), 1)
CG2       = VectorElement("CG", mesh.ufl_cell(), 2)
CG1       = FiniteElement("CG", mesh.ufl_cell(), 1)
# 
CG1v_space        = FunctionSpace(mesh, CG1_v)
CG2_space 	      = FunctionSpace(mesh, CG2)
CG1_space 		  = FunctionSpace(mesh, CG1)
MS                = FunctionSpace(mesh, MixedElement([CG2,CG1,CG1]))
# 
tensor_elem  = TensorElement("CG", mesh.ufl_cell(), degree=1, shape=(3,3))
tensor_space = FunctionSpace(mesh, tensor_elem)

to

# Updated mesh space
updated_mesh_space    = functionspace(mesh, mesh.ufl_domain().ufl_coordinate_element())
# Parameter space
DG0_space = functionspace(mesh, ("DG", 0))
# Mixed Space (R2,R) -> (u,p)
P1_v     = element("P", mesh.topology.cell_name(), degree=1, shape=(mesh.topology.dim,))
P2       = element("P", mesh.topology.cell_name(), degree=2, shape=(mesh.topology.dim,))
P1       = element("P", mesh.topology.cell_name(), degree=1)

# 
P1v_space        = functionspace(mesh, P1_v)
P2_space 	      = functionspace(mesh, P2)
P1_space 		  = functionspace(mesh, P1)
MS                = functionspace(mesh=mesh, element=mixed_element([P2,P1,P1]))
# 
tensor_elem  = element("P", mesh.topology.cell_name(), degree=1, shape=(3,3))
tensor_space = functionspace(mesh, tensor_elem)

This changed performed, the solver is not converging anymore and leads to ‘nan’ values. Please note that the code does work in 0.7.3 with the changes.

The codes / log_outputs and the mesh can be downloaded here.

Thank you for the help.

2

Is there any way you can simplify your problem, for instance with a coarser mesh?

From your log of 0.8.0, every Newton iteration takes 10 minutes, and it’s unlikely anyone will have the patience to run the code more than once when trying to help you debug what is going on.

3

To be fair, a Newton iteration of 0.7.3 takes less than a minute, which to me suggests that you may want to look into this commit Let PETSc choose default LU solver within NewtonSolver (#2803) · FEniCS/dolfinx@6d0720e · GitHub, especially to the change in the demo.

1 Like
4

Dear @francesco-ballarin,

I am unsure to have understood your point.

As you noticed it takes more than 10’ per iteration on the 0.8.0 returning Nan every time. That is my concern. The 0.8.0 code applied in the 0.7.3 as well as the 0.7.3 in the 0.7.3 version are converging without error in less than a minute per iteration. If it passes step one it should work.

You are suggesting that the issue would be the solver. But in both 0.7.3 and 0.8.0 I stick to the default value (around line 400 in the py codes). The demo you are talking about is for demo_cahn-hilliard and I do not understand why it should be introduced here.

Thank you for your time and help.

5

Without a reproducible example it is quite hard to give advice.

If you are using default values, the depend on two things:

  • the version of petsc
  • any hardcoded defaults, such as

Haven’t changed. As you can see here, the newton solver previously had mumps as default, while now it is left to petsc to use its default LU (which is the built in petsc LU).

You can set this to be the mumps by calling something along the lines of

ksp = solver.krylov_solver
opts = PETSc.Options()
option_prefix = ksp.getOptionsPrefix()
opts[f"{option_prefix}ksp_type"] = "preonly"
opts[f"{option_prefix}pc_type"] = "lu"
opts[f"{option_prefix}pc_factor_mat_solver_type"] = "mumps"
ksp.setFromOptions()
7

The changelog in the tutorial tried to list everything, but some things slip under my radar.

1 Like
8

Dear @dokken,

I fully understand your point, so that’s why I wanted to indicate you that it was not present yet if you have time to complete it.

Your tutorials / changelog / help in the discourse are huge ressources to work with, so thank you.