Version:

This documentation is not for the latest stable Salvus version.

This tutorial is presented as Python code running inside a Jupyter Notebook, the recommended way to use Salvus. To run it yourself you can copy/type each individual cell or directly download the full notebook, including all required files.

# Diffusion Equation

We consider a spatial domain $\Omega \subset \mathbf{R}^d$ (d = 2 or 3), a time interval $I = [0, T]$, and a diffusion equation of the following form:

$m_0(\mathbf{x}) \partial_t u(\mathbf{x},t) - \nabla \cdot \big(\mathcal{D}(\mathbf{x}) \nabla u(\mathbf{x},t)\big) = 0,$

with initial conditions

$u(\mathbf{x},0) = u_{\text{init}}(\mathbf{x})$.

Here, $u$ denotes the space- and time-dependent diffusive field and $u_{\text{init}}$ are describes external forces. $\partial_t$ denotes the first time derivative and $\nabla$ the spatial gradient operator. Furthermore, the scalar parameter $m_0$ and the symmetric second-order diffusion tensor $\mathcal{D}$ are space-dependent coefficients.

$\mathcal{D}$ can be related to a Wiener process using the relation

$\mathcal{D} = \frac{1}{2} \sigma \sigma\,^T,$

which direction-dependent smoothing lengths $\sigma_i$.

For the special case of $m_0 = 1$ and $T = 1$, $\sigma$ corresponds to the standard deviation of the Gaussian smoothing in meters.

In the isotropic case, $\mathcal{D}$ simplifies to a scalar value, in which case we may re-write the diffusion equation as

$m_0(\mathbf{x}) \partial_t u(\mathbf{x},t) - \nabla \cdot \big(m_1(\mathbf{x}) \nabla u(\mathbf{x},t)\big) = 0,$

with $m_1 = 0.5 * \sigma^2$ and the isotropic smoothing length $\sigma$.

Copy
%config Completer.use_jedi = False
# Standard Python packages
import toml
import numpy as np
import matplotlib.pyplot as plt

# Salvus imports
from salvus.mesh.structured_grid_2D import StructuredGrid2D
from salvus.mesh.unstructured_mesh import UnstructuredMesh
import salvus.flow.api
import salvus.flow.simple_config as sc
sg = StructuredGrid2D.rectangle(nelem_x=40, nelem_y=60, max_x=4.0, max_y=6.0)
mesh = sg.get_unstructured_mesh()
mesh.find_side_sets("cartesian")
input_mesh = mesh.copy()
input_mesh.attach_field("some_field", np.random.randn(mesh.npoint))
input_mesh.map_nodal_fields_to_element_nodal()
input_mesh.write_h5("initial_values.h5")
input_mesh
/miniconda/envs/salvus/lib/python3.7/_collections_abc.py:702: MatplotlibDeprecationWarning: The global colormaps dictionary is no longer considered public API.
return len(self._mapping)
/miniconda/envs/salvus/lib/python3.7/_collections_abc.py:720: MatplotlibDeprecationWarning: The global colormaps dictionary is no longer considered public API.
yield from self._mapping

When run interactively in a Jupyter Notebook, this output cell contains a widget. Click the button below to load a preview of it. Please note that most interactive functionality does not work without a running Python kernel.
<salvus.mesh.unstructured_mesh.UnstructuredMesh at 0x7f25259ed6d0>
smoothing_length_in_meters = 0.1

mesh.attach_field("M0", np.ones_like(mesh.get_element_nodes()[:, :, 0]))
mesh.attach_field(
"M1",
0.5
* smoothing_length_in_meters ** 2
* np.ones_like(mesh.get_element_nodes()[:, :, 0]),
)
mesh.attach_field("fluid", np.ones(mesh.nelem))
mesh
When run interactively in a Jupyter Notebook, this output cell contains a widget. Click the button below to load a preview of it. Please note that most interactive functionality does not work without a running Python kernel.
<salvus.mesh.unstructured_mesh.UnstructuredMesh at 0x7f25259edb50>
sim = sc.simulation.Diffusion(mesh=mesh)

sim.domain.polynomial_order = 1

sim.physics.diffusion_equation.time_step_in_seconds = 1e-3
sim.physics.diffusion_equation.courant_number = 0.06

sim.physics.diffusion_equation.initial_values.filename = "initial_values.h5"
sim.physics.diffusion_equation.initial_values.format = "hdf5"
sim.physics.diffusion_equation.initial_values.field = "some_field"

sim.physics.diffusion_equation.final_values.filename = "out.h5"

sim.output.volume_data.filename = "diffusion.h5"
sim.output.volume_data.format = "hdf5"
sim.output.volume_data.fields = ["phi"]
sim.output.volume_data.sampling_interval_in_time_steps = 10

sim.validate()
salvus.flow.api.run(
site_name="local",
input_file=sim,
ranks=1,
output_folder="output",
get_all=True,
overwrite=True,
)
Job job_2008112023532204_73b3595308 running on local with 1 rank(s).
Site information:
* Salvus version: 0.11.14
* Floating point size: 32


* Downloaded 30.2 MB of results to output.
* Total run time: 6.67 seconds.
* Pure simulation time: 6.11 seconds.

<salvus.flow.sites.salvus_job.SalvusJob at 0x7f23ff168ad0>
mesh = UnstructuredMesh.from_h5(filename="output/out.h5")
mesh
/miniconda/envs/salvus/lib/python3.7/_collections_abc.py:702: MatplotlibDeprecationWarning: The global colormaps dictionary is no longer considered public API.
return len(self._mapping)
/miniconda/envs/salvus/lib/python3.7/_collections_abc.py:720: MatplotlibDeprecationWarning: The global colormaps dictionary is no longer considered public API.
yield from self._mapping

When run interactively in a Jupyter Notebook, this output cell contains a widget. Click the button below to load a preview of it. Please note that most interactive functionality does not work without a running Python kernel.
<salvus.mesh.unstructured_mesh.UnstructuredMesh at 0x7f23f5b38f50>
PAGE CONTENTS