I deformed my original mesh using ALE.move() . When I deform the mesh the function space of the stress field defined using that mesh changes. For my application. I need the stress field on the deformed shape. so I create a new mesh and deform that new mesh instead of the original mesh. since my tensor function space is changed I cannot get the stress field on the deformed object as the deformation of the object is somewhat large. can some one suggest an option.
import dolfin as df
import numpy as np
from Plot import *
#function definitions
def eps(v):
return df.sym(df.grad(v))
def sigmaE(v, lmbda, mu):
sigE = lmbda * df.tr(eps(v)) * df.Identity(2) + 2.0 * mu * eps(v)
return sigE
def sigmaA(v, sigA, lmbda, mu):
sig = sigA + lmbda * df.tr(eps(v)) * df.Identity(2) + 2.0 * mu * eps(v)
return sig
def btm_Mid_pt(x,on_boundary,xval=5.0,yval=0.0,tolx=0.02,toly=0.03):
return not_(on_boundary) and df.near(x[0],xval,tolx) and df.near(x[1],yval,toly)
def tp_mid_pt(x,on_boundary,xval=5.0,yval=5.0,tolx=0.02,toly=0.03):
return not_(on_boundary) and df.near(x[0],xval,tolx) and df.near(x[1],yval,toly)
def sigA_init(V):
# geometry
#long = np.array(V.tabulate_dof_coordinates())[:, 0].max() - np.array(V.tabulate_dof_coordinates())[:, 0].min()
ep = np.array(V.tabulate_dof_coordinates())[:, 1].max() - np.array(V.tabulate_dof_coordinates())[:, 1].min()
# Residual stress
depth1 = 1.5
sigXX = [1000., 0.]
sigYY = [0., 0.]
sigXY = [0., 0.]
depth = df.Expression('ep - x[1]', degree=1, ep=ep,d=depth1)
sigA_xx = df.Expression('s0 + (s1 - s0) / d1 * d', degree=1, d=depth, d1=depth1, s0=sigXX[0], s1=sigXX[1])
sigA_yy = df.Expression('s0 + (s1 - s0) / d1 * d', degree=1, d=depth, d1=depth1, s0=sigYY[0], s1=sigYY[1])
#sigA_yy = df.Constant(sigYY[0])
sigA_xy = df.Expression('s0 + (s1 - s0) / d1 * d', degree=1, d=depth, d1=depth1, s0=sigXY[0], s1=sigXY[1])
sigA = df.Expression((('d <= d1 ? sigA_xx : 0.', 'd <= d1 ? sigA_xy : 0.'), \
('d <= d1 ? sigA_xy : 0.', 'd <= d1 ? sigA_yy : 0.')), \
degree=1, ep=ep, d=depth, d1=depth1, sigA_xx=sigA_xx, sigA_yy=sigA_yy, sigA_xy=sigA_xy)
return sigA
def resolution(V, sigA, lmbda, mu,bcs_1):
# Problem defintion
du = df.TrialFunction(V)
u_ = df.TestFunction(V)
a = df.inner(sigmaE(du, lmbda, mu), eps(u_)) * df.dx
l = -df.inner(sigA, eps(u_)) * df.dx
# problem solving
u = df.Function(V, name="Displacement")
df.solve(a == l, u, bcs_1)
return u
model = "plane_stress"
# Material
E = df.Constant(110000.)
nu = df.Constant(0.3)
mu = E / (2 * (1 + nu))
lmbda = (E * nu) / ((1 + nu) * (1 - 2 * nu))
if model == "plane_stress":
lmbda = 2 * mu * lmbda / (lmbda + 2 * mu)
mesh = df.RectangleMesh(df.Point(0,0),df.Point(10,5),100,50,diagonal="right")
#creating fuction spaces
V0 = df.VectorFunctionSpace(mesh, "Lagrange", degree=2)
V0sig = df.TensorFunctionSpace(mesh, 'DG', degree=1)
#creating boundqry conditions
bcs_0=df.DirichletBC(V0.sub(0), df.Constant(0.0),lambda x,on_boundary:tp_mid_pt(x,on_boundary,tolx=0.01,toly=0.01),"pointwise")
bcs_1=df.DirichletBC(V0, df.Constant((0.,0.)),lambda x,on_boundary:btm_Mid_pt(x,on_boundary,tolx=0.01,toly=0.01),"pointwise")
bcs = [bcs_0,bcs_1]
#creating new stress field
sig=sigA_init(V0)
u=resolution(V0, sig, lmbda, mu,bcs)
sig_after = sigmaA(u,sig,lmbda,mu)
#exporting stress field to paraview
u.set_allow_extrapolation(True)
file_results = df.XDMFFile('Before_ALE_move' + ".xdmf")
file_results.parameters["flush_output"] = True
file_results.parameters["functions_share_mesh"] = True
sig_ = df.Function(V0sig, name="Sig")
sig_.assign(df.project(sig_after, V0sig))#,mesh=stg_0_msh))
file_results.write(sig_, 0.)
# question - How to get deformed shape without changing the stress field and use it later in program
df.ALE.move(mesh,u)
# stress field after deformation
u.set_allow_extrapolation(True)
file_results = df.XDMFFile('After_ALE_move' + ".xdmf")
file_results.parameters["flush_output"] = True
file_results.parameters["functions_share_mesh"] = True
sig_ = df.Function(V0sig, name="Sig")
sig_.assign(df.project(sig_after, V0sig))#,mesh=stg_0_msh))
file_results.write(sig_, 0.)```
