Local coordinate system in solid mechanics.

8 weeks ago by
Dear all,

I would like to try fenics for solid mechanics of composites, first by modeling a bi-layered laminate with 0 and 90° fiber directions, respectively under uniaxial extension.
My question is,
how can one rotate the strains (into the local coordinate system), operate at an element level and rotate back the stresses in the global coordinate system to perform such an analysis?
Thank you,
Community: FEniCS Project

1 Answer

8 weeks ago by
Here are some routines I wrote for myself a while ago, to do shell structure analysis using local coordinates:
def metricFromBasis(a0,a1):
    return as_matrix(((inner(a0,a0),inner(a0,a1)),

# change of basis of T w/ two lowered indices from a to e
def changeBasis(T,a0,a1,e0,e1):

    # raise indices on from basis
    a = metricFromBasis(a0,a1)
    ac = inv(a)
    a0c = ac[0,0]*a0 + ac[0,1]*a1
    a1c = ac[1,0]*a0 + ac[1,1]*a1

    # change-of-basis matrices
    ea = as_matrix(((inner(e0,a0c),inner(e0,a1c)),
    ae = ea.T

    # apply and return
    return ea*T*ae

If both of your bases are orthonormal, this can be simplified somewhat, since the metric components will always form the identity matrix, and the distinction between raised/lowered indices is unnecessary. 

The basic idea would be to use these routines like:

# Global Cartesian basis:
e0 = Constant((1.0,0.0))
e1 = Constant((0.0,1.0))

# strain in the global basis
eps = ... 

# Local basis vectors expressed in the same basis as e0,e1
a0,a1 = ... 

strainLoc = changeBasis(eps,e0,e1,a0,a1)

# stress in local basis, expressed in terms of strainLoc
stressLoc = ... 

stress = changeBasis(stressLoc,a0,a1,e0,e1)
Thanks for your reply. But, can one apply this to every single element?. One might have laminates with several plies.
written 8 weeks ago by John.  
If the change-of-basis is built into the problem symbolically, in UFL, it will apply to the stress in every element.  The definition of the local basis could be defined to depend on spatial position, which is accessed via

from dolfin import *
mesh = UnitSquareMesh(1,1)
x = SpatialCoordinate(mesh)​

if fiber orientation varies in space.  Or, if plies are being lumped together in a 2D model, the local basis could be defined based on a ply angle, and a loop over plies could be used to add up stress contributions from each ply.  You can also invoke the Jacobian matrix of the mapping from the element parameter space to physical space (which can be used to relate some sort of element-local basis to its pushforward in physical space) via

from dolfin import *
mesh = UnitSquareMesh(1,1)
J = Jacobian(mesh)​

but that should usually only be used in stabilization terms (e.g., SUPG, etc.), and would ideally not be involved in the constitutive model (although I realize that it's common, in practice, for continuum damage models to invoke this mapping for a "length scale" to try and get a mesh-independent fracture energy).
written 8 weeks ago by David Kamensky  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »