2
$\begingroup$

In a elastodynamics problem, I define the stress as

def sigma(v):
    return 2.0*mu*sym(grad(v)) + lmbda*tr(sym(grad(v)))*Identity(v.cell().d) 

I have computed the displacements at the nodes. Now I want to compute the stresses elementwise. Can anyone help?

$\endgroup$

1 Answer 1

4
$\begingroup$

The following compiles:

from dolfin import *
from random import random

lmbda = Constant(1.0)
mu = Constant(1.0)

def sigma(v):
    return 2.0*mu*sym(grad(v)) + lmbda*tr(sym(grad(v)))*Identity(v.cell().d) 

mesh = UnitCubeMesh(8,8,8)
V = VectorFunctionSpace(mesh, "CG", 1)

# Create random displacement vector
displacement = Function(V)
for i in xrange(V.dim()):
    displacement.vector()[i] = random()

T = TensorFunctionSpace(mesh, "DG", 0)
form = inner(sigma(displacement),TestFunction(T))*dx
scaling = 1.0/CellSize(mesh)
stress = Function(T)
assemble(scaling*form, tensor=stress.vector())
$\endgroup$
5
  • $\begingroup$ Hi, thnks for your comments. but when i run your code there is an error Trying to integrate expression of rank 2 with free indices (). Traceback (most recent call last): File "/media/njha/oss/Fenics/demo/NAFEMS.py", line 141, in <module> assemble(sigma(v)*TestFunction(T)*dx(0), tensor=stress.vector()) File "/usr/lib/python2.7/dist-packages/ufl/integral.py", line 194, in rmul $\endgroup$
    – nickrocks
    Commented May 21, 2013 at 10:13
  • $\begingroup$ Oops. You need to change sigma(v)*TestFunction(T) to inner(sigma(v), TestFunction(T)). Edited answer. $\endgroup$ Commented May 21, 2013 at 10:52
  • $\begingroup$ hi i tried assemble(inner(sigma(v),TestFunction(T)*dx(0)), tensor=stress.vector()) but still error is Trying to integrate expression of rank 2 with free indices (). Traceback (most recent call last): File "/media/njha/oss/Fenics/demo/NAFEMS.py", line 146, in <module> assemble(inner(sigma(v),TestFunction(T)*dx(0)), tensor=stress.vector()) $\endgroup$
    – nickrocks
    Commented May 21, 2013 at 11:35
  • $\begingroup$ I'm not entirely sure what you are doing different from me, but I have added a complete working example. $\endgroup$ Commented May 21, 2013 at 12:00
  • $\begingroup$ it works perfectly now ..thnks a ton $\endgroup$
    – nickrocks
    Commented May 21, 2013 at 12:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.