Initial Condition Modification

Objectives

At the end of this example you will know:

  • how to modify GEOS arrays using pygeosx

  • handle parallel communication with pygeosx_tools

Input files

This example requires an input xml and python script located at:

GEOS/examples/pygeosxExamples/sedovWithStressFunction

Description of the case

This example is derived from the sedov integrated test (GEOS/src/coreComponents/physicsSolvers/solidMechanics/integratedTests/sedov.xml), which looks at the propagation of elastic waves due to an initial stress field. The pygeosx interface is used to modify the initial conditions of the problem to something of our choosing.

XML Configuration

As before, the basic sedov input xml file for this example requires some modification in order to work with pygeosx. First, we use the advanced xml input features to include the base problem (this path may need to be updated, depending on where you run the problem).

  <Included>
    <File
      name="../../../inputFiles/solidMechanics/sedov.xml"/>
  </Included>

Next, we add a new entry to the output block Python and an entry in the Events block. Whenever the python event is triggered, GEOS will pause and return to the controlling python script.

  <Events>
    <PeriodicEvent
      name="python"
      cycleFrequency="5"
      target="/Outputs/pythonOutput"/>
  </Events>

  <Outputs>
    <Python
      name="pythonOutput"/>
  </Outputs>

Python Script

Similar to the previous example, the python script begins by importing the required packages, applying the xml preprocessor, GEOS initialization, and key search.

    # Get the MPI rank
    comm = MPI.COMM_WORLD
    rank = comm.Get_rank()

    # Initialize the code and set initial conditions
    args = preprocess_parallel()
    problem = pygeosx.initialize(rank, args)
    pygeosx.apply_initial_conditions()

    # Rather than specifying the wrapper paths explicitly,
    # search for them using a set of filters
    location_key = wrapper.get_matching_wrapper_path(problem, ['Region2', 'elementCenter'])
    stress_key = wrapper.get_matching_wrapper_path(problem, ['Region2', 'shale', 'stress'])
    ghost_key = wrapper.get_matching_wrapper_path(problem, ['Region2', 'cb1', 'ghostRank'])

The next steps rely on a python function that we use to set stress. The argument to this function, x, is assumed to be a numpy array of element centers:

def stress_fn(x):
    """
    Function to set stress values

    Args:
        x (np.ndarray) the element centers

    Returns:
        np.ndarray: stress values
    """
    R = x[:, 0]**2 + x[:, 1]**2 + x[:, 2]**2
    return np.sin(2.0 * np.pi * R / np.amax(R))

In the following section, we zero out the initial stress and then set it based on stress_fn. While doing this, we use wrapper.print_global_value_range to check on the process.

    # Print initial stress
    wrapper.print_global_value_range(problem, stress_key, 'stress')

    # Zero out stress
    wrapper.set_wrapper_to_value(problem, stress_key, 0.0)
    wrapper.print_global_value_range(problem, stress_key, 'stress')

    # Set stress via a function
    wrapper.set_wrapper_with_function(problem, stress_key, location_key, stress_fn, target_index=0)
    wrapper.set_wrapper_with_function(problem, stress_key, location_key, stress_fn, target_index=1)
    wrapper.set_wrapper_with_function(problem, stress_key, location_key, stress_fn, target_index=2)
    wrapper.print_global_value_range(problem, stress_key, 'stress')

Finally, we run the simulation. As an optional step, we extract numpy arrays from the datastructure using different parallel approaches:

    # Run the code
    while pygeosx.run() != pygeosx.COMPLETED:
        wrapper.print_global_value_range(problem, stress_key, 'stress')

        # Gather/allgather tests
        tmp = wrapper.gather_wrapper(problem, stress_key)
        print(wrapper.rank, 'gather', np.shape(tmp), flush=True)

        tmp = wrapper.allgather_wrapper(problem, stress_key)
        print(wrapper.rank, 'allgather', np.shape(tmp), flush=True)

        tmp = wrapper.allgather_wrapper(problem, stress_key, ghost_key=ghost_key)
        print(wrapper.rank, 'allgather_ghost_filtered', np.shape(tmp), flush=True)

Running the Problem

To run the problem, you must use the specific version of python where pygeosx is installed. This is likeley located here:

GEOS/[build_dir]/lib/PYGEOSX/bin/python

Note that you may need to manually install the pygeosx_tools package (and its pre-requisites) into this python distribution. To do so, you can run the following:

cd GEOS/[build_dir]/lib/PYGEOSX/bin
pip install --upgrade ../../../../src/coreComponents/python/modules/pygeosx_tools_package/

To run the code, you will call the pygeosx run script with python, and supply the typical geosx command-line arguments and any parallel arguments. For example:

# Load the correct python environment
# If you are not using a bash shell, you may need to target one of
# the other activation scripts
source GEOS/[build_dir]/lib/PYGEOSX/bin/activate

# Move to the correct directory and run
cd /path/to/problem
python run_sedov_problem.py -i modified_sedov.xml -o results

To go further

Feedback on this example

For any feedback on this example, please submit a GitHub issue on the project’s GitHub page.

For more details