check and update solution variable at new time step

Question

Hi, I am solving a time marching problem. I have to compare and update variable from previous time step and new time step (Conditional statement): I am attaching minimal working code. any idea to fix this issue.

Thanks in advance.

 from dolfin import *
 mesh = UnitSquareMesh(4,4)
 plot(mesh, title = 'mesh', interactive = True)
 def eps(v):
     return sym(grad(v))
 def eps_dev(du):
     return dev(eps(du))
 def sigma(u):
     return 2.0*mu*eps(u) + lmbda*tr(eps(u))*Identity(2)
 V = VectorFunctionSpace(mesh, 'Lagrange', 1)
 u = TrialFunction(V)
 v = TestFunction(V)
 du = Function(V)

 mu    = 80.77e3
 lmbda = 121.15e3
 EDef = inner(sigma(u),grad(v))*dx


 class top(SubDomain):
    def inside(self,x,on_boundary):
         tol = 1e-10
         return abs(x[1]-1.0) < tol and on_boundary

 class bottom(SubDomain):
     def inside(self,x,on_boundary):
         tol = 1e-10
         return abs(x[1]) < tol and on_boundary
 Top = top()
 Bottom = bottom()
 bclx= DirichletBC(V.sub(0), Constant(0.0), Bottom)
 bcly = DirichletBC(V.sub(1), Constant(0.0), Bottom)
 bcty = DirichletBC(V.sub(1), Constant(1e-4), Top)
 bc_u = [bcly,bclx,  bcty]
 solve(lhs(EDef)==rhs(EDef),du,bc_u)

 #---------------------
 # Definition for energy
 #--------------------
 def en_dens(u):
     str_ele = 0.5*(grad(u) + grad(u).T)
     IC = tr(str_ele)
     ICC = tr(str_ele * str_ele)
     return (0.5*lmbda*IC**2) + mu*ICC



 WW = FunctionSpace(mesh, 'CG',1)
 Histold = Function(WW)
 HistNew = Function(WW)
 uold = Function(V)

 Hist = en_dens(du)
 Hisarray = Hist.vector().array()
 Hisarray[Hisarray < Histold ] = en_dens(du)
 updateHist = Hist.vector().set_local(Hisarray)

Comments

Hi there,

as far as I can see, there is no time stepping in the code snippet you posted. Please clarify. Thanks.

---

This is Just part of main code. In the last section I have here previous time step solution and current time step solution. I have check and update

Thanks for the reply

---

Please include the full code. Thanks.

---

I am facing difficulty in the last section only. Last two line.

Thanks

Answer 1

From what I can see, en_dens(du) is a scalar and hence it does not make sense to use the vector() method. Do you want to compute en_dens as a function of the coordinates on the whole two-dimensional domain? If so, see e.g. here or here for how to deal with functional involving traces and projections.