cuIBM/io_8cu_source.html

 #include <sys/stat.h>

 #include "io.h"
 #include "utilities/types.h"
 #include "utilities/boundaryCondition.h"


 template <typename T>
 T toNumber(std::string str)
 {
     T num;
     std::stringstream ss(str); //turn the string into a stream
     ss >> num; //convert
     return num;
 } // toNumber


 namespace io
 {

 std::vector<std::string> &split(const std::string &s, char delim, std::vector<std::string> &elems)
 {
     std::stringstream ss(s);
     std::string item;
     while (std::getline(ss, item, delim))
     {
         elems.push_back(item);
     }
     return elems;
 } // split


 std::vector<std::string> split(const std::string &s, char delim)
 {
     std::vector<std::string> elems;
     split(s, delim, elems);
     return elems;
 } // split


 void readInputs(int argc, char **argv, parameterDB &DB, domain &D)
 {
     // get a default database
     initialiseDefaultDB(DB);

     // first pass of command line arguments
     commandLineParse1(argc, argv, DB);

     // case folder
     std::string folder = DB["inputs"]["caseFolder"].get<std::string>();
     std::string fname;

     // read the simulation file
     fname = folder + "/simParams.yaml";
     parseSimulationFile(fname, DB);

     // read the flow file
     fname = folder + "/flow.yaml";
     parseFlowFile(fname, DB);

     // read the domain file
     fname = folder + "/domain.yaml";
     parseDomainFile(fname, D);

     // read the body file
     fname = folder + "/bodies.yaml";
     parseBodiesFile(fname, DB);

     // second pass of command line -- overwrite values in DB
     commandLineParse2(argc, argv, DB);
 } // readInputs


 void initialiseDefaultDB(parameterDB &DB)
 {
     DB["inputs"] = componentParameter();
     DB["flow"] = componentParameter();
     DB["simulation"] = componentParameter();
     DB["velocitySolve"] = componentParameter();
     DB["PoissonSolve"] = componentParameter();

     // default input files
     std::string inputs = "inputs";
     DB[inputs]["caseFolder"].set<std::string>(".");
     DB[inputs]["deviceNumber"].set<int>(0);

     // flow parameters
     std::string flow = "flow";
     DB[flow]["nu"].set<real>(0.01);
     DB[flow]["uInitial"].set<real>(1.0);
     DB[flow]["vInitial"].set<real>(0.0);
     DB[flow]["numBodies"].set<int>(0);
     std::vector<body> *bodyVec = new std::vector<body>;
     DB[flow]["bodies"].set<std::vector<body> *>(bodyVec);

     // boundary conditions
     boundaryCondition **bc = new boundaryCondition*[4];
     for (int i=0; i<4; i++)
         bc[i] = new boundaryCondition[2];
     DB[flow]["boundaryConditions"].set<boundaryCondition **>(bc);

     // simulation parameters
     std::string sim = "simulation";
     DB[sim]["dt"].set<real>(0.02);
     DB[sim]["nt"].set<int>(100);
     DB[sim]["nsave"].set<int>(100);
     DB[sim]["startStep"].set<bool>(0);
     DB[sim]["convTimeScheme"].set<timeScheme>(EULER_EXPLICIT);
     DB[sim]["diffTimeScheme"].set<timeScheme>(EULER_IMPLICIT);
     DB[sim]["ibmScheme"].set<ibmScheme>(TAIRA_COLONIUS);
     DB[sim]["interpolationType"].set<interpolationType>(LINEAR);

     // velocity solver
     std::string solver = "velocitySolve";
     DB[solver]["solver"].set<std::string>("CG");
     DB[solver]["preconditioner"].set<preconditionerType>(DIAGONAL);
     DB[solver]["rTol"].set<real>(1.0E-05);
     DB[solver]["aTol"].set<real>(1.0E-50);
     DB[solver]["maxIterations"].set<int>(10000);

     // Poisson solver
     solver = "PoissonSolve";
     DB[solver]["solver"].set<std::string>("CG");
     DB[solver]["preconditioner"].set<preconditionerType>(DIAGONAL);
     DB[solver]["rTol"].set<real>(1.0E-05);
     DB[solver]["aTol"].set<real>(1.0E-50);
     DB[solver]["maxIterations"].set<int>(20000);
 } // initialiseDefaultDB


 void commandLineParse1(int argc, char **argv, parameterDB &DB)
 {
     for (int i=1; i<argc; i++)
     {
         if (strcmp(argv[i],"-directory")==0)
         {
             i++;
             DB["inputs"]["caseFolder"].set<std::string>(std::string(argv[i]));
         }
         else if (strcmp(argv[i],"-deviceNumber")==0)
         {
             i++;
             int devNum = toNumber<int>(std::string(argv[i]));
             DB["inputs"]["deviceNumber"].set<int>(devNum);
             // sets devNum as the current device for the calling host thread
             cudaSetDevice(devNum);
         }
     }
 } // commandLineParse1


 void commandLineParse2(int argc, char **argv, parameterDB &DB)
 {
     for (int i=1; i<argc; i++)
     {
         // kinematic viscosity
         if ( strcmp(argv[i],"-nu")==0 )
         {
             i++;
             DB["flow"]["nu"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // perturbation in the x-velocity
         if ( strcmp(argv[i],"-uPerturb")==0 )
         {
             i++;
             DB["flow"]["uPerturb"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // perturbation in the y-velocity
         if ( strcmp(argv[i],"-vPerturb")==0 )
         {
             i++;
             DB["flow"]["vPerturb"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // scale the CV with respect to the body
         if ( strcmp(argv[i],"-scaleCV")==0 )
         {
             i++;
             DB["simulation"]["scaleCV"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // frequency of saving the data
         if ( strcmp(argv[i],"-nsave")==0 )
         {
             i++;
             DB["simulation"]["nsave"].set<int>(toNumber<int>(std::string(argv[i])));
         }
         // total number of time steps
         if ( strcmp(argv[i],"-nt")==0 )
         {
             i++;
             DB["simulation"]["nt"].set<int>(toNumber<int>(std::string(argv[i])));
         }
         // size of time increment
         if ( strcmp(argv[i],"-dt")==0 )
         {
             i++;
             DB["simulation"]["dt"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // relative tolerance for the velocity solve
         if ( strcmp(argv[i],"-velocity-rtol")==0 )
         {
             i++;
             DB["velocitySolve"]["rTol"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // absolute tolerance for the velocity solve
         if ( strcmp(argv[i],"-velocity-atol")==0 )
         {
             i++;
             DB["velocitySolve"]["aTol"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // relative tolerance for the Poisson solve
         if ( strcmp(argv[i],"-poisson-rtol")==0 )
         {
             i++;
             DB["PoissonSolve"]["rTol"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // absolute tolerance for the Poisson solve
         if ( strcmp(argv[i],"-poisson-atol")==0 )
         {
             i++;
             DB["PoissonSolve"]["aTol"].set<real>(toNumber<real>(std::string(argv[i])));
         }
         // IBM Scheme
         if ( strcmp(argv[i],"-ibmScheme")==0 )
         {
             i++;
             if ( strcmp(argv[i],"NavierStokes")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(NAVIER_STOKES);
             else
             if ( strcmp(argv[i],"TairaColonius")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(TAIRA_COLONIUS);
             else
             if ( strcmp(argv[i],"DirectForcing")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(DIRECT_FORCING);
             else
             if ( strcmp(argv[i],"FadlunEtAl")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(FADLUN_ET_AL);
             else
             if ( strcmp(argv[i],"Diffusion")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(DIFFUSION);
             else
             if ( strcmp(argv[i],"DFModified")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(DF_MODIFIED);
             else
             if ( strcmp(argv[i],"FEAModified")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(FEA_MODIFIED);
             else
             if ( strcmp(argv[i],"DFImproved")==0 )
                 DB["simulation"]["ibmScheme"].set<ibmScheme>(DF_IMPROVED);
         }
         // interpolation type for Eulerian direct forcing methods
         if ( strcmp(argv[i],"-interpolationType")==0 )
         {
             i++;
             if ( strcmp(argv[i],"constant")==0 )
                 DB["simulation"]["interpolationType"].set<interpolationType>(CONSTANT);
             else
             if ( strcmp(argv[i],"linear")==0 )
                 DB["simulation"]["interpolationType"].set<interpolationType>(LINEAR);
             else
             if ( strcmp(argv[i],"quadratic")==0 )
                 DB["simulation"]["interpolationType"].set<interpolationType>(QUADRATIC);
         }
     }
 } // commandLineParse2


 std::string stringFromPreconditionerType(preconditionerType s)
 {
     if (s == NONE)
         return "None";
     else if (s == DIAGONAL)
         return "Diagonal";
     else if (s == SMOOTHED_AGGREGATION)
         return "Smoothed Aggregation";
     else if (s == AINV)
         return "Approximate Inverse";
     else
     {
         printf("Error: Unknown preconditionerType.\n");
         exit(-1);
     }
 } // stringFromPreconditionerType


 std::string stringFromTimeScheme(timeScheme s)
 {
     if (s == EULER_EXPLICIT)
         return "Explicit Euler Method";
     else if (s == EULER_IMPLICIT)
         return "Implicit Euler Method";
     else if (s == ADAMS_BASHFORTH_2)
         return "2nd Order Adams-Bashforth";
     else if (s == CRANK_NICOLSON)
         return "Crank-Nicolson";
     else
     {
         printf("Error: Unknown timeScheme!\n");
         exit(-1);
     }
 } // stringFromTimeScheme


 void printSimulationInfo(parameterDB &DB, domain &D)
 {
     real dt = DB["simulation"]["dt"].get<real>(),
          scaleCV = DB["simulation"]["scaleCV"].get<real>();
     int  nt = DB["simulation"]["nt"].get<int>(),
          nsave = DB["simulation"]["nsave"].get<int>(),
          startStep = DB["simulation"]["startStep"].get<int>();
     interpolationType interpType = DB["simulation"]["interpolationType"].get<interpolationType>();
     ibmScheme ibmSchm = DB["simulation"]["ibmScheme"].get<ibmScheme>();


     std::cout << '\n';

     std::cout << "\nFlow parameters" << '\n';
     std::cout << "---------------" << '\n';
     std::cout << "nu = " << DB["flow"]["nu"].get<real>() << '\n';

     std::cout << "\nDomain" << '\n';
     std::cout << "------" << '\n';
     std::cout << D.nx << " x " << D.ny << '\n';

     std::cout << "\nSimulation parameters" << '\n';
     std::cout << "---------------------" << '\n';
     std::cout << "dt = " << dt << '\n';
     std::cout << "scaleCV = " << scaleCV << '\n';
     std::cout << "startStep = " << startStep << '\n';
     std::cout << "nt = " << nt << '\n';
     std::cout << "nsave = " << nsave << '\n';
     std::cout << "Convection time scheme = " << stringFromTimeScheme(DB["simulation"]["convTimeScheme"].get<timeScheme>()) << '\n';
     std::cout << "Diffusion time scheme  = " << stringFromTimeScheme(DB["simulation"]["diffTimeScheme"].get<timeScheme>()) << '\n';
     if (ibmSchm == FADLUN_ET_AL ||
         ibmSchm == DIRECT_FORCING ||
         ibmSchm == DIFFUSION ||
         ibmSchm == DF_IMPROVED ||
         ibmSchm == DF_MODIFIED ||
         ibmSchm == FEA_MODIFIED)
     {
         std::cout << "Interpolation type: ";
         switch(interpType)
         {
             case CONSTANT : std::cout << "Constant\n"; break;
             case LINEAR   : std::cout << "Linear\n"; break;
             case QUADRATIC: std::cout << "Quadratic\n"; break;
             default : std::cout << "Unknown\n"; break;
         }
     }

     std::cout << "\nVelocity Solve" << '\n';
     std::cout << "--------------" << '\n';
     std::cout << "Solver = " << DB["velocitySolve"]["solver"].get<std::string>() << '\n';
     std::cout << "Preconditioner = " << stringFromPreconditionerType(DB["velocitySolve"]["preconditioner"].get<preconditionerType>()) << '\n';
     std::cout << "Relative tolerance = " << DB["velocitySolve"]["rTol"].get<real>() << '\n';
     std::cout << "Absolute tolerance = " << DB["velocitySolve"]["aTol"].get<real>() << '\n';

     std::cout << "\nPoisson Solve" << '\n';
     std::cout << "-------------" << '\n';
     std::cout << "Solver = " << DB["PoissonSolve"]["solver"].get<std::string>() << '\n';
     std::cout << "Preconditioner = " << stringFromPreconditionerType(DB["PoissonSolve"]["preconditioner"].get<preconditionerType>()) << '\n';
     std::cout << "Relative tolerance = " << DB["PoissonSolve"]["rTol"].get<real>() << '\n';
     std::cout << "Absolute tolerance = " << DB["PoissonSolve"]["aTol"].get<real>() << '\n';

     std::cout << "\nOutput parameters" << '\n';
     std::cout << "-----------------" << '\n';
     std::cout << "Output folder = " << DB["inputs"]["caseFolder"].get<std::string>() << '\n';
     std::cout << "nsave = " << DB["simulation"]["nsave"].get<int>() << '\n';

     cudaDeviceProp deviceProp;
     int gpu = DB["inputs"]["deviceNumber"].get<int>();
     cudaGetDeviceProperties(&deviceProp, gpu);
     std::cout << "\nDevice Properties" << '\n';
     std::cout << "-----------------" << '\n';
     std::cout << "Name = " << deviceProp.name << '\n';
     std::cout << "Number = " << gpu << '\n';
     std::string ecc = deviceProp.ECCEnabled ? "yes" : "no";
     std::cout << "Compute capability = " << deviceProp.major << "." << deviceProp.minor << '\n';
     std::cout << "ECC Enabled = " << ecc << std::endl;
 } // printSimulationInfo


 void printTimingInfo(Logger &logger)
 {
     logger.printAllTime();
     std::cout << std::endl;
 } // printTimingInfo


 void writeGrid(std::string &caseFolder, domain &D)
 {
     std::stringstream out;
     out << caseFolder << "/grid";
     std::ofstream file(out.str().c_str(), std::ios::binary);
     file.write((char*)(&D.nx), sizeof(int));
     file.write((char*)(&D.x[0]), (D.nx+1)*sizeof(real));
     file.write((char*)(&D.ny), sizeof(int));
     file.write((char*)(&D.y[0]), (D.ny+1)*sizeof(real));
     file.close();
 } // writeGrid


 template <>
 void writeData<vecH>(std::string &caseFolder, int n, vecH &q, vecH &lambda, domain &D)//, bodies &B)
 {
     std::string path;
     std::stringstream out;
     int N;

     out << caseFolder << '/' << std::setfill('0') << std::setw(7) << n;
     path = out.str();

     mkdir(path.c_str(), S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);

     out.str("");
     out << path << "/q";
     std::ofstream file(out.str().c_str(), std::ios::binary);
     N = q.size();
     file.write((char*)(&N), sizeof(int));
     file.write((char*)(&q[0]), N*sizeof(real));
     file.close();

     out.str("");
     out << path << "/lambda";
     file.open(out.str().c_str(), std::ios::binary);
     N = lambda.size();
     file.write((char*)(&N), sizeof(int));
     file.write((char*)(&lambda[0]), N*sizeof(real));
     file.close();

     std::cout << "Data saved to folder " << path << std::endl;
 } // writeData


 template <>
 void writeData<vecD>(std::string &caseFolder, int n, vecD &q, vecD &lambda, domain &D)//, bodies &B)
 {
     vecH qH = q,
          lambdaH = lambda;

     writeData(caseFolder, n, qH, lambdaH, D);
 } // writeData


 void readData(std::string &caseFolder, int timeStep, real *x, std::string name)
 {
     std::stringstream in;
     std::string inFilePath;
     int n;

     in << caseFolder << "/" << std::setfill('0') << std::setw(7) << timeStep << "/" << name;
     inFilePath = in.str();
     std::cout << "Reading fluxes from " << inFilePath << " ... ";
     std::ifstream inFile(inFilePath.c_str(), std::ifstream::binary);
     inFile.read((char*)(&n), sizeof(int));
     inFile.read((char*)(&x[0]), n*sizeof(real));
     inFile.close();
     std::cout << "done" << std::endl;
 } // readData


 void printDeviceMemoryUsage(std::string label)
 {
     size_t _free, _total;
     cudaMemGetInfo(&_free, &_total);
     std::cout << label << ": Memory Usage " << std::setprecision(3) << (_total-_free)/(1024.0*1024*1024) \
               << " / " << std::setprecision(3) << _total/(1024.0*1024*1024) << " GB" << std::setprecision(6) << '\n' << std::endl;
 } // printDeviceMemoryUsage

 } // End of namespace io
Logger
Monitors the time spent to achieve a certain task.
Definition: logger.h:33

domain::x
vecH x
x-coordinates of the nodes
Definition: domain.h:22

io::readInputs
void readInputs(int argc, char **argv, parameterDB &DB, domain &D)
Reads data inputs from the command-line and the simulation files.
Definition: io.cu:83

FEA_MODIFIED
Fadlun et al. modified.
Definition: types.h:83

types.h
Definition of custom types required by the code.

io::writeData< vecH >
void writeData< vecH >(std::string &caseFolder, int n, vecH &q, vecH &lambda, domain &D)
Writes numerical data at a given time-step (on the host).
Definition: io.cu:509

CONSTANT
constant
Definition: types.h:94

io::parseFlowFile
void parseFlowFile(std::string &flowFile, parameterDB &DB)
Parses the flow file and stores the parameters in the database.
Definition: parseFlowFile.cu:113

io::split
std::vector< std::string > & split(const std::string &s, char delim, std::vector< std::string > &elems)
Splits a string given a delimiter.
Definition: io.cu:47

Logger::printAllTime
void printAllTime()
Prints time spent for each event as well as the total time.
Definition: logger.h:185

boundaryCondition.h
Definition of the class boundaryCondition.

real
double real
Is a float or a double depending on the machine precision.
Definition: types.h:116

io::readData
void readData(std::string &caseFolder, int timeStep, real *x, std::string name)
Reads numerical data at a given time-step.
Definition: io.cu:571

DIRECT_FORCING
fully discrete direct forcing method
Definition: types.h:78

io::stringFromTimeScheme
std::string stringFromTimeScheme(timeScheme s)
Converts a timeScheme to a std::string.
Definition: io.cu:361

io::commandLineParse1
void commandLineParse1(int argc, char **argv, parameterDB &DB)
Parses the command-line to get the case folder name and the device number.
Definition: io.cu:186

ibmScheme
ibmScheme
Specifies the immersed boundary method used to solve the flow.
Definition: types.h:74

EULER_IMPLICIT
implicit Euler method (first order)
Definition: types.h:63

io.h
Declaration of the functions of the namespace io.

timeScheme
timeScheme
Specifies the numerical scheme used for time-integration.
Definition: types.h:60

io::parseSimulationFile
void parseSimulationFile(std::string &simFile, parameterDB &DB)
Parses simParams.yaml and stores the simulation parameters.
Definition: parseSimulationFile.cu:180

preconditionerType
preconditionerType
Specifies the type of preconditioner.
Definition: types.h:104

parameterDB
std::map< std::string, componentParameter > parameterDB
Map from a string to a componentParameter.
Definition: parameterDB.h:64

toNumber
T toNumber(std::string str)
Converts a string to a number.
Definition: io.cu:22

SMOOTHED_AGGREGATION
smoothed aggregation preconditioner
Definition: types.h:108

io::printTimingInfo
void printTimingInfo(Logger &logger)
Prints the time spent to execute tasks.
Definition: io.cu:469

vecH
cusp::array1d< real, host_memory > vecH
Cusp 1D array stored in the host.
Definition: types.h:154

io
Contains functions related to I/O tasks.

domain
Stores information about the computational grid.
Definition: domain.h:16

NAVIER_STOKES
no immersed bodies - Perot (1993)
Definition: types.h:76

DF_IMPROVED
Direct Forcing Improved.
Definition: types.h:84

FADLUN_ET_AL
Fadlun et al (2000)
Definition: types.h:79

io::writeData< vecD >
void writeData< vecD >(std::string &caseFolder, int n, vecD &q, vecD &lambda, domain &D)
Writes numerical data at a given time-step (on the device).
Definition: io.cu:554

QUADRATIC
quadratic
Definition: types.h:96

io::stringFromPreconditionerType
std::string stringFromPreconditionerType(preconditionerType s)
Converts a preconditionerType to a std::string.
Definition: io.cu:336

interpolationType
interpolationType
Specifies the type of interpolation.
Definition: types.h:92

io::writeGrid
void writeGrid(std::string &caseFolder, domain &D)
Writes grid-points coordinates into the file grid.
Definition: io.cu:482

TAIRA_COLONIUS
Taira & Colonius (2007)
Definition: types.h:80

domain::ny
int ny
number of cells in the y-direction
Definition: domain.h:19

io::initialiseDefaultDB
void initialiseDefaultDB(parameterDB &DB)
Initializes the database with default values.
Definition: io.cu:121

NONE
no preconditioner
Definition: types.h:106

domain::nx
int nx
number of cells in the x-direction
Definition: domain.h:19

io::parseBodiesFile
void parseBodiesFile(std::string &bodiesFile, parameterDB &DB)
Parses the bodies.yaml file and stores information about the immersed bodies.
Definition: parseBodiesFile.cu:117

AINV
approximate inverse preconditioner
Definition: types.h:109

io::printSimulationInfo
void printSimulationInfo(parameterDB &DB, domain &D)
Prints the parameters of the simulation.
Definition: io.cu:385

io::printDeviceMemoryUsage
void printDeviceMemoryUsage(std::string label)
Prints device memory usage.
Definition: io.cu:593

ADAMS_BASHFORTH_2
second-order Adams-Bashforth scheme
Definition: types.h:64

CRANK_NICOLSON
Crank-Nicolson scheme (second order)
Definition: types.h:66

LINEAR
linear
Definition: types.h:95

DIAGONAL
diagonal preconditioner
Definition: types.h:107

componentParameter
std::map< std::string, property > componentParameter
Map from a string to a property object.
Definition: parameterDB.h:57

DIFFUSION
Diffusion.
Definition: types.h:81

io::parseDomainFile
void parseDomainFile(std::string &domFile, domain &D)
Parses the domain file and generates the computational grid.
Definition: parseDomainFile.cu:135

domain::y
vecH y
y-coordinates of the nodes
Definition: domain.h:22

EULER_EXPLICIT
explicit Euler method (first order)
Definition: types.h:62

boundaryCondition
Stores the type of boundary condition and its value.
Definition: boundaryCondition.h:19

DF_MODIFIED
Direct Forcing modified.
Definition: types.h:82

io::commandLineParse2
void commandLineParse2(int argc, char **argv, parameterDB &DB)
Overwrites parameters with additional arguments of the command-line.
Definition: io.cu:214

io::writeData
void writeData(std::string &caseFolder, int n, Vector &q, Vector &lambda, domain &D)

vecD
cusp::array1d< real, device_memory > vecD
Cusp 1D array stored in the device.
Definition: types.h:161