Version:

Volume output in a selected region

Volumetric wavefield visualizations are a powerful tool, allowing us to observe and understand wave behavior in two and three dimensions. However, these visualizations come with a significant challenge: they can generate massive amounts of data that strain computational resources and storage systems.
This tutorial introduces a practical solution in Salvus: recording data only in regions of interest (ROI). Think of an ROI as a spotlight that focuses on just the most important part of your simulation. While the simulation runs across your entire model, you'll only capture and store the wavefield data within your defined ROI. This targeted approach dramatically reduces data volumes. While we'll demonstrate this concept using a 2-D example for clarity, the exact same principles and workflow apply to 3-D simulations.
Key benefits you'll learn about:
  • Efficient data management through selective storage
  • Full-quality visualization of critical model areas
  • Significant reduction in storage requirements
Let's explore how to implement this feature in your wave propagation studies.
Copy
# Let's start with a few imports.
import dataclasses
import matplotlib.pyplot as plt
import numpy as np
import os

import salvus.namespace as sn
from salvus.toolbox.helpers.wavefield_output import wavefield_output_to_xarray

SALVUS_FLOW_SITE_NAME = os.environ.get("SITE_NAME", "local")
# We'll use a simple 2-D box.
domain = sn.domain.dim2.BoxDomain(x0=0.0, x1=10.0, y0=0.0, y1=5.0)

p = sn.Project.from_domain(domain=domain, path="project", load_if_exists=True)
To implement regions of interest in Salvus, we'll need to use an UnstructuredMeshSimulationConfiguration. This means we'll first create our mesh manually and then define our region of interest within it. Let's walk through this process step by step.
mesh = sn.layered_meshing.mesh_from_domain(
    domain=domain,
    model=sn.material.from_params(VP=1450.0, RHO=1000.0),
    mesh_resolution=sn.MeshResolution(
        reference_frequency=10000.0, elements_per_wavelength=1.5
    ),
)
The shape and size of your region of interest depends on what you want to study. In this tutorial, we'll keep things simple by creating a spherical region around our source location. However, you might choose different shapes for various applications – for example, you could create a spotlight region to focus on specific features, or define thin slices through the mesh to study wave propagation along particular planes.
# Create an elemental boolean field. False in this array means that Salvus Compute,
# i.e. the simulation, won't save the wavefield in this element.
roi = np.zeros(mesh.nelem, dtype=bool)

# Compute the centroids and set the roi to True in spherical regions around the
# source.
centroids = mesh.get_element_centroid()
roi[
    np.sqrt((centroids[:, 0] - 5.0) ** 2 + (centroids[:, 1] - 2.5) ** 2) < 1.5
] = True

# Lastly, this field is attached to the mesh. For Salvus Compute to interpret the
# ROI correctly, it needs to be named exactly `region_of_interest`.
mesh.attach_field("region_of_interest", roi)

# Visualize the mesh and its ROI. The lines below create a copy of the mesh
# with just one elemental field instead of all the ones present, which is
# just a trick to not visualize the other fields on the mesh. As this method
# returns another mesh, we don't capture it, and it will simply display
# the copy with only the ROI, and not save it to any variable.
dataclasses.replace(mesh, elemental_fields={"region_of_interest": roi})
<salvus.mesh.data_structures.unstructured_mesh.unstructured_mesh.UnstructuredMesh object at 0x7d5934ef2710>
Now let's bring all these elements together: we'll configure our complete simulation by integrating our mesh (with attached region of interest) and event parameters. We'll add a source, and then launch the simulation to see our limited wavefield storage in action.
p += sn.UnstructuredMeshSimulationConfiguration(
    name="simulation",
    unstructured_mesh=mesh,
    event_configuration=sn.EventConfiguration(
        wavelet=sn.simple_config.stf.Ricker(center_frequency=5000.0),
        waveform_simulation_configuration=sn.WaveformSimulationConfiguration(
            end_time_in_seconds=0.005
        ),
    ),
)

p += sn.Event(
    event_name="event_a",
    sources=[
        sn.simple_config.source.cartesian.ScalarPoint2D(x=5.0, y=2.5, f=10.0)
    ],
    receivers=[],
)

p.simulations.launch(
    simulation_configuration="simulation",
    events=p.events.list(),
    site_name=SALVUS_FLOW_SITE_NAME,
    ranks_per_job=1,
    extra_output_configuration={
        "volume_data": {
            "fields": ["phi"],
            "sampling_interval_in_time_steps": 100,
        }
    },
)
p.simulations.query(block=True)
[2024-11-20 13:22:08,286] INFO: Submitting job ...
Uploading 1 files...

🚀  Submitted job_2411201322369024_0f2a8b77a1@local
True

Visualizing the Results

Let's examine our simulation results by displaying a single time-slice of the wavefield within our region of interest. While we'll use Python for this visualization, you can also explore the data interactively using ParaView for more detailed analysis and different visualization options.
ed = p.waveforms.get("simulation", p.events.list()[0])[0]
volume_output = wavefield_output_to_xarray(
    ed.get_wavefield_output(output_type="volume", field="phi"),
    [
        np.linspace(0.0, 10.0, 200),
        np.linspace(0.0, 5.0, 100),
    ],
).T
ax = volume_output.sel(t=0.003, method="nearest").plot(
    shading="gouraud", infer_intervals=False
)
ax.axes.set_aspect("equal")
plt.show()
[2024-11-20 13:22:15,959] WARNING - salvus.mesh.algorithms.unstructured_mesh.utils: 17164 points were not claimed by enclosing elements. Depending on your use case, this may not be an issue.
PAGE CONTENTS