Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update all FEniCS participants to be compatible with precice:develop (v3) #333

Merged
merged 18 commits into from
Sep 28, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
18 commits
Select commit Hold shift + click to select a range
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -27,9 +27,8 @@ def inside(self, x, on_boundary):
# Initialize preCICE
precice = fenicsprecice.Adapter(
adapter_config_filename="chemical-reaction-advection-diffusion.json")
precice_dt = precice.initialize(
coupling_subdomain=CouplingDomain(),
read_function_space=W)
precice.initialize(coupling_subdomain=CouplingDomain(), read_function_space=W)
precice_dt = precice.get_max_time_step_size()

flow_expr = precice.create_coupling_expression()

Expand Down Expand Up @@ -95,13 +94,14 @@ def inside(self, x, on_boundary):
# No implicit coupling
while precice.is_coupling_ongoing():

read_data = precice.read_data()
precice_dt = precice.get_max_time_step_size()
dt = np.min([default_dt, precice_dt])
read_data = precice.read_data(dt)
precice.update_coupling_expression(flow_expr, read_data)
flow_old.assign(flow)
flow.interpolate(flow_expr)
# If we add writing, do it here

dt = np.min([default_dt, precice_dt])
k.assign(1. / dt)

t += dt
Expand All @@ -128,6 +128,6 @@ def inside(self, x, on_boundary):
vtkfileC << u_C, t
vtkfileFlow << flow, t

precice_dt = precice.advance(dt)
precice.advance(dt)

precice.finalize()
9 changes: 4 additions & 5 deletions channel-transport-reaction/fluid-fenics/fluid.py
Original file line number Diff line number Diff line change
Expand Up @@ -118,10 +118,8 @@ def inside(self, x, on_boundary):


precice = fenicsprecice.Adapter(adapter_config_filename="fluid-config.json")
precice_dt = precice.initialize(
coupling_subdomain=CouplingDomain(),
write_object=u_)

precice.initialize(coupling_subdomain=CouplingDomain(), write_object=u_)
precice_dt = precice.get_max_time_step_size()

dt = np.min([default_dt, precice_dt])
k.assign(dt)
Expand Down Expand Up @@ -155,7 +153,8 @@ def inside(self, x, on_boundary):

vtkfile << u_

precice_dt = precice.advance(dt)
precice.advance(dt)
precice_dt = precice.get_max_time_step_size()
dt = np.min([default_dt, precice_dt])
k.assign(dt)

Expand Down
20 changes: 11 additions & 9 deletions elastic-tube-3d/solid-fenics/solid.py
Original file line number Diff line number Diff line change
Expand Up @@ -65,7 +65,8 @@ def neumann_boundary(x, on_boundary):
precice = Adapter(adapter_config_filename="precice-adapter-config-fsi-s.json")

# Initialize the coupling interface
precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=V, fixed_boundary=fixed_boundary)
precice.initialize(coupling_boundary, read_function_space=V, write_object=V, fixed_boundary=fixed_boundary)
precice_dt = precice.get_max_time_step_size()

fenics_dt = precice_dt # if fenics_dt == precice_dt, no subcycling is applied
dt = Constant(np.min([precice_dt, fenics_dt]))
Expand Down Expand Up @@ -168,15 +169,18 @@ def avg(x_old, x_new, alpha):

u_n.rename("Displacement", "")
u_np1.rename("Displacement", "")
displacement_out << u_n
displacement_out << (u_n, t)

while precice.is_coupling_ongoing():

if precice.is_action_required(precice.action_write_iteration_checkpoint()): # write checkpoint
if precice.requires_writing_checkpoint(): # write checkpoint
precice.store_checkpoint(u_n, t, n)

precice_dt = precice.get_max_time_step_size()
dt = Constant(np.min([precice_dt, fenics_dt]))

# read data from preCICE and get a new coupling expression
read_data = precice.read_data()
read_data = precice.read_data(dt)

# Update the point sources on the coupling boundary with the new read data
forces_x, forces_y, forces_z = precice.get_point_sources(read_data)
Expand All @@ -195,17 +199,15 @@ def avg(x_old, x_new, alpha):
assert (b is not b_forces)
solve(A, u_np1.vector(), b_forces)

dt = Constant(np.min([precice_dt, fenics_dt]))

# Write relative displacements to preCICE
u_delta.vector()[:] = u_np1.vector()[:] - u_n.vector()[:]
precice.write_data(u_delta)

# Call to advance coupling, also returns the optimum time step value
precice_dt = precice.advance(dt(0))
precice.advance(dt(0))

# Either revert to old step if timestep has not converged or move to next timestep
if precice.is_action_required(precice.action_read_iteration_checkpoint()): # roll back to checkpoint
if precice.requires_reading_checkpoint(): # roll back to checkpoint
u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
u_n.assign(u_cp)
t = t_cp
Expand All @@ -221,6 +223,6 @@ def avg(x_old, x_new, alpha):
displacement_out << (u_n, t)

