This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

DOLFIN Function for L2 average on some elements - interpolate a defined function

–1 vote

I want to implement Perzyna viscoplastic model on FEniCS for a plane stress problem. I need to find average stress (which is a function of both elastic and viscoplastic strains) on each element to check for yield criterion. I assumed that viscoplastic strains are constant on elements (DG 0 space). Is there any built in function to find average of a defined function (for e.g. sigma11 in the following code) not a given expression on another space or to interpolate a defined function?

def eps(v):
    return sym(grad(v))
def sigma11(u,e_vp11,e_vp22):
    return (2*mu+lmbda)*(eps(u)[0,0]-e_vp11)+lmbda*(eps(u)[1,1]-e_vp22)
def sigma22(u,e_vp11,e_vp22):
    return lmbda*(eps(u)[0,0]-e_vp11)+(2*mu+lmbda)*(eps(u)[1,1]-e_vp22)
def sigma12(u,e_vp12):
    return 2*mu*(eps(u)[0,1]-e_vp12)
V_u = VectorFunctionSpace(mesh, 'Lagrange', 1)
u = Function(V_u)
v = TestFunction(V_u)
V_vp = FunctionSpace(mesh, 'DG', 0)
e_vp11 = Function(V_vp)
e_vp22 = Function(V_vp)
e_vp12 = Function(V_vp)
f=Constant((0.0,0.0))
# Linear and bilinear forms
a= (sigma11(u,e_vp11,e_vp22)*eps(v)[0,0]+ sigma22(u,e_vp11,e_vp22)*eps(v)[1,1]+ sigma12(u,e_vp12)*eps(v)[0,1])*dx
L=dot(f,v)*dx
# bcs ...
solve(a == L, u, bcs)
Vol = assemble(Constant(1.0)*dx(), mesh=mesh)
s11_avg = Function(V_vp)
s11_avg.vector()[:] = assemble(sigma11(u,e_vp11,e_vp22)*dx())/Vol
s22_avg = Function(V_vp)
s22_avg.vector()[:] = assemble(sigma22(u,e_vp11,e_vp22)*dx())/Vol 
s12_avg = Function(V_vp)
s12_avg.vector()[:] = assemble(sigma12(u,e_vp11,e_vp22)*dx())/Vol
# Loop over elements
for k in range(mesh.num_cells()):
# calculate yield function over each element
    f = s11_avg.vector()[k]*s11_avg.vector()[k] +s22_avg.vector()[k]*s22_avg.vector()[k]+s12_avg.vector()[k]*s12_avg.vector()[k]
    if f > 0.0 :
        e_vp11.vector().array()[k] = Gamma*dt*((sign(f)/Sigma_0)**N)*(0.5*f)*s11_avg.vector()[k]
asked Feb 3, 2015 by Navid Mozaffari FEniCS Novice (510 points)

You probably need to make your question a bit more clear. From what I guess you are asking, it sounds like a typical job for a projection into a DG0 space. This can be done by either using the project function or by solving auxiliary problems (which project also does internally).

Thanks for your reply, as you can see in the code, I want to find s11_avg (which is a function of DG 0 space) and it is computed based on the defined function sigma11(u,e_vp11,e_vp22), which is a function of both CG 1 and DG 0 spaces. That's why I defined Vol = assemble(Constant(1.0)*dx(), mesh=mesh) to compute integration of test function and then compute the average. I would appreciate if you explain how to solve the auxiliary problem or use project function to find this average.

...