Tutorial

This page outlines how to run several demos included with Ascent.

These demos assume you have an Ascent install. If you have access to Docker, the easiest way to test the waters is via the alpinedav/ascent Docker image. For more details about using this image see Running Demos using Docker. You can also build Ascent with Spack. For more details see Building with Spack.

Demo 1: First Light

Render a sample dataset using Ascent from C++ and Python

For this demo, we run some of the “First Light” examples are installed with Ascent to enable users to quickly test ascent in their build system.

C++ Example (examples/ascent/using-with-make/ascent_render_example.cpp):

//-----------------------------------------------------------------------------
///
/// file: ascent_render_example.cpp
///
//-----------------------------------------------------------------------------

#include <iostream>

#include "ascent.hpp"

#include "conduit_blueprint.hpp"

using namespace ascent;
using namespace conduit;


int main(int argc, char **argv)
{
    std::cout << ascent::about() << std::endl;

    Ascent a;

    // open ascent
    a.open();

    // create example mesh using conduit blueprint
    Node n_mesh;
    conduit::blueprint::mesh::examples::braid("hexs",
                                              10,
                                              10,
                                              10,
                                              n_mesh);
    // publish mesh to ascent
    a.publish(n_mesh);

    // declare a scene to render the dataset
    Node scenes;
    scenes["s1/plots/p1/type"] = "pseudocolor";
    scenes["s1/plots/p1/field"] = "braid";
    // Set the output file name (ascent will add ".png")
    scenes["s1/image_prefix"] = "out_ascent_render_3d";

    // setup actions
    Node actions;
    Node &add_act = actions.append();
    add_act["action"] = "add_scenes";
    add_act["scenes"] = scenes;

    actions.append()["action"] = "execute";

    // execute
    a.execute(actions);

    // close alpine
    a.close();
}

Python Example (examples/ascent/python/ascent_python_render_example.py):


"""
 file: ascent_python_render_example.py

 description:
   Demonstrates using ascent to render a pseudocolor plot.

"""

import conduit
import conduit.blueprint
import ascent

# print details about ascent
print(ascent.about())


# open ascent
a = ascent.Ascent()
a.open()


# create example mesh using conduit blueprint
n_mesh = conduit.Node()
conduit.blueprint.mesh.examples.braid("hexs",
                                      10,
                                      10,
                                      10,
                                      n_mesh)
# publish mesh to ascent
a.publish(n_mesh)

# declare a scene to render the dataset
scenes  = conduit.Node()
scenes["s1/plots/p1/type"] = "pseudocolor"
scenes["s1/plots/p1/params/field"] = "braid"
# Set the output file name (ascent will add ".png")
scenes["s1/image_prefix"] = "out_ascent_render_3d"

# setup actions to 
actions = conduit.Node()
add_act =actions.append()
add_act["action"] = "add_scenes"
add_act["scenes"] = scenes

actions.append()["action"] = "execute"

# execute
a.execute(actions)

# close alpine
a.close()




These examples render the same example data set using ray casting to create a pseudocolor plot. The data set is one of the built-in Conduit Mesh Blueprint examples, in this case an unstructured mesh composed of hexagons. The Conduit C++ and Python interfaces are very similar, with the C++ interface heavily influenced by the ease of use of Python.

Demo 2: Using Pipelines and Scenes

Use Ascent’s pipelines in Cloverleaf3D to transform data and render in situ

For this demo, we will use Ascent to process data from Cloverleaf3D in situ.

To begin, make sure you are in the examples/ascent/proxies/cloverleaf3d directory of your Ascent install.

The default integration example for Cloverleaf3D sets up Ascent to volume render the energy field. Here is a snippet of the related Fortran code that specifies the default actions:

      execute_act = conduit_node_append(sim_actions)
      CALL conduit_node_set_path_char8_str(execute_act,"action", "execute")

      reset_act = conduit_node_append(sim_actions)
      CALL conduit_node_set_path_char8_str(reset_act,"action", "reset")

      CALL ascent_publish(my_ascent, sim_data)
      CALL ascent_execute(my_ascent, sim_actions)

      CALL conduit_node_destroy(sim_actions)
      CALL conduit_node_destroy(sim_data)

      DEALLOCATE(xvel, yvel, zvel)
      DEALLOCATE(density, energy, pressure)
      DEALLOCATE(xcoords, ycoords, zcoords)

Ascent also allows you to override compiled in actions with a ascent_actions.json file. In this case, the file we provide with Cloverleaf3D mirrors the compiled in actions:

Cloverleaf3D default ascent_actions.json file (examples/ascent/proxies/cloverleaf3d-ref/ascent_actions.json):

[
  {
    "action": "add_scenes",
    "scenes":
    {
      "s1":
      {
        "plots":
        {
          "p1":
          {
            "type": "volume",
            "field": "energy"
          }
        }
      }
    }
  },

  {
    "action": "execute"
  },

  {
    "action": "reset"
  }
]

