### state variables

250
views
0
9 months ago by
Hi, I would like to save a variable "Igama" which is a scaler function of solution variable u1. I want to save its previous values in time to perform some operations. I defined "updatedamage" function but this does not work. Igama_n.assign(Igama) also does not work and gives an error message as underlined in the last row. How can I assign Igama to Igama_n correctly? Thank you !

mesh = UnitCubeMesh(1, 1, 1)
V = VectorFunctionSpace(mesh, "Lagrange", 1)
VF = FunctionSpace(mesh, "Lagrange", 1)

d = len(u1)
#d = u1.geometric_dimension()
I = Identity(d) # Identity tensor
C = F.T*F # Right Cauchy-Green tensor
Cinv=inv(C)
Finv=inv(F)
J = det(F)
I1 = tr(C)
I2 = 0.5*(tr(C)**2-tr(dot(C,C)))
I3 = J**2.0
Jinv=J**(-1.0)
I1bar=I1/(I3**(1.0/3.0))
I2bar=I2/(I3**(2.0/3.0))
# Compute Igama
Igama = (2.0*(I1bar**2)-6.0*I2bar)
tolV=project(Constant("1E-14"),VF)
Igama=conditional(le(Igama,tolV), tolV, Igama) # required to avoid minus values in sqrt()
Igama = sqrt(Igama)/6.0
# Values of Igama at previous time
Igama_n=Function(VF)

# Update variables for t=t(n)
def updatedamage(Igama):
Igama_n = Igama

# Set up PDE
problem = NonlinearVariationalProblem(Fsol, u1, bcs=bc, J=dF)
solver = NonlinearVariationalSolver(problem)
#solver.parameters["linear_solver"] = "cg"
#solver.parameters["preconditioner"] = "ilu"
#prm = solver.parameters["krylov_solver"]
prm = solver.parameters["newton_solver"]
prm["absolute_tolerance"] = 1E-7
prm["relative_tolerance"] = 1E-7
prm["maximum_iterations"] = 100
#prm["newton_solver"]["relaxation_parameter"] = 0.1

# Solve
while t < T-dt:
t += dt
h.t=t
imposed.t=t
solver.solve()
updatedamage(Igama) # This definition does not update Igama_n
Igama_n.assign(Igama) # I get this error: Expects only linear combinations of Functions in the same FunctionSpaces.

Community: FEniCS Project
Hi, please include a minimal working example and use the proper syntax highlighting functions.

At least the quantities u1, t, dt, T, h and imposed are not defined in your code (I stopped looking then). The import of fenics is not declared either.

Igama_n.assign(project(Igama,VF))
written 9 months ago by klunkean

0
9 months ago by
Thank you , I tried the following:

def updatedamage(Igama, void, g, P):
Igama_np = project(Igama,VF)
Igama_n.vector()[:] = Igama_np.vector()

both worked well.

Regards,
Birkan