Any strategies for speading up compilation time?

variational formulation
1

I’ve been coding some gimmicks to parameterize my solvers and to avoid error due to missing settings (I’m using dolfin). To this aim, I’m using common python features, such as classes and dictionaries. I notice an increase in compilation time. Is it expected behavior? Would you have any tips and tricks to optimize code to shrink compilation time?

Some examples: I’m saving the traction vectors inside a dictionary; the stress tensors given by a hyperelastic constitutive model are stored in classes.

Thank you in advance

2

It is extremely difficult to give advice without a minimal working example.

The typical complaint regarding “long compilation times” arises from incorrectly compiling finite element formulations dolfinx.fem.form() inside for-loops.

3

There follows a minimal working example. I would like to check with you if this approach is suitable. Someone might advise to create a field of material properties; I made this way because my application requires multiple constitutive models in different physical groups.


from dolfin import *
import ufl_legacy as ufl

L, H, W = 1.0, 0.2, 0.3
mesh = BoxMesh(Point(0,0,0), Point(W,H,L), 10,10,10)

lower_facet = CompiledSubDomain("near(x[1], 0)")
left_facet = CompiledSubDomain("near(x[2], 0)")
right_facet = CompiledSubDomain("near(x[2], L)", L=L)
volume_1 = CompiledSubDomain("x[2]<=L*0.5", L=L)
volume_2 = CompiledSubDomain("x[2]>L*0.5", L=L)

boundary_markers = MeshFunction("size_t", mesh, mesh.topology().dim()-1)
boundary_markers.set_all(0)
lower_facet.mark(boundary_markers, 2)
left_facet.mark(boundary_markers, 3)
right_facet.mark(boundary_markers, 6)

volume_markers = MeshFunction("size_t", mesh, mesh.topology().dim())
volume_markers.set_all(0)
volume_1.mark(volume_markers, 7)
volume_2.mark(volume_markers, 8)

dx = Measure("dx", domain=mesh, subdomain_data=volume_markers, metadata={
"quadrature_degree": 2})
ds = Measure("ds", domain=mesh, subdomain_data=boundary_markers)

# Constitutive model

class Neo_Hookean:

    def __init__(self, E, v):

        self.E = E
        self.v = v
        self.mu = Constant(self.E/(2*(1+self.v)))
        self.lmbda = Constant(self.v*self.E/((1+self.v)*(1-2*self.v)))

    def strain_energy(self, C):

        I1_C = ufl.tr(C)
        I2_C = ufl.det(C)
        J = ufl.sqrt(I2_C)
        
        psi_1 = (self.mu/2)*(I1_C - 3)

        ln_J = ufl.ln(J)
        psi_2 = -(self.mu*ln_J)+((self.lmbda*0.5)*((ln_J)**2))

        return psi_1+psi_2

    def first_piolaStress(self, u):

        I = Identity(3)
        F = grad(u)+I
        C = (F.T)*F
        C = variable(C)
        
        W = self.strain_energy(C)
        S = 2*diff(W,C)

        return F*S

constitutive_models = {7: Neo_Hookean(1E6,0.3), 8: Neo_Hookean(0.5E6, 0.4)}
traction_vectors = {2: Constant([0.0, 0.0, 0.0]), 6: Constant([0.0, 0.0, 0.0])}

# Sets the forms

V = VectorFunctionSpace(mesh, "Lagrange", 2)
u_trial = TrialFunction(V)
v = TestFunction(V)
u_solution = Function(V)

bc = [DirichletBC(V, Constant((0.0, 0.0, 0.0)), boundary_markers, 3)]

internal_work = 0.0
external_work = 0.0

for physical_group, model in constitutive_models.items():
    internal_work += inner(model.first_piolaStress(u_solution),grad(v))*dx(physical_group)

for physical_group, traction in traction_vectors.items():
    external_work += dot(traction,v)*ds(physical_group)

# Solver
residual_form = internal_work-external_work
residual_derivative = derivative(residual_form , u_solution, u_trial)

Res = NonlinearVariationalProblem(residual_form, u_solution, bc, J=residual_derivative)

solver = NonlinearVariationalSolver(Res)
4

It’s been a very long time since I’ve used legacy DOLFIN, but here are some simple suggestions.

Hyperelasticity models are pretty notorious for having long JIT compile times due to their complicated formulations.

You can try a more efficient representation of the UFL form for FFC to parse using UFLACS. But I think at some point this was made the default so this might not achieve much.

from dolfin import *
import ufl_legacy as ufl

parameters['form_compiler']['representation'] = 'uflacs'
#...

Rather than defining a formulation for each constant used in your hyperelasticity model, just use a material coefficient defined piecewise constant (E.g., DG0) with its values set accordingly for the subdomains. Thereby the lines

for physical_group, model in constitutive_models.items():
    internal_work += inner(model.first_piolaStress(u_solution),grad(v))*dx(physical_group)

will be reduced to

E = dolfin.Function(DG0)
E.x.array[:] = #... set constant values according to each cell appropriately

nu = dolfin.Function(DG0)
nu.x.array[:] = #... set constant values according to each cell appropriately

internal_work = inner(first_piolaStress(u_solution, E, nu),grad(v))*dx

which is much simpler for the JIT compiler to handle (forgive my syntax I’ve forgotten most of the legacy DOLFIN interface).

Lastly, I’ve noticed large speedups after moving to DOLFINx. Particularly with my own work with the entropy formulation of the compressible Navier Stokes formulation. This used to take a good number of minutes to JIT compile with legacy DOLFIN, but the improvements of DOLFINx reduced this to seconds.

1 Like
5

Thanks for your answer. The problem with material parameters as functions is that they do not really solve my problem, because I need different material models in different physical groups of the mesh. I built this sort of dynamically programmed framework because my application has different material models assigned to physical groups defined in gmsh. I’ve supposed that the problem was all the fancy structure using classes and dictionaries, but the minimal working example also took a while to compile. Picking it up from your answer, I believe the ultimate performance jump would be to set for dolfinx anyway. Thank you again!