We will override the default actions to compute contours of the input data and render the result. To do this we use a pipeline. Ascent’s pipelines allow you add transforms that modify the mesh data published to Ascent. Pipeline results can be rendered in a scene or used as input to extracts.

Edit the ascent_actions.json to create a pipeline that computes contours and renders them using a pseudocolor plot:

[
   {
      "action": "add_pipelines",
      "pipelines":
      {
        "pl1":
        {
          "f1":
          {
            "type": "contour",
            "params":
            {
              "field": "velocity_x",
              "iso_values": 1.0
            }
          }
       }
     }
    },
   {
     "action": "add_scenes",
     "scenes":
     {
       "scene1":
       {
         "plots":
         {
           "plt1":
           {
             "type": "pseudocolor",
             "pipeline": "pl1",
             "field": "pressure"
           }
         }
       }
     }
   },
   {
     "action": "execute"
   },

   {
     "action": "reset"
   }
]

(Also available in install directory: examples/ascent/tutorial/demo_2/contour.json)

You can also compose more complex scenes that render both pipeline results and the published data. To demonstrate this, we combine the pseudocolor rendering of the contour results with a volume rendering of the entire mesh:

[
   {
      "action": "add_pipelines",
      "pipelines":
      {
        "pl1":
        {
          "f1":
          {
            "type": "contour",
            "params":
            {
              "field": "velocity_x",
              "iso_values": 1.0
            }
          }
       }
     }
    },
   {
     "action": "add_scenes",
     "scenes":
     {
       "scene1":
       {
         "plots":
         {
           "plt1":
           {
             "type": "pseudocolor",
             "pipeline": "pl1",
             "field": "pressure"
           },
           "plt2":
           {
             "type": "volume",
             "field": "energy"
           }
         }
       }
     }
   },
   {
     "action": "execute"
   },

   {
     "action": "reset"
   }
]

(Also available in install directory: examples/ascent/tutorial/demo_2/volume_contour.json)

Demo 3: Creating Cinema Extracts

Use Ascent to create a Cinema database from Cloverleaf3D that can be explored after the simulation finishes

In this demo, we use Ascent to render out a Cinema database of a plot in Cloverleaf3D.

Make sure you are in examples/ascent/proxies/cloverleaf3d directory of your Ascent install. Edit the Cloverleaf3D ascent_actions.json file to direct Ascent to render out a scene using the Cinema Astaire specification (Spec-A):

[
  {
    "action": "add_scenes",
    "scenes":
    {
      "s1":
      {
        "plots":
        {
          "p1":
          {
            "type": "volume",
            "field": "energy",
            "min_value" : 1,
            "max_value" : 3
          }
        },
        "renders":
        {
          "r1" :
          {
            "type" : "cinema",
            "phi" : 4,
            "theta" : 4,
            "db_name" : "clover_db"
          }
        }
      }
    }
  },

  {
    "action": "execute"
  },

  {
    "action": "reset"
  }
]

(Also available in install directory: examples/ascent/tutorial/demo_3/ascent_actions.json)

Run Cloverleaf3D with this setup and it will render several viewpoints and construct Cinema database in the current directory. You can then open this database with a Cinema viewer and interactively explore views of data set after the simulation finishes.

Demo 4: Custom Python Extracts

Use Ascent to run a Python script which computes a histogram of Cloverleaf3D’s energy field in situ

Ascent’s Python extract provides a simple path to run Python scripts for custom analysis. Ascent provides the Python environment, so Python extracts can for any host code (even those without a Python interface).

For this demo we use numpy and mpi4py to compute a histogram of Cloverleaf3D’s energy field.

Again, since we will use the Cloverleaf3D Ascent integration, make sure you are in examples/ascent/proxies/cloverleaf3d directory of your Ascent install. Then edit the ascent_actions.json file to define a single python extract that runs a script file:

[
  {
    "action": "add_extracts",
    "extracts":
    {
      "e1":
      {
        "type": "python",
         "params":
         {
           "file": "ascent_tutorial_demo_4_histogram.py"
         }
       }
    }
  },

  {
    "action": "execute"
  },

  {
    "action": "reset"
  }
]

(Also available in install directory: examples/ascent/tutorial/demo_4/ascent_actions.json)

This requests a python extract that will use an embedded python interpreter to execute ascent_tutorial_demo_4_histogram.py, which is specified using the file parameter. The python extract also supports a source parameter that allows you to pass a python script as a string.

Next, create our analysis script ascent_tutorial_demo_4_histogram.py:

###############################################################################
#
# file: ascent_tutorial_demo_4_histogram.py
#
###############################################################################

import numpy as np
from mpi4py import MPI

# obtain a mpi4py mpi comm object
comm = MPI.Comm.f2py(ascent_mpi_comm_id())

