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

New tutorial: Volume coupled diffusion #219

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
41 changes: 41 additions & 0 deletions volume-coupled-diffusion/README.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,41 @@
---
title: Volume coupled diffusion
permalink: tutorials-volume-coupled-diffusion.html
keywords: FEniCS, Diffusion, Volume Coupling
summary: This tutorial illustrates volume coupling with a simple example.
---

{% include note.html content="Get the [case files of this tutorial](https://github.com/precice/tutorials/tree/develop/volume-coupled-diffusion). Read how in the [tutorials introduction](https://www.precice.org/tutorials.html)." %}

## Setup

This case illustrates how to implement volume coupling in a simple toy problem. Two diffusion problems are coupled via volume terms. One domain (the source) has constant non-zero Dirichlet boundary conditions. The other domain (the drain) has Neumann boundary conditions and a zero Dirichlet boundary condition at the right edge of the domain. The quantity u flows from the source to the drain.

![Case setup of volume-coupled-diffusion case](images/tutorials-volume-coupled-diffusion-setup.png)

## Available solvers and dependencies

* FEniCS. Install [FEniCS](https://fenicsproject.org/download/) and the [FEniCS-adapter](https://github.com/precice/fenics-adapter). Additionally, you will need to have preCICE and the python bindings installed on your system.

## Running the simulation

This tutorial is for FEniCS. You can find the corresponding `run.sh` script in the folder `fenics`.

To choose whether you want to run the source or the drain solver, please provide the following command line input:

* `-s` flag will create a source.
* `-d` flag will create a drain.

For running the case, open two terminals run:

```bash
cd fenics
./run.sh -s
```

and

```bash
cd fenics
./run.sh -d
```
1 change: 1 addition & 0 deletions volume-coupled-diffusion/clean-tutorial.sh
4 changes: 4 additions & 0 deletions volume-coupled-diffusion/fenics/.gitignore
Original file line number Diff line number Diff line change
@@ -0,0 +1,4 @@
venv
*.pyc
*.log
out
6 changes: 6 additions & 0 deletions volume-coupled-diffusion/fenics/clean.sh
Original file line number Diff line number Diff line change
@@ -0,0 +1,6 @@
#!/bin/sh
set -e -u

. ../../tools/cleaning-tools.sh

clean_fenics .
Original file line number Diff line number Diff line change
@@ -0,0 +1,9 @@
{
"participant_name": "Drain",
"config_file_name": "../precice-config.xml",
"interface": {
"coupling_mesh_name": "Drain-Mesh",
"write_data_name": "Drain-Data",
"read_data_name": "Source-Data"
}
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,9 @@
{
"participant_name": "Source",
"config_file_name": "../precice-config.xml",
"interface": {
"coupling_mesh_name": "Source-Mesh",
"write_data_name": "Source-Data",
"read_data_name": "Drain-Data"
}
}
24 changes: 24 additions & 0 deletions volume-coupled-diffusion/fenics/run.sh
Original file line number Diff line number Diff line change
@@ -0,0 +1,24 @@
#!/bin/sh
set -e -u

usage() { echo "Usage: cmd [-s] [-d]" 1>&2; exit 1; }

# Check if no input argument was provided
if [ -z "$*" ] ; then
usage
fi

# Select appropriate case
while getopts ":sd" opt; do
case ${opt} in
s)
python3 volume-coupled-diffusion.py --source
;;
d)
python3 volume-coupled-diffusion.py --drain
;;
*)
usage
;;
esac
done
131 changes: 131 additions & 0 deletions volume-coupled-diffusion/fenics/volume-coupled-diffusion.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,131 @@
from fenics import Function, FunctionSpace, Expression, Constant, DirichletBC, TrialFunction, TestFunction, File, \
solve, lhs, rhs, dx, UnitSquareMesh, SubDomain, inner, grad, MeshFunction, MPI, interpolate
from fenicsprecice import Adapter
import numpy as np
import argparse


class AllDomain(SubDomain):
def inside(self, x, on_boundary):
return True


class AllBoundary(SubDomain):
def inside(self, x, on_boundary):
return on_boundary


class RightBoundary(SubDomain):
def inside(self, x, on_boundary):
return on_boundary and x[0] == 1


parser = argparse.ArgumentParser(description="Solving a volume coupled problem")
command_group = parser.add_mutually_exclusive_group(required=True)
command_group.add_argument("-s", "--source", help="create a source", dest="source", action="store_true")
command_group.add_argument("-d", "--drain", help="create a drain", dest="drain", action="store_true")
args = parser.parse_args()

if args.source:
precice = Adapter(adapter_config_filename="precice-adapter-config-source.json")
elif args.drain:
precice = Adapter(adapter_config_filename="precice-adapter-config-drain.json")

mesh = UnitSquareMesh(10, 10)
V = FunctionSpace(mesh, "P", 1)

u = TrialFunction(V)
v = TestFunction(V)
u_n = Function(V)
if args.source:
u_ini = Expression("1", degree=1)
bc = DirichletBC(V, u_ini, AllBoundary())
elif args.drain:
u_ini = Expression("0", degree=1)
bc = DirichletBC(V, u_ini, RightBoundary())


