Mercurial > fem-fenics-eugenio
view inst/example/Elasticity/HyperElasticity.ufl @ 159:d03d998e64af
Move examples in the ./inst directory.
author | gedeone-octave <marcovass89@hotmail.it> |
---|---|
date | Thu, 12 Sep 2013 15:22:19 +0200 |
parents | example/Elasticity/HyperElasticity.ufl@38b31ad6ab4e |
children |
line wrap: on
line source
# Copyright (C) 2009-2010 Harish Narayanan and Garth N. Wells # # This file is part of DOLFIN. # # DOLFIN is free software: you can redistribute it and/or modify # it under the terms of the GNU Lesser General Public License as published by # the Free Software Foundation, either version 3 of the License, or # (at your option) any later version. # # DOLFIN is distributed in the hope that it will be useful, # but WITHOUT ANY WARRANTY; without even the implied warranty of # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # GNU Lesser General Public License for more details. # # You should have received a copy of the GNU Lesser General Public License # along with DOLFIN. If not, see <http://www.gnu.org/licenses/>. # # Modified by Anders Logg 2011 # # First added: 2009-09-29 # Last changed: 2011-06-28 # # The bilinear form a(du, v) and linear form L(v) for # a hyperelastic model. # # Compile this form with FFC: ffc -l dolfin -feliminate_zeros -fprecompute_basis_const -fprecompute_ip_const HyperElasticity.ufl # Function spaces element = VectorElement("Lagrange", tetrahedron, 1) # Trial and test functions du = TrialFunction(element) # Incremental displacement v = TestFunction(element) # Test function # Functions u = Coefficient(element) # Displacement from previous iteration B = Coefficient(element) # Body force per unit volume T = Coefficient(element) # Traction force on the boundary # Kinematics I = Identity(element.cell().d) # Identity tensor F = I + grad(u) # Deformation gradient C = F.T*F # Right Cauchy-Green tensor # Invariants of deformation tensors Ic = tr(C) J = det(F) # Elasticity parameters mu = Constant(tetrahedron) lmbda = Constant(tetrahedron) # Stored strain energy density (compressible neo-Hookean model) psi = (mu/2)*(Ic - 3) - mu*ln(J) + (lmbda/2)*(ln(J))**2 # Total potential energy Pi = psi*dx - inner(B, u)*dx - inner(T, u)*ds # First variation of Pi (directional derivative about u in the direction of v) F = derivative(Pi, u, v) # Compute Jacobian of F J = derivative(F, u, du)