# get this MPI task's published blueprint data
mesh_data = ascent_data().child(0)

# fetch the numpy array for the energy field values
e_vals = mesh_data["fields/energy/values"]

# find the data extents of the energy field using mpi

# first get local extents
e_min, e_max = e_vals.min(), e_vals.max()

# declare vars for reduce results
e_min_all = np.zeros(1)
e_max_all = np.zeros(1)

# reduce to get global extents
comm.Allreduce(e_min, e_min_all, op=MPI.MIN)
comm.Allreduce(e_max, e_max_all, op=MPI.MAX)

# compute bins on global extents 
bins = np.linspace(e_min_all, e_max_all)

# get histogram counts for local data
hist, bin_edges = np.histogram(e_vals, bins = bins)

# declare var for reduce results
hist_all = np.zeros_like(hist)

# sum histogram counts with MPI to get final histogram
comm.Allreduce(hist, hist_all, op=MPI.SUM)

# print result on mpi task 0
if comm.Get_rank() == 0:
    print("\nEnergy extents: {} {}\n".format(e_min_all[0], e_max_all[0]))
    print("Histogram of Energy:\n")
    print("Counts:")
    print(hist_all)
    print("\nBin Edges:")
    print(bin_edges)
    print("")

(Also available in install directory: examples/ascent/tutorial/demo_4/ascent_tutorial_demo_4_histogram.py)

This script computes a basic histogram counting the number of energy field elements that fall into a set of uniform bins. It uses numpy’s histogram function and mpi4py to handle distributed-memory coordination.

Note, there are only two functions provided by ascent:

  • ascent_data()

    Returns a Conduit tree with the data published to this MPI Task.

    Conduit’s Python API mirrors its C++ API, but with leaves returned as numpy.ndarrays. For examples of how to use Conduit’s Python API, see the Conduit Python Tutorial. In this script, we simply fetch a ndarray that points to the values of a known field, energy.

  • ascent_mpi_comm_id()

    Returns the Fortran style MPI handle (an integer) of the MPI Communicator Ascent is using.

    In this script, we use this handle and mpi4py.MPI.Comm.f2py() to obtain a mpi4py Communicator.

Finally, run Cloverleaf3D:

mpiexec -n 2 ./cloverleaf3d_par

With the default clover.in settings, Ascent execute the python script every 10th cycle. The script computes the histogram of the energy field and prints a summary like the following:

Energy extents: 1.0 2.93874088025

Histogram of Energy:

Counts:
[159308   4041   1763   2441   2044   1516   1780   1712   1804   1299
   1366   1959   1668   2176   1287   1066    962   2218   1282   1006
   1606   2236   1115   1420   1185   1293   2495   1255   1191   1062
   1435   1329   2371   1619   1067   2513   3066   2124   2755   3779
   3955   4933   2666   3279   3318   3854   3123   4798   2604]

Bin Edges:
[ 1.          1.03956614  1.07913228  1.11869842  1.15826456  1.1978307
  1.23739684  1.27696298  1.31652912  1.35609526  1.3956614   1.43522754
  1.47479368  1.51435983  1.55392597  1.59349211  1.63305825  1.67262439
  1.71219053  1.75175667  1.79132281  1.83088895  1.87045509  1.91002123
  1.94958737  1.98915351  2.02871965  2.06828579  2.10785193  2.14741807
  2.18698421  2.22655035  2.26611649  2.30568263  2.34524877  2.38481491
  2.42438105  2.4639472   2.50351334  2.54307948  2.58264562  2.62221176
  2.6617779   2.70134404  2.74091018  2.78047632  2.82004246  2.8596086
  2.89917474  2.93874088]

Running Demos using Docker

If you have Docker installed you can obtain a Docker image with a ready-to-use ascent install from Docker Hub.

Fetch the latest Ascent image:

docker pull alpinedav/ascent

After the download completes, create and run a container using this image:

docker run -p 8000:8000 -p 10000:10000 -t -i alpinedav/ascent

(The -p is used to forward ports between the container and your host machine, we use these ports to allow web servers on the container to serve data to the host.)

You will now be at a bash prompt in you container.

To add the proper paths to Python and MPI to your environment run:

source ascent_docker_setup.sh

The ascent source code is at /ascent/src/, and the install is at /ascent/install-debug.

Next, try running the python example mentioned in Demo 1: First Light:

cd /ascent/install-debug/examples/ascent/python
python ascent_python_render_example.py

You should see some verbose output and out_ascent_render_3d.png will be created.

To view output files you can use a simple Python web server to expose files from the container to your host machine:

python -m SimpleHTTPServer 10000

With this server running, open up a web browser on your host machine and view localhost:10000. You should be able to click on out_ascent_render_3d.png and view the rendered result in your web browser.

You should now be ready to run the other demos, remember to use the Python web server to browse results.