# Plot tip displacement evolution
displacement_out << u_n
displacement_out << (u_n, t)

precice.finalize()
19 changes: 10 additions & 9 deletions flow-over-heated-plate/solid-fenics/solid.py
Original file line number Diff line number Diff line change
Expand Up @@ -115,13 +115,14 @@ def determine_heat_flux(V_g, u, k, flux):
# Adapter definition and initialization
precice = Adapter(adapter_config_filename="precice-adapter-config.json")

precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=V_flux_y)
precice.initialize(coupling_boundary, read_function_space=V, write_object=V_flux_y)

# Create a FEniCS Expression to define and control the coupling boundary values
coupling_expression = precice.create_coupling_expression()

# Assigning appropriate dt
dt = Constant(0)
precice_dt = precice.get_max_time_step_size()
dt.assign(np.min([fenics_dt, precice_dt]))

# Define variational problem
Expand Down Expand Up @@ -149,7 +150,7 @@ def determine_heat_flux(V_g, u, k, flux):

# Create output file
file_out = File("output/%s.pvd" % precice.get_participant_name())
file_out << u_n
file_out << (u_n, t)

print("output vtk for time = {}".format(float(t)))
n = 0
Expand All @@ -159,16 +160,16 @@ def determine_heat_flux(V_g, u, k, flux):

while precice.is_coupling_ongoing():

if precice.is_action_required(precice.action_write_iteration_checkpoint()): # write checkpoint
if precice.requires_writing_checkpoint(): # write checkpoint
precice.store_checkpoint(u_n, t, n)

read_data = precice.read_data()
precice_dt = precice.get_max_time_step_size()
dt.assign(np.min([fenics_dt, precice_dt]))
read_data = precice.read_data(dt)

# Update the coupling expression with the new read data
precice.update_coupling_expression(coupling_expression, read_data)

dt.assign(np.min([fenics_dt, precice_dt]))

# Compute solution
solve(a == L, u_np1, bcs)

Expand All @@ -177,9 +178,9 @@ def determine_heat_flux(V_g, u, k, flux):
fluxes_y = fluxes.sub(1) # only exchange y component of flux.
precice.write_data(fluxes_y)

precice_dt = precice.advance(dt(0))
precice.advance(dt(0))

if precice.is_action_required(precice.action_read_iteration_checkpoint()): # roll back to checkpoint
if precice.requires_reading_checkpoint(): # roll back to checkpoint
u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
u_n.assign(u_cp)
t = t_cp
Expand All @@ -193,7 +194,7 @@ def determine_heat_flux(V_g, u, k, flux):
tol = 10e-5 # we need some tolerance, since otherwise output might be skipped.
if abs((t + tol) % dt_out) < 2 * tol: # output if t is a multiple of dt_out
print("output vtk for time = {}".format(float(t)))
file_out << u_n
file_out << (u_n, t)

# Update dirichlet BC
u_D.t = t + float(dt)
Expand Down
21 changes: 11 additions & 10 deletions partitioned-heat-conduction-complex/fenics/heat.py
Original file line number Diff line number Diff line change
Expand Up @@ -103,14 +103,15 @@ def determine_gradient(V_g, u, flux):
# Initialize the adapter according to the specific participant
if problem is ProblemType.DIRICHLET:
precice = Adapter(adapter_config_filename="precice-adapter-config-D.json")
precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N_function)
precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N_function)
elif problem is ProblemType.NEUMANN:
precice = Adapter(adapter_config_filename="precice-adapter-config-N.json")
precice_dt = precice.initialize(coupling_boundary, read_function_space=V_g, write_object=u_D_function)
precice.initialize(coupling_boundary, read_function_space=V_g, write_object=u_D_function)

boundary_marker = False

dt = Constant(0)
precice_dt = precice.get_max_time_step_size()
dt.assign(np.min([fenics_dt, precice_dt]))

# Define variational problem
Expand Down Expand Up @@ -173,7 +174,7 @@ def determine_gradient(V_g, u, flux):
# output solution and reference solution at t=0, n=0
n = 0
print('output u^%d and u_ref^%d' % (n, n))
temperature_out << u_n
temperature_out << (u_n, t)
ref_out << u_ref
ranks << mesh_rank

Expand All @@ -191,10 +192,12 @@ def determine_gradient(V_g, u, flux):
while precice.is_coupling_ongoing():

# write checkpoint
if precice.is_action_required(precice.action_write_iteration_checkpoint()):
if precice.requires_writing_checkpoint():
precice.store_checkpoint(u_n, t, n)