u_n = interpolate(u_ini, V)


dt = precice.initialize(AllDomain(), read_function_space=V, write_object=u_n)
volume_term = precice.create_coupling_expression()
f = Function(V)


dt_inv = Constant(1 / dt)

diffusion_source = 1
diffusion_drain = 1
if args.source:
F = dt_inv * (u - u_n) * v * dx - (f - u_ini) * v * dx + diffusion_source * inner(grad(u), grad(v)) * dx
elif args.drain:
F = dt_inv * (u - u_n) * v * dx - (f - u) * v * dx + diffusion_drain * inner(grad(u), grad(v)) * dx

# Time-stepping
u_np1 = Function(V)
if args.source:
u_n.rename("Source-Data", "")
u_np1.rename("Source-Data", "")
elif args.drain:
u_n.rename("Drain-Data", "")
u_np1.rename("Drain-Data", "")

t = 0

mesh_rank = MeshFunction("size_t", mesh, mesh.topology().dim())
if args.source:
mesh_rank.set_all(MPI.rank(MPI.comm_world) + 4)
else:
mesh_rank.set_all(MPI.rank(MPI.comm_world) + 0)
mesh_rank.rename("myRank", "")

# Generating output files
solution_out = File("out/%s.pvd" % precice.get_participant_name())
ranks = File("out/ranks%s.pvd" % precice.get_participant_name())

# output solution and reference solution at t=0, n=0
n = 0
print('output u^%d and u_ref^%d' % (n, n))
solution_out << u_n
ranks << mesh_rank

while precice.is_coupling_ongoing():

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

read_data = precice.read_data()

# Update the coupling expression with the new read data
precice.update_coupling_expression(volume_term, read_data)
f.assign(interpolate(volume_term, V))

dt_inv.assign(1 / dt)

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

# Write data to preCICE according to which problem is being solved
precice.write_data(u_np1)

dt = precice.advance(dt)

# roll back to checkpoint
if precice.is_action_required(precice.action_read_iteration_checkpoint()):
u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
u_n.assign(u_cp)
t = t_cp
n = n_cp
else: # update solution
u_n.assign(u_np1)
t += float(dt)
n += 1

if precice.is_time_window_complete():
solution_out << u_n

# Hold plot
precice.finalize()
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
61 changes: 61 additions & 0 deletions volume-coupled-diffusion/precice-config.xml
Original file line number Diff line number Diff line change
@@ -0,0 +1,61 @@
<?xml version="1.0"?>

<precice-configuration>

<log>
<sink filter="%Severity% > debug and %Rank% = 0" format="---[precice] %ColorizedSeverity% %Message%" enabled="true"/>
</log>

<solver-interface dimensions="2">

<data:scalar name="Drain-Data"/>
<data:scalar name="Source-Data"/>

<mesh name="Source-Mesh">
<use-data name="Drain-Data"/>
<use-data name="Source-Data"/>
</mesh>

<mesh name="Drain-Mesh">
<use-data name="Drain-Data"/>
<use-data name="Source-Data"/>
</mesh>

<participant name="Source">
<use-mesh name="Source-Mesh" provide="yes"/>
<use-mesh name="Drain-Mesh" from="Drain"/>
<write-data name="Source-Data" mesh="Source-Mesh"/>
<read-data name="Drain-Data" mesh="Source-Mesh"/>
<mapping:rbf-thin-plate-splines direction="read" from="Drain-Mesh" to="Source-Mesh" constraint="consistent"/>
</participant>

<participant name="Drain">
<use-mesh name="Drain-Mesh" provide="yes"/>
<use-mesh name="Source-Mesh" from="Source"/>
<write-data name="Drain-Data" mesh="Drain-Mesh"/>
<read-data name="Source-Data" mesh="Drain-Mesh"/>
<mapping:rbf-thin-plate-splines direction="read" from="Source-Mesh" to="Drain-Mesh" constraint="consistent"/>
</participant>

<m2n:sockets from="Source" to="Drain" exchange-directory=".."/>

<coupling-scheme:serial-implicit>
<participants first="Source" second="Drain"/>
<max-time value="10.0"/>
<time-window-size value="0.1"/>
<max-iterations value="100"/>
<exchange data="Source-Data" mesh="Source-Mesh" from="Source" to="Drain" />
<exchange data="Drain-Data" mesh="Drain-Mesh" from="Drain" to="Source" initialize="true"/>
<relative-convergence-measure data="Source-Data" mesh="Source-Mesh" limit="1e-5"/>
<relative-convergence-measure data="Drain-Data" mesh="Drain-Mesh" limit="1e-5"/>
<acceleration:IQN-ILS>
<data name="Drain-Data" mesh="Drain-Mesh"/>
<initial-relaxation value="0.1"/>
<max-used-iterations value="10"/>
<time-windows-reused value="5"/>
<filter type="QR2" limit="1e-3"/>
</acceleration:IQN-ILS>
</coupling-scheme:serial-implicit>

</solver-interface>
</precice-configuration>