__author__ = "Anders Logg <logg@simula.no>"
__date__ = "2012-01-17"
__copyright__ = "Copyright (C) 2012 Anders Logg"
__license__  = "GNU LGPL version 3 or any later version"

from cbc.twist import *

class LEGOBeam(StaticHyperelasticity):

    def mesh(self):
        return Mesh("lego_beam.xml")

    def dirichlet_boundaries(self):
        left = "x[0] < 0.1*0.001 + DOLFIN_EPS"
        return [left]

    def dirichlet_values(self):
        clamp = Constant((0, 0, 0))
        return [clamp]

    def neumann_boundaries(self):
        right = "x[0] > 10.0*0.008 - 0.1*0.001 - DOLFIN_EPS"
        return [right]

    def neumann_conditions(self):
        traction = Constant((0, 0, -5e3))
        return [traction]

    def reference_density(self):
        # Density of PVC plastic
        return 1.45*1e3

    def body_force(self):
        # Acceleration of gravity
        return Constant((0, 0, -9.81*self.reference_density()))

    def material_model(self):
        # Material parameters for PVC plastic
        mu = 0.00230 * 1e9
        lmbda = 0.0105 * 1e9
        material = StVenantKirchhoff([mu, lmbda])
        return material

    def __str__(self):
        return "A LEGO beam"

# Setup problem
problem = LEGOBeam()

# Solve problem
u = problem.solve()

# Compute average displacement in z-direction
mesh = u.function_space().mesh()
displacement = assemble(u[2]*dx) / assemble(Constant(1)*dx, mesh=mesh)
print
print "Average z-displacement: %.16g" % displacement