1 Basic API

Download Jupyter Notebook

This page is automatically generated from a Jupyter Notebook. You can view and download it here and here.

Note that everything here is fully automatically tested against - and thus guaranteed to work - only for the latest repository versions of all the Salvus packages. So please update if anything works differently on your machine.

The most basic mesh to run with Salvus

This notebooks shows the most basic mesh with all necessary commands to put material properties and boundaries ready to run with Salvus.

# initialize notebook
%matplotlib inline
from __future__ import print_function
import matplotlib.pyplot as plt
import numpy as np
plt.style.use('ggplot')
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.     # maximum frequency in Hz
elements_per_wavelength = 2.  # 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)

Build a structured Grid and plot it.

from pymesher.structured_grid_2D import StructuredGrid2D

sg = StructuredGrid2D.rectangle(nelem_x=nelem_x, nelem_y=nelem_y, max_x=max_x, max_y=max_y)
sg.plot(show=True)

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.npoint))
m.attach_field('VS', vs * np.ones(m.npoint))
m.attach_field('RHO', rho * np.ones(m.npoint))
# this is currently necessary to tell salvus that this is elastic material everywhere. Set to 1 to 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')

Write to an exodus file that can be read by paraview and salvus.

m.write_exodus('basic2D.e')

Open with paraview. You propably need to enable the advanced mode clicking on the gear symbol to see the boundaries.

!~/ParaView/bin/paraview basic2D.e

3D

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.     # maximum frequency in Hz
elements_per_wavelength = 2.  # 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 pymesher.structured_grid_3D import StructuredGrid3D

sg = StructuredGrid3D.cube(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.npoint))
m.attach_field('VS', vs * np.ones(m.npoint))
m.attach_field('RHO', rho * np.ones(m.npoint))
m.find_side_sets(mode='cartesian')
m.write_exodus('basic3D.e')
!~/ParaView/bin/paraview basic3D.e