OOPS 2D Tutorial: Scalar Wave Equation

OOPS-2D Structure

The most basic OOPS-2D program contains the following components:

• A Domain object, which defines the space for your problem and contains communication information for how the physical domain is split across processors.
• A Grid object, which is a subdivision of the Domain into discrete points to be solved on a particular processor.
• A Solver object, such as RK4, which will numerically solve any PDE fed into it.
• An ODE object, which provides the right-hand side for a system of time-dependent ODEs (such as PDEs discretized with finite-difference methods), boundary conditions, and initial conditions. More complex codes will also implement Parameters and ParamParser classes, but these are not much more difficult to add in. As mentioned before, OOPS-2D works as hard as possible to hide all the gory details behind the scenes.

Setting Up a New Project

The build system for OOPS automatically searches all subdirectories for CMake build lists. Therefore, it's recommended that you set up your projects inside the OOPS directory itself. Begin by building a new directory tree for your project:

mkdir WaveEquation
cd WaveEquation
mkdir include src pars scripts

The include directory is where we will store all our header files (*.h), and src will contain all source files (*.cpp). The scripts and pars directories will not be used for the first section of this tutorial, but they will be used when we add custom parameters to the code.

Go ahead at this time and add a new file, CMakeLists.txt, into OOPS-2D/WaveEquation/. Add the following code to it:

cmake_minimum_required(VERSION 3.0)
project(WaveEquation)

set(WAVE_INCLUDE_FILES
     include/wave.h
   )
set(WAVE_SOURCE_FILES
    src/wave.cpp
    src/main.cpp
   )

set(SOURCE_FILES ${WAVE_INCLUDE_FILES} ${WAVE_SOURCE_FILES})
add_executable(Wave ${SOURCE_FILES})
target_include_directories(Wave PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/include)
target_include_directories(Wave PRIVATE ${CMAKE_SOURCE_DIR}/include)
target_link_libraries(Wave oops2d ${EXTRA_LIBS})

This file contains the information CMake needs to build our new project. The set lines tell CMake exactly which files our project needs to compile. The add_executable command specifies links the source files we've set to Wave, which will be the name of our final executable. Both target_include_directories commands specify additional locations CMake can look for header files. Finally, the target_link_libraries command links our Wave project to the core OOPS-2D library (which is compiled separately) and any additional libraries, such as MPI, that are specified in the primary build list for OOPS-2D.

The Scalar Wave Equation Code

The Scalar Wave Equation

Before we start writing actual code, let's talk about the problem we're solving. Because it's easy, we'll start with a simple wave equation on a two-dimensional sheet with outflow/advection boundaries. In two dimensions, the scalar wave equation takes the form

.

In order for OOPS-2D to solve this equation, we need to reduce the order of the time derivative. So, we introduce a new variable, , which we define as

.

Therefore, we now have a system of two first-order PDEs rather than a single second-order PDE:

,

.

We can now introduce our boundary conditions:

,

where the sign of the right-hand side depends on the boundary itself (positive for left, negative for right). The boundaries also apply to the y-direction, simply substituting y for x, down for left, and up for right. This guarantees that any data that hits the edge of the grid will simply propagate off.

The Wave Class

Now we're in a position to start writing code. In the include directory, create a new file, wave.h, and add the following code.

#ifndef WAVE_H
#define WAVE_H

#include <ode.h>

class Wave : public ODE{
  private:
    // Variable labels
    static const unsigned int U_PHI = 0;
    static const unsigned int U_PI  = 1;

  protected:
    virtual void rhs(std::shared_ptr<FieldMap>& fieldMap);

  public:
    Wave(Domain* d, Solver* s);
    virtual ~Wave();

    virtual void initData();
};
#endif

The ode.h file contains the declaration of the ODE class, which we must inherit from. ODE is an abstract class which provides a useful framework for things like our right-hand side, boundary conditions, and initial data. When OOPS-2D is run with multiple cores, ODE also takes care of transferring data between neighboring grids, and it also takes care of formatting data for output.

Underneath the private label, we define two constant variables, U_PHI and U_PI, which are simply indices. These are not strictly necessary, but they are good practice because they make for more readable code. For larger systems of variables, some users may prefer to declare an enumerator instead.

The protected label contains (for now) one function, rhs. This is a virtual function declared in ODE, which is used by Solver to evolve the system of equations forward in time. The std::shared_ptr<FieldMap> type is a smart pointer (specifically, an std::shared_ptr) to FieldMap, which contains the data for all the various fields our ODE will evolve. (Users unfamiliar with the concept of a "smart pointer" may want to read more about them. They are a useful feature introduced in C++11 which, used properly, offer no noticeable performance decrease over traditional pointers while greatly reducing the risk of introducing memory leaks into the code.)

Inside the public function, we define a constructor (which requires Domain and Solver objects), a destructor (which, as Wave is a derived class, is declared as virtual), and initData, which is an inherited function we can override to set up our problem.

Following this, we need to set up main.cpp inside WaveEquation/src/, where we will define the Wave class's methods. We'll start with the constructor and destructor:

#include <wave.h>
#include <iostream>
#include <cmath>

// Constructor
Wave::Wave(Domain* d, Solver* s) : ODE(2, d, s){
  addField("Evolution", nEqs, true, true);

  reallocateData();
}

// Destructor
Wave::~Wave(){

}

We make sure to include the wave.h class we just created, as well as some utility classes from the Standard Template Library. The destructor is completely empty, as any memory we allocate is handled by ODE or through a smart pointer.

The constructor itself is a little more interesting. Right off the bat, we have to make a call to the ODE constructor, which has the following prototype:

ODE::ODE(unsigned int nEqs, Domain* d, Solver* s)

The first argument is just the number of equations in our ODE system. With the more modern FieldMap interface functionality implemented, this is somewhat obsolete, but it can be a useful piece of information to hold onto for information about the primary evolved field. The second and third arguments provide pointers to the Domain and Solver objects, respectively, that will be used by the ODE object.

The next line makes a call to addField:

Result addField(std::string name, unsigned int nEqs, bool isEvolved, bool isComm)

The addField function creates a new set of variables for the ODE object. The first parameter, name, assigns a unique name to these variables so they can be accessed via FieldMap. The second parameter, nEqs, defines how many equations this particular set of fields has. The third parameter, isEvolved, is a boolean value that tells ODE whether or not this particular field needs to be evolved or not. This can be useful for defining utility fields (such as constraint equations) which shouldn't be evolved by the Solver object and will be handled manually by the user. The final parameter, isComm, tells ODE whether or not this field needs to be communicated between different processors. Again, this adds flexibility for fields whose values may not be independent of a field which is communicated.

In our particular case, we call our field "Evolution", state that it has two variables, and that it will be evolved by the Solver and communicated across processors.