Saint Venant model of hyperelasticity in 2D

Asked by minak

My question is dedicated to solving silmpe 2D system for Saint Venant–Kirchhoff model of hyperelasticity

-div(E) + grad(p) = 0
div u = 0

where
E = 0.5*(F.T*F - I)
F = grad(u) + I

Any suggestions how to solve it in fenics?

I have tried Mixed fem formulation with NonlinearVariationalProblem. Following code and my comments

 from dolfin import *
 mesh = Rectangle(0, 0, 1, 1, 10, 10, 'right')
 # velocity field
 V = VectorFunctionSpace(mesh, "CG", 2)
 # pressure field
 Q = FunctionSpace(mesh, "CG", 1)
 W = V * Q
 (u, p) = TrialFunctions(W)
 (v, q) = TestFunctions(W)
 def E(v):
     return sym(grad(v)) + grad(v)*grad(v).T
 F = inner(E(u), grad(v))*dx + div(v)*p*dx + q*div(u)*dx
 U = Function(W)
 dw = Function(W)
 J = derivative(F,dw)
 solve(F == 0, U, bcs=None, J=J, solver_parameters={"linear_solver": "lu"}, form_compiler_parameters={"optimize": False})
 u, p = U.split()

The Error code:

 Traceback (most recent call last):
   File "/home/minak/workspace/fenics tutorial/fenics tutorial/Working codes/ready.py", line 25, in <module>
     solve(F == 0, U, bcs=None, J=J, solver_parameters={"linear_solver": "lu"}, form_compiler_parameters={"optimize": False})
   File "/usr/lib/python2.6/dist-packages/dolfin/fem/solving.py", line 257, in solve
     _solve_varproblem(*args, **kwargs)
   File "/usr/lib/python2.6/dist-packages/dolfin/fem/solving.py", line 294, in _solve_varproblem
     form_compiler_parameters=form_compiler_parameters)
   File "/usr/lib/python2.6/dist-packages/dolfin/fem/solving.py", line 116, in __init__
     F = Form(F, form_compiler_parameters=form_compiler_parameters)
   File "/usr/lib/python2.6/dist-packages/dolfin/fem/form.py", line 54, in __init__
     common_cell)
   File "/usr/lib/python2.6/dist-packages/dolfin/compilemodules/jit.py", line 66, in mpi_jit
     return local_jit(*args, **kwargs)
   File "/usr/lib/python2.6/dist-packages/dolfin/compilemodules/jit.py", line 148, in jit
     return jit_compile(form, parameters=p, common_cell=common_cell)
   File "/usr/lib/python2.6/dist-packages/ffc/jitcompiler.py", line 73, in jit
     return jit_form(object, parameters, common_cell)
   File "/usr/lib/python2.6/dist-packages/ffc/jitcompiler.py", line 130, in jit_form
     common_cell=common_cell)
   File "/usr/lib/python2.6/dist-packages/ffc/compiler.py", line 155, in compile_form
     ir = compute_ir(analysis, parameters)
   File "/usr/lib/python2.6/dist-packages/ffc/representation.py", line 83, in compute_ir
     for (i, fd) in enumerate(form_datas)]
   File "/usr/lib/python2.6/dist-packages/ffc/representation.py", line 213, in _compute_integral_ir
     parameters)
   File "/usr/lib/python2.6/dist-packages/ffc/quadrature/quadraturerepresentation.py", line 131, in compute_integral_ir
     ir["trans_integrals"] = _transform_integrals(transformer, integrals_dict, domain_type)
   File "/usr/lib/python2.6/dist-packages/ffc/quadrature/quadraturerepresentation.py", line 298, in _transform_integrals
     terms = transformer.generate_terms(integrand)
   File "/usr/lib/python2.6/dist-packages/ffc/quadrature/quadraturetransformerbase.py", line 715, in generate_terms
     loop, entry = self._create_loop_entry(key)
   File "/usr/lib/python2.6/dist-packages/ffc/quadrature/quadraturetransformerbase.py", line 768, in _create_loop_entry
     error("Only rank 0, 1 and 2 tensors are currently supported: " + repr(key))
   File "<string>", line 1, in <lambda>
   File "/usr/lib/python2.6/dist-packages/ufl/log.py", line 148, in error
     raise self._exception_type(self._format_raw(*message))
 Exception: Only rank 0, 1 and 2 tensors are currently supported: ((0, 'j', 15, 15), (1, 'k', 15, 15), (1, 'k', 15, 15))

I have checked the ranks of used components and it looks like they are right, but there is something wrong with rank(F) i.e.
 ufl.log.UFLException: Invalid type conversion: <class 'ufl.form.Form'> can not be converted to any UFL type.

Thank you in advance.
minak

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
Praveen C
Solved:
Last query:
Last reply:
Revision history for this message
Praveen C (cpraveen) said :
#1

For non-linear problem, dont use TrialFunction. Your F does not depend on U at all. Try something like this

 U = Function(W)
(u,p) = (as_vector((U[0],U[1])), U[2])

Revision history for this message
Best Praveen C (cpraveen) said :
#2

This seems to run

from dolfin import *
mesh = Rectangle(0, 0, 1, 1, 10, 10, 'right')
# velocity field
V = VectorFunctionSpace(mesh, "CG", 2)
# pressure field
Q = FunctionSpace(mesh, "CG", 1)
W = V * Q
U = Function(W)
(u, p) = (as_vector((U[0],U[1])), U[2])
(v, q) = TestFunctions(W)
def E(v):
   return sym(grad(v)) + grad(v)*grad(v).T
F = inner(E(u), grad(v))*dx + div(v)*p*dx + q*div(u)*dx
dw = TrialFunction(W)
J = derivative(F,U,dw)
solve(F == 0, U, bcs=None, J=J, solver_parameters={"linear_solver": "lu"}, form_compiler_parameters={"optimize": False})
u, p = U.split()

Revision history for this message
minak (piotrekminak) said :
#3

Thanks Praveen C, that solved my question.

Revision history for this message
Anders Logg (logg) said :
#4

On Thu, Nov 24, 2011 at 04:25:51PM -0000, minak wrote:
> Question #179868 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/179868
>
> Status: Answered => Solved
>
> minak confirmed that the question is solved:
> Thanks Praveen C, that solved my question.

You can also try out CBC.Twist which has St-Venant and many other
models built-in.

--
Anders