--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Elasticity/HyperElasticity.m	Sun Aug 11 22:10:47 2013 +0200
@@ -0,0 +1,35 @@
+clear all
+pkg load msh
+pkg load fem-fenics
+
+fem_init_env ();
+problem = 'HyperElasticity';
+fem_create_all (problem);
+
+x = y = z = linspace (0, 1, 16);
+msho = msh3m_structured_mesh (x, y, z, 1, 1:6);
+mshd = Mesh (msho);
+
+global V  = HyperElasticity_FunctionSpace (mshd);
+
+global bc2 = DirichletBC (V, @(x,y,z) [0; ...
+ 0.5*(0.5 + (y - 0.5)*cos(pi/3) - (z-0.5)*sin(pi/3) - y);...
+ 0.5*(0.5 + (y - 0.5)*sin(pi/3) + (z-0.5)*cos(pi/3) - z)], 2);
+
+global bc1 = DirichletBC (V, @(x,y,z) [0; 0; 0], 1);
+
+global B = Expression ('B', @(x,y,z) [0.0; -0.5; 0.0]);
+global T = Expression ('T', @(x,y,z) [0.1; 0.0; 0.0]);
+global E = 10.0;
+global nu = 0.3;
+global mu = Expression ('mu', @(x,y,z) E./(2*(1+nu)));
+global lmbda = Expression ('lmbda', @(x,y,z) E*nu./((1+nu)*(1-2*nu)));
+u = Expression ('u', @(x,y,z) [0; 0; 0]);
+
+L = HyperElasticity_LinearForm (V, mu, lmbda, B, T, u);
+u0 = assemble (L, bc1, bc2);
+
+[x, fval, info] = fsolve (@f, u0, optimset ("jacobian", "on"));
+
+func = Function ('u', V, x);
+fem_plot (func);

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Elasticity/HyperElasticity.ufl	Sun Aug 11 22:10:47 2013 +0200
@@ -0,0 +1,63 @@
+# 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)

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Elasticity/f.m	Sun Aug 11 22:10:47 2013 +0200
@@ -0,0 +1,22 @@
+function [y, jac] = f (xx)
+
+  global V;
+  global mu;
+  global lmbda;
+  global T;
+  global B;
+  global bc1;
+  global bc2;
+
+  u = Function ('u', V, xx);
+
+  a = HyperElasticity_BilinearForm (V, mu, lmbda, u);
+  L = HyperElasticity_LinearForm (V, mu, lmbda, B, T, u);
+
+  if (nargout == 1)
+    [y, xx] = assemble (L, xx, bc1, bc2);
+  elseif (nargout == 2)
+    [jac, y, xx] = assemble_system (a, L, xx, bc1, bc2);
+  endif
+
+endfunction