Just in time compilation error

1

Hi,

In the following code I have a growing material that grows with a growth tensor g_i. Now my code is supposed to minimize the energy with respect to displacement and growth tensor entries. But it breaks with a jit.py error. Here is a simplified version of my code:

from dolfin import *

#Upload the mesh
mesh = UnitCubeMesh(8, 8, 8) #imports the converted mesh
###############################################################
#define function space     
V1 = VectorElement('P', tetrahedron,1)         #Vector elements for the mixed space
V2 = FiniteElement('R', tetrahedron,0)         #Real line elements for the mixed space
element = MixedElement([V1,V2,V2])
Z = FunctionSpace(mesh,element)             #mixed space
######################################################################
#Define functions
du = TrialFunction(Z)            # Incremental displacement
v  = TestFunction(Z)             # Test function
Sol  = Function(Z)               # Solution function
######################################################################
#Breaking the solution
u , alpha, beta = split(Sol)
######################################################################
#Elasticity parameter for the Hopzafel cylinder (intima)
mu_i, nu_i, beta_i, eta_i, rho_i= 27.9, 0.49, 170.88, 263.66, 0.51
######################################################################
#Fibers direction 
b_i = Constant((1,0,0))
######################################################################
#looping prameters
epsilon = 0.001
deps= 0.001
epsmax = 3.0
######################################################################

#Loop
while epsilon<=3.0:

 #Building the Growth tensor
   gamma = (((1+2.7*epsilon)/((1+alpha*epsilon)*(1+beta*epsilon))-1)/epsilon)
   PARS2 = as_tensor([[alpha,0,0],[0,beta,0],[0,0,gamma]])
   g_i = Identity(3)+epsilon*PARS2          #this is y growth tensor with unknowns alpha and beta 
   G_i = det(g_i)
   ginv_i = inv(g_i)    
   #####################################################################
 
   # Kinematics 
   II = Identity(3)            # Identity tensor
   F =II + grad(u)             # Deformation gradient 
   J = det(F)
   Fe_i= F*ginv_i
   Fe_i = Fe_i
   Ce_i= Fe_i.T*Fe_i           #elastic response right cauchy-green tensor
   ######################################################################

   # Invariants of deformation tensors 
   Ie_i=tr(Ce_i)
   J_i = det(Fe_i)   #Jacobian Fe
   ######################################################################

   #Stiffness
   I4_i = conditional(lt(dot(b_i, Ce_i*b_i),Constant(1)),1,dot(b_i, Ce_i*b_i))
   ######################################################################
 
   #printing the current loop variables
   print(epsilon)
   ######################################################################

   # Stored strain energy density (neo-Hookean Hopzafel model for growing material)
   psi_i = G_i*((mu_i/2)*(Ie_i - 3)+(eta_i/beta_i)*(exp(beta_i*(rho_i*(I4_i-1)**2+(1-rho_i)*(Ie_i-3)**2))-  1)+((mu_i*nu_i)/(1-2*nu_i))*(J_i-1)**2- mu_i*ln(J_i))
   ######################################################################
 
   #The energy form
   Pi = psi_i*dx
   #####################################################################
 
   #Gateaux derivatives
   Fac = derivative(Pi, Sol, v)
   Jac = derivative(Fac, Sol, du)
   ######################################################################
 
   ##In here we have no boundary conditions 
   bcs= []
   ######################################################################
 
   #solve
   solve(Fac == 0, Sol, J=Jac,
         solver_parameters={"newton_solver": {"relative_tolerance": 9e-10,   
         "absolute_tolerance": 9e-10,"maximum_iterations": 80}})

   #Break the solution into components
   u, alpha  = split(Sol)
   ######################################################################

   epsilon+=deps

I use Dolfin version 2017.2.0 on linux mint and below is the error I receive:
I have done this on a hpc cluster and still get the same error. I use Dolfin version 2017.2.0 on Linux mint. Also this code works fine in 2D.

*** Error:   Unable to perform just-in-time compilation of form.
*** Reason:  ffc.jit failed with message:
Traceback (most recent call last):
 File "/usr/lib/python2.7/dist-packages/dolfin/compilemodules/jit.py", line 142, in jit
    result = ffc.jit(ufl_object, parameters=p)
  File "/usr/lib/python2.7/dist-packages/ffc/jitcompiler.py", line 218, in jit
    module = jit_build(ufl_object, module_name, parameters)
  File "/usr/lib/python2.7/dist-packages/ffc/jitcompiler.py", line 134, in jit_build
    generate=jit_generate)
  File "/usr/lib/python2.7/dist-packages/dijitso/jit.py", line 167, in jit
header, source, dependencies = generate(jitable, name, signature, params["generator"])
 File "/usr/lib/python2.7/dist-packages/ffc/jitcompiler.py", line 67, in jit_generate
   prefix=module_name, parameters=parameters, jit=True)
 File "/usr/lib/python2.7/dist-packages/ffc/compiler.py", line 143, in compile_form
   prefix, parameters, jit)
File "/usr/lib/python2.7/dist-packages/ffc/compiler.py", line 190, in compile_ufl_objects
 ir = compute_ir(analysis, prefix, parameters, jit)
File "/usr/lib/python2.7/dist-packages/ffc/representation.py", line 186, in compute_ir
 for (form_id, fd) in enumerate(form_datas)]
File "/usr/lib/python2.7/dist-packages/ffc/representation.py", line 463, in _compute_integral_ir
 parameters)
File "/usr/lib/python2.7/dist-packages/ffc/uflacs/uflacsrepresentation.py", line 139, in     compute_integral_ir
 parameters)
 File "/usr/lib/python2.7/dist-packages/ffc/uflacs/build_uflacs_ir.py", line 353, in build_uflacs_ir
   rtol=p["table_rtol"], atol=p["table_atol"])
  File "/usr/lib/python2.7/dist-packages/ffc/uflacs/elementtables.py", line 565, in     build_optimized_tables
    unique_table_ttypes = analyse_table_types(unique_tables, rtol=rtol, atol=atol)
  File "/usr/lib/python2.7/dist-packages/ffc/uflacs/elementtables.py", line 544, in analyse_table_types
    for uname, table in unique_tables.items() }
  File "/usr/lib/python2.7/dist-packages/ffc/uflacs/elementtables.py", line 544, in <dictcomp>
    for uname, table in unique_tables.items() }
  File "/usr/lib/python2.7/dist-packages/ffc/uflacs/elementtables.py", line 517, in analyse_table_type
   elif is_quadrature_table(table, rtol=rtol, atol=atol):
 File "/usr/lib/python2.7/dist-packages/ffc/uflacs/elementtables.py", line 493, in is_quadrature_table
   Id = numpy.eye(num_points)
 File "/usr/lib/python2.7/dist-packages/numpy/lib/twodim_base.py", line 233, in eye
m = zeros((N, M), dtype=dtype)
MemoryError

Thank you in advance for your help.

2

Please take some time to carefully format your question such that code snippets and mathematical formulae are correctly presented.

Furthermore construct a minimal working example as described here.

The more succinct your question, the more likely it is that people will be able to offer help.

1 Like
3

Thank you for the tip. I edited my post with the most minimal example I could come up with.

4

So I figure out how to fix it. Apparently the problem occurs because of the inv function. I defined the inverse manually since it is the inverse of a diagonal matrix and it is not that difficult. But I still appreciate it if someone could explain why this happens.