Parameter study in 1D parameter space

Problem description

This tutorial shows how to do a parameter study by using the autonumber function in DPMBase. This tutorial treats the case in which only a 1D parameter study has to be performed, for example one wants to study the same setup but for different particle sizes. In the image below a schematic overview is given for a 1D parameter study, This tutorial explains only how to setup a parameter study, for (more) realistic geomtries or complicated setups one is referred to the other tutorials.
Definition of words:
Parameter study: The complete set of simulations that needs to be run in which parameters change value based on the dimension set. This is done by creating a exectuable of the complete parameter study that is capable of starting many simulations. Parameter study is sometimes abreviated to study.
Simulation: A single simulation that is part of the parameter study. The simulation is not a distinct executable, but run through the parameter study executable.

Schematic overview of a 1D parameter study. Definitions of i and dim1 will be given later.

Headers

The following headers are included:

#include "Mercury3D.h"
#include "Math/Helpers.h"
#include "Species/LinearViscoelasticFrictionSpecies.h"
#include <iostream>
#include <vector>

These are the headers needed for this tutorial.

Before the main function

The main class inherits from Mercury3D.

class ParameterStudy1DDemo: public Mercury3D 

All steps needed in order to create this class are explained below.
Step 1: (Private) Member variables.
The private member variables are needed in order to keep track of which simulation has to be run next, these are stored in studyNum. The total amount of simulations that need to be performed are stored in dim1_. studyNum is build up as follows:
The first entry is storing an integer giving information on the completeness of the study. If the value is 0 the parameter study is incomplete, if the value is 1 (or larger) the study is complete.
The second entry is storing the current run in the first dimension of the parameter study.

    int dim1_ = 0;
    std::vector<int> studyNum;

Step 2: Set- and Get-functions for private member variables.
The set- and get-functions are needed to access the private member variables in the public member functions of ParameterStudy1DDemo.

    void setStudyDimensions(int dim1)
    {
        dim1_ = dim1;
    }
 
    std::vector<int> getCurrentStudyNum()
    {
        return DPMBase::get1DParametersFromRunNumber(dim1_);
    }

Step 3: Create the particle species.
In this case it has been chosen to create a separate function for the particle species. In principle any contact model can be chosen here, but for the sake of simplicity a simple LinearViscoelasticFrictionSpecies is implemented.

    void createSpecies()
    {
        // Setup the contact model for the species
        LinearViscoelasticFrictionSpecies species;
        species.setDensity(2000);
        //specify contact properties
        //normal forces
        species.setStiffness(1e5);
        species.setDissipation(0.63);
        //tangential (sliding) forces
        species.setSlidingFrictionCoefficient(0.5);
        species.setSlidingStiffness(1.2e4);
        species.setSlidingDissipation(0.16);
        //tangential (rolling) torques
        species.setRollingFrictionCoefficient(0.2);
        species.setRollingStiffness(1.2e4);
        species.setRollingDissipation(6.3e-2);
        //normal (torsion/spin) torques
        species.setTorsionFrictionCoefficient(0.1);
        species.setTorsionStiffness(1.2e4);
        species.setSlidingDissipation(6.3e-2);
        speciesHandler.copyAndAddObject(species);
    }

Step 4: Setup the initial conditions.
The initial conditions are run inside the DPMBase::solve() routine. This means that before the actual simulation is started, this function is run first. This makes us able to do the following things. First we get the studyNum, as we need this information to setup the current simulation. Also the createSpecies function is called (see Step 3). Afterwards two particles are inserted. Note that in these lines the studyNum is used to change the particle radius of particle 1. the variable studyNum gets assigned integer values, starting from 1 and ranging to dim1_ (which will be set through main()).

    void setupInitialConditions() override
    {
        studyNum = getCurrentStudyNum();
        createSpecies();
 
        SphericalParticle P;
        P.setSpecies(speciesHandler.getObject(0));
        // Create 2 particles
        particleHandler.copyAndAddObject(P);
        particleHandler.copyAndAddObject(P);
        //Set up the stuff that is the same for all runs
        particleHandler.getObject(0)->setPosition(Vec3D(0.006,0.005,0.0));
        particleHandler.getObject(1)->setPosition(Vec3D(0.004,0.005,0.0));
        particleHandler.getObject(0)->setVelocity(Vec3D(-0.1,0.0,0.0));
        particleHandler.getObject(1)->setVelocity(Vec3D( 0.1,0.0,0.0));
        // We have two different sized particles for the left and right particle
        // The value gets updated based on what studyNum is, these come from the automatic counter
        particleHandler.getObject(0)->setRadius(0.0001*studyNum[1]);
        particleHandler.getObject(1)->setRadius(0.0001);
    }

