Mondaic
PricingAbout UsDocumentation
Mondaic
Version:
No results found.
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 Ω⊂Rd\Omega \subset \mathbf{R}^d (d = 2 or 3), a time interval I=[0,T]I = [0, T], and a diffusion equation of the following form:
m0(x)∂tu(x,t)−∇⋅(D(x)∇u(x,t))=0, 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(x,0)=uinit(x)u(\mathbf{x},0) = u_{\text{init}}(\mathbf{x}).
Here, uu denotes the space- and time-dependent diffusive field and uinitu_{\text{init}} are describes external forces. ∂t\partial_t denotes the first time derivative and ∇\nabla the spatial gradient operator. Furthermore, the scalar parameter m0m_0 and the symmetric second-order diffusion tensor D\mathcal{D} are space-dependent coefficients.
D\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 σi\sigma_i.
For the special case of m0=1m_0 = 1 and T=1T = 1, σ\sigma corresponds to the standard deviation of the Gaussian smoothing in meters.
In the isotropic case, D\mathcal{D} simplifies to a scalar value, in which case we may re-write the diffusion equation as
m0(x)∂tu(x,t)−∇⋅(m1(x)∇u(x,t))=0, 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 m1=0.5∗σ2m_1 = 0.5 * \sigma^2 and the isotropic smoothing length σ\sigma.
Copy
%config Completer.use_jedi = False
# Standard Python packages
import matplotlib.pyplot as plt
import numpy as np
import os
import toml

# Salvus imports
from salvus.mesh import mesh_block
from salvus.mesh.unstructured_mesh import UnstructuredMesh
import salvus.flow.api
import salvus.flow.simple_config as sc

SALVUS_FLOW_SITE_NAME = os.environ.get("SITE_NAME", "token")
sg = mesh_block.generators.cartesian.rectangle_2d(
    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