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.

# Basic Meshing API

## The most basic mesh to run with Salvus

This exercise shows the most basic mesh with all necessary commands to put material properties and boundaries ready to run with Salvus.
Copy
# initialize notebook
%matplotlib inline
%config Completer.use_jedi = False
import matplotlib.pyplot as plt
import numpy as np

plt.rcParams["figure.figsize"] = (10, 8)

### 2D

Set basic input parameters for a rectangular domain with homogeneous material properties.
vp = 1700  # p-wave velocity in m/s
vs = 1000  # s-wave velocity in m/s
rho = 2700  # density in kg / m^3

max_x = 10000  # domain size in x direction
max_y = 5000  # domain size in y direction

fmax = 1.0  # maximum frequency in Hz
elements_per_wavelength = 2.0  # resolution criterion
Compute element size and number of elements needed according to the resolution criterion.
hmax = vs / fmax / elements_per_wavelength
nelem_x = int(max_x / hmax)
nelem_y = int(max_y / hmax)
In the mesher, we build unstructured meshes by combining multiple structured grids. For the simplest case here, we build a single structured grid and plot it.
from salvus.mesh import mesh_block

sg = mesh_block.generators.cartesian.rectangle_2d(
nelem_x=nelem_x, nelem_y=nelem_y, max_x=max_x, max_y=max_y
)
sg.plot(show=True)
Then, we convert the structured grid to an unstructured mesh (trivial in this special case, but uses the same general routines).
m = sg.get_unstructured_mesh()
Attach the material properties, constant here, but using the general array storage on element nodes.
m.attach_field("VP", vp * np.ones([m.nelem, m.nodes_per_element]))
m.attach_field("VS", vs * np.ones([m.nelem, m.nodes_per_element]))
m.attach_field("RHO", rho * np.ones([m.nelem, m.nodes_per_element]))

# this is necessary to tell salvus that this is elastic material everywhere.
# Setting it to 1 would make it fluid.
m.attach_field("fluid", np.zeros(m.nelem))
Find outer Boundaries of the domain assuming the domain is rectangular.
m.find_side_sets(mode="cartesian")
Open a preview inline, showing all parameters stored on the mesh.
m
<salvus.mesh.unstructured_mesh.UnstructuredMesh at 0x7f35e0248610>

### 3D

The interfaces are very much the same as in 2D, so less comments.
vp = 1700  # p-wave velocity in m/s
vs = 1000  # s-wave velocity in m/s
rho = 2700  # density in kg / m^3

max_x = 10000  # domain size in x direction
max_y = 5000  # domain size in y direction
max_z = 3000  # domain size in y direction

fmax = 1.0  # maximum frequency in Hz
elements_per_wavelength = 2.0  # resolution criterion
hmax = vs / fmax / elements_per_wavelength
nelem_x = int(max_x / hmax)
nelem_y = int(max_y / hmax)
nelem_z = int(max_z / hmax)
Build a structured Grid. Don't plot it, because in 3D matplotlib is to slow.
from salvus.mesh import mesh_block

sg = mesh_block.generators.cartesian.cube_3d(
nelem_x=nelem_x,
nelem_y=nelem_y,
nelem_z=nelem_z,
max_x=max_x,
max_y=max_y,
max_z=max_z,
)
m = sg.get_unstructured_mesh()
m.attach_field("VP", vp * np.ones([m.nelem, m.nodes_per_element]))
m.attach_field("VS", vs * np.ones([m.nelem, m.nodes_per_element]))
m.attach_field("RHO", rho * np.ones([m.nelem, m.nodes_per_element]))
m.attach_field("fluid", np.zeros(m.nelem))
m.find_side_sets(mode="cartesian")
m