Running the models

Contents

Running the models#

We understand plume evolution and melt generation as it ascends using the accompanying parameter file, plume_with_melt.prm:

In the cells below, it is assumed that the ASPECT directory is installed system-wide or that you are using the HubZero ASPECT installation. If not, modify the ASPECT executable to the location where it is installed.

ASPECT can be configured in either “release” mode or “debug” mode. Release mode is more optimized than debug mode, and will usually run much faster but lacks debugging symbols and more descriptive error messages. If you are making changes to ASPECT or altering parameter files, it might be a good idea to try running in debug mode in order to resolve any problems. The code below is calling the release version using aspect-release . If you want to try running debug mode, call aspect-debug instead.

! aspect-release ./prm_files/plume_with_melt.prm
Note: The provided .prm files are configured to run on a single node in a reasonable time, but a finer mesh is recommended that can resolve the melting process without introducing temperature oscillations. To do this, increase `set Initial global refinement` and `set Minimum refinement level` to a higher value. The plots below use 4 as the minimum refinement level.
import os
os.environ['VTK_DEFAULT_OPENGL_WINDOW'] = 'vtkEGLRenderWindow' # Suppresses a warning with Pyvista
# Load the relevant libraries
import pyvista as pv
import glob
from IPython.display import HTML

pv.global_theme.allow_empty_mesh = True
Note: If you are using the HubZero environment, please uncomment the following line to spin up a virtual X server to run the pyvista plotting commands.
# pv.start_xvfb()
def plot_deformation(output_dir, field, gif_file, cmap='RdBu', ncontours=5):
    '''
    Process a series of VTU files in the given output directory, returning the 'field'
    data for plotting.
    Inputs:
    -----------
    output_dir : str
        Path to the folder containing solution VTU files (expected in solution/solution-*.vtu)
    field : str
        Name of the field we want to look at in the solution
    label : str
        Label for the field being plotted
    '''

    solution_file_names = sorted(glob.glob(output_dir + 'solution/solution-*.pvtu'))

    mesh = pv.read(solution_file_names[0])

    plotter = pv.Plotter()

    # Create the first frame
    plotter.add_mesh(mesh, scalars=field, show_scalar_bar=False)
    plotter.add_axes()

    plotter.hide_axes()

    # We know that the model is XY plane so the following
    # lines just sets up the camera
    plotter.view_xy()
    plotter.camera.roll = 0
    plotter.camera.parallel_projection = True

    sargs_1 = dict(height=0.5, vertical=True, position_x=0.05, position_y=0.05)
    sargs_2 = dict(height=0.5, vertical=True, position_x=0.85, position_y=0.05)

    plotter.open_gif(gif_file, fps=4)

    # Update function for each frame
    def update(frame):
        frame_int = int(frame)
        new_mesh = pv.read(solution_file_names[frame_int])

        # Generate contour lines
        contours = new_mesh.contour(isosurfaces=ncontours, scalars='porosity')

        plotter.add_mesh(new_mesh, clim=[0.0, 1.0], scalars=field, cmap=cmap, \
                         scalar_bar_args=sargs_1)

        if contours.n_points > 0:
            plotter.add_mesh(contours, clim=[0.0, 0.5], line_width=2, cmap='viridis_r', scalar_bar_args=sargs_2)

        plotter.add_text(f"Time step: {frame}", name="tstep")
        plotter.add_text(field,  position='upper_right', font_size=14)

    n_frames = len(solution_file_names)

    for i in range(n_frames):
        update(i)
        plotter.write_frame()

    plotter.close()
plot_deformation('plume-with-melt/', 'peridotite', 'images/peridotite.gif')

Now, you would have the gifs saved in the current directory this cell visualizes the model outputs side by side for comparison

../../../_images/peridotite.gif

Exercises#

  • How does the plume evolve if you increase/decrease the shear viscosity? How is the melt generated affected?

  • The melt volume generated in the above example is ~ 20 \(km^3\), much smaller than what is expected for continental flood basalts. What parameters could you change to increase the melt volume? Try changing them and see how the plume evolves.

Appendix#

The example file provided takes about 20 minutes on 6 cores. When the model is complete, you should expect the following output on your terminal:

../../../_images/terminal_output1.png icon