read_data = precice.read_data()
precice_dt = precice.get_max_time_step_size()
dt.assign(np.min([fenics_dt, precice_dt]))
read_data = precice.read_data(dt)
if problem is ProblemType.DIRICHLET and (domain_part is DomainPart.CIRCULAR or domain_part is DomainPart.RECTANGLE):
# We have to data for an arbitrary point that is not on the circle, to obtain exact solution.
# See https://github.com/precice/fenics-adapter/issues/113 for details.
Expand All @@ -203,8 +206,6 @@ def determine_gradient(V_g, u, flux):
# Update the coupling expression with the new read data
precice.update_coupling_expression(coupling_expression, read_data)

dt.assign(np.min([fenics_dt, precice_dt]))

# Compute solution u^n+1, use bcs u_D^n+1, u^n and coupling bcs
solve(a == L, u_np1, bcs)

Expand All @@ -217,10 +218,10 @@ def determine_gradient(V_g, u, flux):
# Neumann problem reads flux and writes temperature on boundary to Dirichlet problem
precice.write_data(u_np1)

precice_dt = precice.advance(dt(0))
precice.advance(dt(0))

# roll back to checkpoint
if precice.is_action_required(precice.action_read_iteration_checkpoint()):
if precice.requires_reading_checkpoint():
u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
u_n.assign(u_cp)
t = t_cp
Expand All @@ -237,7 +238,7 @@ def determine_gradient(V_g, u, flux):
print('n = %d, t = %.2f: L2 error on domain = %.3g' % (n, t, error))
# output solution and reference solution at t_n+1
print('output u^%d and u_ref^%d' % (n, n))
temperature_out << u_n
temperature_out << (u_n, t)
ref_out << u_ref
error_out << error_pointwise

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@
"config_file_name": "../precice-config.xml",
"interface": {
"coupling_mesh_name": "Dirichlet-Mesh",
"write_data_name": "Flux",
"write_data_name": "Heat-Flux",
"read_data_name": "Temperature"
}
}
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,6 @@
"interface": {
"coupling_mesh_name": "Neumann-Mesh",
"write_data_name": "Temperature",
"read_data_name": "Flux"
"read_data_name": "Heat-Flux"
}
}
14 changes: 7 additions & 7 deletions partitioned-heat-conduction-complex/precice-config.xml
Original file line number Diff line number Diff line change
Expand Up @@ -6,22 +6,22 @@
</log>

<data:scalar name="Temperature"/>
<data:vector name="Flux"/>
<data:vector name="Heat-Flux"/>

<mesh name="Dirichlet-Mesh" dimensions="2">
<use-data name="Temperature"/>
<use-data name="Flux"/>
<use-data name="Heat-Flux"/>
</mesh>

<mesh name="Neumann-Mesh" dimensions="2">
<use-data name="Temperature"/>
<use-data name="Flux"/>
<use-data name="Heat-Flux"/>
</mesh>

<participant name="Dirichlet">
<provide-mesh name="Dirichlet-Mesh" />
<receive-mesh name="Neumann-Mesh" from="Neumann" />
<write-data name="Flux" mesh="Dirichlet-Mesh" />
<write-data name="Heat-Flux" mesh="Dirichlet-Mesh" />
<read-data name="Temperature" mesh="Dirichlet-Mesh" />
<mapping:nearest-projection direction="read" from="Neumann-Mesh" to="Dirichlet-Mesh" constraint="consistent" />
</participant>
Expand All @@ -30,7 +30,7 @@
<provide-mesh name="Neumann-Mesh" />
<receive-mesh name="Dirichlet-Mesh" from="Dirichlet"/>
<write-data name="Temperature" mesh="Neumann-Mesh"/>
<read-data name="Flux" mesh="Neumann-Mesh"/>
<read-data name="Heat-Flux" mesh="Neumann-Mesh"/>
<mapping:nearest-projection direction="read" from="Dirichlet-Mesh" to="Neumann-Mesh" constraint="consistent" />
</participant>

Expand All @@ -41,9 +41,9 @@
<max-time value="1.0"/>
<time-window-size value="0.1"/>
<max-iterations value="100"/>
<exchange data="Flux" mesh="Dirichlet-Mesh" from="Dirichlet" to="Neumann" />
<exchange data="Heat-Flux" mesh="Dirichlet-Mesh" from="Dirichlet" to="Neumann" />
<exchange data="Temperature" mesh="Neumann-Mesh" from="Neumann" to="Dirichlet" initialize="true"/>
<relative-convergence-measure data="Flux" mesh="Dirichlet-Mesh" limit="1e-5"/>
<relative-convergence-measure data="Heat-Flux" mesh="Dirichlet-Mesh" limit="1e-5"/>
<relative-convergence-measure data="Temperature" mesh="Neumann-Mesh" limit="1e-5"/>
<acceleration:IQN-ILS>
<data name="Temperature" mesh="Neumann-Mesh"/>
Expand Down
Loading
Loading