Step 5: Seting up the parameter study.
This is probably the most difficult part of the tutorial, as the previous steps are similar to other tutorials that can be found. The desire is to create a smart algorithm to decide what simulation runs have to be performed. It might happen that one runs the executable twice in the same directory, or that you want to extend the previous dimension of your parameter study. In this case the previous simulations do not need to be rerun, as all the data is already known and stored. Hence in actionsBeforeTimeLoop this is checked by looking for data-files of a certain name as well as checking the studyNum. If the study is complete the parameter study is done, and hence the process is terminated by the exit code. If the data-file has been found, but the study is not yet complete the program is taking a short-cut by skipping the current simulation (without rerunning) and starting the next.
Please Note: Not all systems use the same syntax to start simulations. For this reason the entry in the functions "launchNewRun(..)" might be slightly different depending on the operating system you use.

    void actionsBeforeTimeLoop() override
    {
        std::string fileName = this->dataFile.getName();
        if (studyNum[0] > 0)
        {
            logger(INFO, "Study is complete.");
            std::exit(0);
        }
        else if (helpers::fileExists(fileName))
        {
            //If it has move on to the next run immediately
            logger(INFO,"This run has been done: %", this->dataFile.getName());
            logger(INFO,"Launching next run \n\n");
            launchNewRun("./ParameterStudy1DDemo");
            std::exit(0);
        }
    }

If none ot the above options are triggered, in other words the simulation has to be performed, the code is simply run. At the end of the simulation it is once again checked if the parameter study is completed (the case if this simulation is the last one of the study) or if a new simulation has to be started.

    void actionsAfterSolve() override
    {
        if (studyNum[0] > 0)
        {
            logger(ERROR, "Study was already completed, no new simulations ran.");
        }
        else
        {
            //If the study is not complete save the data to disk and move on
            writeRestartFile();
            logger(INFO,"Launching next run \n\n");
            launchNewRun("./ParameterStudy1DDemo");
            std::exit(0);
        }
    }

Step 6: Creating the main() function The main function contains some core function calls that setup the parameter study, but also have some information that is the same for each simulation (e.g. the system dimensions). In principle one could place everything that is equal in all simulations here. Typically one wants to perform a parameter study in which particle properties or contact models are changed, hence these can be found in the ParameterStudy1DDemo class, it is less likely that the system dimensions need to be changed. If however this is the case for you, simply move the lines setXMax(), etc. to the setupInitialConditions(). The key line in main() is problem.autonumber(), this creates the file COUNTERDONOTDEL and stores the information needed for the parameter study. This file should only be deleted in the cases: The data is backed up in a different directory (as then the executable is no longer able to make changes).
One wants to rerun the parameter study completely.
The last thing we need to do is to set the dimension os the 1D parameter study, in this case it is set to 4. This means that in total we want 4 simulations.

int main(int argc UNUSED, char *argv[] UNUSED)
{
    //Problem constructor, using class definition above
    ParameterStudy1DDemo problem;
    // Setup the study dimensions that you want to perform
    problem.setStudyDimensions(4);
    // Setup the automatic numbering system (create/read COUNTER_DO_NOT_DEL file)
    problem.autoNumber();
    // Setup output file name
    problem.setName("ParameterStudy1DDemo");
    // General simulation settings
    problem.setTimeStep(1e-5);
    problem.setTimeMax(0.1);
    problem.setSaveCount(5000);
    problem.setXMin(0.0);
    problem.setXMax(1.0);
    problem.setYMin(0.0);
    problem.setYMax(1.0);
    problem.setZMin(0.0);
    problem.setZMax(1.0);
 
    // Solve the problem
    problem.solve();
    // Give a successful exit code to the terminal
    return 0;
}

If you would like to continue to the 2D or 3D parameter study demos:
ParameterStudy2DDemo
ParameterStudy3DDemo
Or return to Overview of advanced tutorials