Hello, I am trying to model 3D metal printing process with moving Heat source. I have Robin bcs in boundary facets of working area and base plate except for bottom facets where i add Dirichlet bcs. After evaluation of Heat source as a dolfinx.fem.Function at some point using Heat.interpolate(lambda x: self.HeatSource(current_point, x)) i update dolfinx.mesh.meshtags to make cells near heat source alive
MWE: (just algorithm due to complexity)
points = [[0.,0.,0.],[0.,0.5,0.],[0.,1.,0.]]
num_steps = len(points)
cnt = 0
dt = fem.Constant(domain,0.1)
V = fem.functionspace(domain, ("Lagrange", 1))
facet_tag = get_meshtags_facet()
# to integrate specific boundaries ds(1) or ds(2)
ds = ufl.Measure("ds", domain=domain, subdomain_data=facet_tag)
# solution function
u_next = dolfinx.fem.Function(V)
u_next.name = "u_next"
# initial condition
u_n = dolfinx.fem.Function(V)
u_n.x.array[:] = 20
#Heat source
Heat = dolfinx.fem.Function(V)
u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)
F = ro * c * u * v * ufl.dx + dt * kappa * ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx\
- (ro * c * u_n + dt * Heat) * v * ufl.dx
boundary_conditions = [ BoundaryCondition("Robin" , 1, (s, H_c)), # our domain
BoundaryCondition("Dirichlet", 3, s) , # # bottom of plate
BoundaryCondition("Robin" , 5, (s, H_c))] # base platform
bcs = []
for condition in boundary_conditions:
if condition.type == "Dirichlet":
bcs.append(condition.bc)
else:
F += condition.bc
bilinear_form = fem.form(ufl.lhs(F))
linear_form = fem.form(ufl.rhs(F))
A = assemble_matrix(bilinear_form, bcs)
A.assemble()
B = create_vector(linear_form)
ksp = PETSc.KSP().create(domain.comm)
ksp.setOperators(A)
ksp.setType(PETSc.KSP.Type.PREONLY)
ksp.getPC().setType(PETSc.PC.Type.LU)
ksp.getPC().setFactorSolverType("mumps")
while cnt < num_steps:
if cnt < len(points):
current_point = point_at(cnt)
Heat.interpolate(lambda x: HeatSource(current_point, x))
# update alive cells
DACells(current_point,HeatSource.size)
# kill cells outside visited area
u_n.x.array[DACells.dead_tags] = 0
u_next.x.array[DACells.dead_tags] = 0
# Update the right hand side reusing the initial vector
with B.localForm() as loc_b:
loc_b.set(0)
assemble_vector(B, linear_form)
# Apply Dirichlet boundary condition to the vector
apply_lifting(B, [bilinear_form], [bcs])
B.ghostUpdate(
addv=PETSc.InsertMode.ADD_VALUES,
mode=PETSc.ScatterMode.REVERSE)
set_bc(B, bcs)
# Solve linear problem
ksp.solve(B, u_next.x.petsc_vec)
# to update values in different ranks
if cnt < len(points):
u_next.x.array[DACells.dead_tags] = 0
u_next.x.scatter_forward()
# Update solution at previous time step (u_n)
u_n.x.array[:] = u_next.x.array
cnt += 1
And thats how i update meshtags
num_cells = domain.topology.index_map(domain.topology.dim).size_local
midpoints = mesh.compute_midpoints(domain, domain.topology.dim, np.arange(num_cells, dtype=np.int32))
tags = np.zeros(len(domain.geometry.x),dtype = int)
def __call__(self, point,sizes):
def in_cube(x):
return np.array((((abs(x.T[0] - point[0]) <= sizes[0] / 2)
& (abs(x.T[1] - point[1]) <= sizes[1] / 2)
& (abs(x.T[2] - point[2]) <= 0.003 / 2))), dtype=np.int32)
cell_tags = dolfinx.mesh.meshtags(domain, domain.topology.dim
, np.arange(num_cells),in_cube(midpoints) )
temp = dolfinx.fem.locate_dofs_topological(V,3, cell_tags.find(1))
tags[temp] = 1
def dead_tags(self):
return np.logical_not(self.tags)
My question is how can i add time dependent boundary condition for facets of alive cells only?
How can i make sure that dead cells wonât affect my computational area? Because now i canât replicate results of other models due to strong heat losses. Are there any other suggestions to better approach for dead-alive elements method in dolfinx?
