First method, I try to solve my problem by following (1) Calculate sensitivity when the control is a fem function
The objective function is, which epsilon is the design variable
, the adjoint equation is
and the sensitivity
First, i solve the concentration in the function space V and the distribution seems good
import dolfinx
import ufl
from dolfinx.plot import vtk_mesh
import pyvista
from ufl import Measure, inner, grad, div, TestFunction, TrialFunction, adjoint, derivative
from dolfinx.mesh import CellType, create_rectangle
from dolfinx import fem, nls, io
from mpi4py import MPI
from petsc4py import PETSc
import numpy as np
import time
from matplotlib import*
from ufl import SpatialCoordinate, TestFunction, TrialFunction, dot, ds, dx, grad, inner
import dolfinx.nls.petsc
import dolfinx.fem.petsc
from dolfinx.fem import (Constant, Function, FunctionSpace,
assemble_scalar, dirichletbc, form, locate_dofs_geometrical)
from dolfinx import default_scalar_type
from dolfinx import *
#from dolfinx import has_patsc_complex
#Parameter
theta0 = 0.5;
theta = theta0;
D0 = 1.0;
k0 = -1.0;
eta = 2.0;
gamma = 2.0;
Rfl = 0.0018;
URF = 0.1;
# mesh and function space
mesh = create_rectangle(MPI.COMM_WORLD, [[0, 0], [1, 1]], [20, 20], CellType.triangle)
#domainU = mesh.create_unit_square(MPI.COMM_WORLD, 100, 100, mesh.CellType.quadrilateral)
V = fem.FunctionSpace(mesh, ("CG", 1))
c = fem.Function(V)
phi = TestFunction(V)
eps = fem.Function(V) # The eps is the design variable
eps.x.array[:] = 0.5
dx = Measure("dx")
x = SpatialCoordinate(mesh)
gu = Constant(mesh, default_scalar_type(0))
# Residual of the variational form of Poisson problem
R = inner(D0*eps**eta*grad(c), grad(phi))*dx - inner(k0*(1-eps)**gamma*c,phi)*dx - inner(gu,phi) *ds
##Dirichlet
def on_boundary(x):
return np.isclose(x[0], 0)
dofs_L = fem.locate_dofs_geometrical(V , on_boundary)
bc_L = fem.dirichletbc(default_scalar_type(1), dofs_L, V)
bcs = [bc_L]
#Forming and solving the linear system
problem_u = fem.petsc.NonlinearProblem(R, c, bcs)
solver = nls.petsc.NewtonSolver(MPI.COMM_WORLD,problem_u)
solver.rtol = 1e-16
ch = solver.solve(c)
but in the second part of solving adjoint equation has an error ,when I apply this step that the link I following above is error
#Define the objective function
J = ufl.inner(k0*(1-eps)**gamma,c)*dx
J_form = fem.form(J)
J_old = fem.assemble_scalar(J_form)
#Bilinear and linear form of the adjoint problem
lhs = ufl.adjoint(ufl.derivative(R,c))
rhs = -ufl.derivative(J,c)
#Solve adjoint equation
bc = []
problem_adj = dolfinx.fem.petsc.LinearProblem(lhs, rhs, bcs=bc, petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
lambda = problem_adj.solve()
#compute sensitivity
dJde_form = fem.form(ufl.derivative(J, eps))
dJde = fem.petsc.assemble_vector(dJde_form)
dJde.assemble()
# partial derivative of R w.r.t. eps
dRde_form = fem.form(adjoint(derivative(R, eps)))
dRde = fem.petsc.assemble_matrix(dRde_form)
dRde.assemble()
# Calculate the gradient of J with respect to eps (the epsilon)
dJde += dRde*lambda.vector
dJde.ghostUpdate(addv=PETSc.InsertMode.INSERT, mode=PETSc.ScatterMode.FORWARD)
error has occured at this line #compute sensitivity
dJde_form = fem.form(ufl.derivative(J, eps))
"TimeoutError: JIT compilation timed out, probably due to a failed previous compile.Try cleaning cache (e.g. remove /root/.cache/fenics/libffcx_forms_5eb899f22d9b05f26c47e06ab454b4561b19ac5a.c) or increase timeout option."
then i try following from this link How to increase the timeout parameter in JIT compilation in FEniCSX? but still has error.
Second method, I try to substitute variable values instead because i already have the sensitivity equation as this pic
dj_0 = fem.Function(V) # The eps is the design variable
dj_0.x.array[:] = 0
dj_0 = dolfinx.fem.Expression((k0*gamma*(1-eps)**(gamma-1)*c)+((D0*eta*eps**(eta-1))*dot(grad(c),grad(lamda)))+(k0*gamma*c*lamda*(1-eps)**(gamma-1)),V.element.interpolation_points())
print(dj_0)
the result is, which i don’t know what does it mean
<dolfinx.fem.function.Expression object at 0x7f52f46480a0>
Now i try 2 methods
- Following the first link Calculate sensitivity when the control is a fem function
- Substituting variable values in the sensitivity equation directly
My questions are
- Which methods should i follow?
- If i choose the first method how could i figure out the error?
- The result when i substitute variable values in the sensitivity equation, what does it mean?
, but can you help me fix this error again. Thank you.