|
moab
|
#include "moab/Core.hpp"#include "moab/CpuTimer.hpp"#include "DataCoupler.hpp"#include "ElemUtil.hpp"#include "iMesh.h"#include "MBiMesh.hpp"#include <iostream>#include <iomanip>#include <sstream>#include <assert.h>Go to the source code of this file.
Defines | |
| #define | STRINGIFY_(A) #A |
| #define | STRINGIFY(A) STRINGIFY_(A) |
| #define | RRA(a) |
| #define | PRINT_LAST_ERROR |
| #define | PRINT_LAST_ERR |
Functions | |
| std::string | TestDir (".") |
| void | print_usage (char **argv) |
| ErrorCode | get_file_options (int argc, char **argv, std::vector< std::string > &meshFiles, DataCoupler::Method &method, std::string &interpTag, std::string &gNormTag, std::string &ssNormTag, std::vector< const char * > &ssTagNames, std::vector< const char * > &ssTagValues, std::string &readOpts, std::string &outFile, std::string &writeOpts, std::string &dbgFile, bool &help, double &epsilon) |
| ErrorCode | test_interpolation (Interface *mbImpl, DataCoupler::Method method, std::string &interpTag, std::string &gNormTag, std::string &ssNormTag, std::vector< const char * > &ssTagNames, std::vector< const char * > &ssTagValues, iBase_EntitySetHandle *roots, std::vector< ParallelComm * > &pcs, double &instant_time, double &pointloc_time, double &interp_time, double &gnorm_time, double &ssnorm_time, double &toler) |
| void | reduceMax (double &v) |
| int | main (int argc, char **argv) |
| bool | check_for_flag (const char *str) |
| #define PRINT_LAST_ERR |
if (iBase_SUCCESS != err) {\ std::string tmp_str;\ std::cout << "Failure; message:" << std::endl;\ mbImpl->get_last_error(tmp_str);\ std::cout << tmp_str << std::endl;\ MPI_Abort(MPI_COMM_WORLD, result); \ return result;\ }
Definition at line 44 of file datacoupler_test.cpp.
| #define PRINT_LAST_ERROR |
if (MB_SUCCESS != result) {\ std::string tmp_str;\ std::cout << "Failure; message:" << std::endl;\ mbImpl->get_last_error(tmp_str);\ std::cout << tmp_str << std::endl;\ MPI_Abort(MPI_COMM_WORLD, result); \ return result;\ }
Definition at line 34 of file datacoupler_test.cpp.
| #define RRA | ( | a | ) |
if (MB_SUCCESS != result) {\ std::string tmp_str; mbImpl->get_last_error(tmp_str);\ tmp_str.append("\n"); tmp_str.append(a);\ dynamic_cast<Core*>(mbImpl)->get_error_handler()->set_last_error(tmp_str.c_str()); \ return result;}
Definition at line 28 of file datacoupler_test.cpp.
| #define STRINGIFY | ( | A | ) | STRINGIFY_(A) |
Definition at line 21 of file datacoupler_test.cpp.
| #define STRINGIFY_ | ( | A | ) | #A |
Definition at line 20 of file datacoupler_test.cpp.
| bool check_for_flag | ( | const char * | str | ) |
Definition at line 314 of file datacoupler_test.cpp.
{
if (str[0] == '-')
return true;
else
return false;
}
| ErrorCode get_file_options | ( | int | argc, |
| char ** | argv, | ||
| std::vector< std::string > & | meshFiles, | ||
| DataCoupler::Method & | method, | ||
| std::string & | interpTag, | ||
| std::string & | gNormTag, | ||
| std::string & | ssNormTag, | ||
| std::vector< const char * > & | ssTagNames, | ||
| std::vector< const char * > & | ssTagValues, | ||
| std::string & | readOpts, | ||
| std::string & | outFile, | ||
| std::string & | writeOpts, | ||
| std::string & | dbgFile, | ||
| bool & | help, | ||
| double & | epsilon | ||
| ) |
Definition at line 322 of file datacoupler_test.cpp.
{
// Initialize some of the outputs to null values indicating not present
// in the argument list.
gNormTag = "";
ssNormTag = "";
readOpts = "PARALLEL=READ_PART;PARTITION=PARALLEL_PARTITION;PARTITION_DISTRIBUTE;PARALLEL_RESOLVE_SHARED_ENTS;PARALLEL_GHOSTS=3.0.1;CPUTIME";
outFile = "";
writeOpts = "PARALLEL=WRITE_PART;CPUTIME";
dbgFile = "";
std::string defaultDbgFile = argv[0]; // The executable name will be the default debug output file.
// These will indicate if we've gotten our required parameters at the end of parsing.
bool haveMeshes = false;
bool haveInterpTag = false;
// Loop over the values in argv pulling out an parsing each one
int npos = 1;
if (argc > 1 && argv[1] == std::string("-h")) {
help = true;
return MB_SUCCESS;
}
while (npos < argc) {
if (argv[npos] == std::string("-meshes")) {
// Parse out the mesh filenames
npos++;
int numFiles = 2;
meshFiles.resize(numFiles);
for (int i = 0; i < numFiles; i++) {
if ((npos < argc) && (!check_for_flag(argv[npos])))
meshFiles[i] = argv[npos++];
else {
std::cerr << " ERROR - missing correct number of mesh filenames" << std::endl;
return MB_FAILURE;
}
}
haveMeshes = true;
}
else if (argv[npos] == std::string("-itag")) {
// Parse out the interpolation tag
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
interpTag = argv[npos++];
else {
std::cerr << " ERROR - missing <interp_tag>" << std::endl;
return MB_FAILURE;
}
haveInterpTag = true;
}
else if (argv[npos] == std::string("-meth")) {
// Parse out the interpolation tag
npos++;
if (argv[npos][0] == '0') method = DataCoupler::CONSTANT;
else if (argv[npos][0] == '1') method = DataCoupler::LINEAR_FE;
else if (argv[npos][0] == '2') method = DataCoupler::QUADRATIC_FE;
else if (argv[npos][0] == '3') method = DataCoupler::SPECTRAL;
else {
std::cerr << " ERROR - unrecognized method number " << method << std::endl;
return MB_FAILURE;
}
npos++;
}
else if (argv[npos] == std::string("-eps")) {
// Parse out the tolerance
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
epsilon = atof(argv[npos++]);
else {
std::cerr << " ERROR - missing <epsilon>" << std::endl;
return MB_FAILURE;
}
}
else if (argv[npos] == std::string("-gnorm")) {
// Parse out the global normalization tag
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
gNormTag = argv[npos++];
else {
std::cerr << " ERROR - missing <gnorm_tag>" << std::endl;
return MB_FAILURE;
}
}
else if (argv[npos] == std::string("-ssnorm")) {
// Parse out the subset normalization tag and selection criteria
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
ssNormTag = argv[npos++];
else {
std::cerr << " ERROR - missing <ssnorm_tag>" << std::endl;
return MB_FAILURE;
}
if ((npos < argc) && (!check_for_flag(argv[npos]))) {
char* opts = argv[npos++];
char sep1[1] = {';'};
char sep2[1] = {'='};
bool end_vals_seen = false;
std::vector<char *> tmpTagOpts;
// first get the options
for (char* i = strtok(opts, sep1); i; i = strtok(0, sep1)) {
tmpTagOpts.push_back(i);
}
// parse out the name and val or just name.
for (unsigned int j = 0; j < tmpTagOpts.size(); j++) {
char* e = strtok(tmpTagOpts[j], sep2);
ssTagNames.push_back(e);
e = strtok(0, sep2);
if (e != NULL) {
// We have a value
if (end_vals_seen) {
// ERROR we should not have a value after none are seen
std::cerr << " ERROR - new value seen after end of values in <ssnorm_selection>" << std::endl;
return MB_FAILURE;
}
// Otherwise get the value string from e and convert it to an int
int *valp = new int;
*valp = atoi(e);
ssTagValues.push_back((const char *) valp);
}
else {
// Otherwise there is no '=' so push a null on the list
end_vals_seen = true;
ssTagValues.push_back((const char *) 0);
}
}
}
else {
std::cerr << " ERROR - missing <ssnorm_selection>" << std::endl;
return MB_FAILURE;
}
}
else if (argv[npos] == std::string("-ropts")) {
// Parse out the mesh file read options
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
readOpts = argv[npos++];
else {
std::cerr << " ERROR - missing <roptions>" << std::endl;
return MB_FAILURE;
}
}
else if (argv[npos] == std::string("-outfile")) {
// Parse out the output file name
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
outFile = argv[npos++];
else {
std::cerr << " ERROR - missing <out_file>" << std::endl;
return MB_FAILURE;
}
}
else if (argv[npos] == std::string("-wopts")) {
// Parse out the output file write options
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
writeOpts = argv[npos++];
else {
std::cerr << " ERROR - missing <woptions>" << std::endl;
return MB_FAILURE;
}
}
else if (argv[npos] == std::string("-dbgout")) {
// Parse out the debug output file name.
// If no name then use the default.
npos++;
if ((npos < argc) && (!check_for_flag(argv[npos])))
dbgFile = argv[npos++];
else
dbgFile = defaultDbgFile;
}
else {
// Unrecognized parameter. Skip it and move along.
std::cerr << " ERROR - Unrecognized parameter:" << argv[npos] << std::endl;
std::cerr << " Skipping..." << std::endl;
npos++;
}
}
if (!haveMeshes) {
meshFiles.resize(2);
meshFiles[0] = std::string(TestDir + "/64bricks_1khex.h5m");
meshFiles[1] = std::string(TestDir + "/64bricks_12ktet.h5m");
std::cout << "Mesh files not entered; using default files "
<< meshFiles[0] << " and " << meshFiles[1] << std::endl;
}
if (!haveInterpTag) {
interpTag = "vertex_field";
std::cout << "Interpolation field name not given, using default of " << interpTag << std::endl;
}
#ifdef HDF5_FILE
if (1 == argc) {
std::cout << "No arguments given; using output file dum.h5m." << std::endl;
outFile = "dum.h5m";
}
#endif
return MB_SUCCESS;
}
| int main | ( | int | argc, |
| char ** | argv | ||
| ) |
Definition at line 109 of file datacoupler_test.cpp.
{
int err;
#ifdef USE_MPI
// need to init MPI first, to tell how many procs and rank
err = MPI_Init(&argc, &argv);
#endif
std::vector<const char *> ssTagNames, ssTagValues;
std::vector<std::string> meshFiles;
std::string interpTag, gNormTag, ssNormTag, readOpts, outFile, writeOpts, dbgFile;
DataCoupler::Method method = DataCoupler::CONSTANT;
ErrorCode result = MB_SUCCESS;
bool help = false;
double toler = 5.e-10;
result = get_file_options(argc, argv, meshFiles, method, interpTag,
gNormTag, ssNormTag, ssTagNames, ssTagValues,
readOpts, outFile, writeOpts, dbgFile, help, toler);
if (result != MB_SUCCESS || help) {
print_usage(argv);
#ifdef USE_MPI
err = MPI_Finalize();
#endif
return 1;
}
int nprocs = 1, rank = 0;
#ifdef USE_MPI
err = MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
err = MPI_Comm_rank(MPI_COMM_WORLD, &rank);
std::vector<ParallelComm *> pcs(meshFiles.size());
#endif
// redirect stdout and stderr if dbgFile is not null
if (!dbgFile.empty()) {
std::stringstream dfname;
dfname << dbgFile << rank << ".txt";
if (!std::freopen(dfname.str().c_str(), "a", stdout)) return false;
if (!std::freopen(dfname.str().c_str(), "a", stderr)) return false;
}
// create MOAB instance based on that
Interface *mbImpl = new Core();
if (NULL == mbImpl) return 1;
iMesh_Instance iMeshInst = reinterpret_cast<iMesh_Instance>( new MBiMesh(mbImpl) );
// read in mesh(es)
// Create root sets for each mesh using the iMesh API. Then pass these
// to the load_file functions to be populated.
iBase_EntitySetHandle *roots = (iBase_EntitySetHandle *) malloc(sizeof(iBase_EntitySetHandle) * meshFiles.size());
for (unsigned int i = 0; i < meshFiles.size(); i++) {
std::string newReadopts;
std::ostringstream extraOpt;
#ifdef USE_MPI
pcs[i] = new ParallelComm(mbImpl, MPI_COMM_WORLD);
int index = pcs[i]->get_id();
extraOpt << ";PARALLEL_COMM=" << index;
newReadopts = readOpts+extraOpt.str();
#endif
iMesh_createEntSet(iMeshInst, 0, &(roots[i]), &err);
result = mbImpl->load_file( meshFiles[i].c_str(), (EntityHandle *)&roots[i], newReadopts.c_str() );
PRINT_LAST_ERROR;
}
#ifdef USE_MPI
result = report_iface_ents(mbImpl, pcs, true);
PRINT_LAST_ERROR;
#endif
double instant_time=0.0, pointloc_time=0.0, interp_time=0.0, gnorm_time=0.0, ssnorm_time=0.0;
// test interpolation and global normalization and subset normalization
result = test_interpolation(mbImpl, method, interpTag, gNormTag, ssNormTag,
ssTagNames, ssTagValues, roots, pcs,
instant_time, pointloc_time, interp_time,
gnorm_time, ssnorm_time, toler);
PRINT_LAST_ERROR;
reduceMax(instant_time);
reduceMax(pointloc_time);
reduceMax(interp_time);
if(0==rank)
printf("\nMax time : %g %g %g (inst loc interp -- %d procs )\n", instant_time,
pointloc_time, interp_time, nprocs);
// output mesh
if (!outFile.empty()) {
Range partSets;
// only save the target mesh
partSets.insert((EntityHandle)roots[1]);
std::string newwriteOpts = writeOpts;
std::ostringstream extraOpt;
#ifdef USE_MPI
extraOpt << ";PARALLEL_COMM=" << 1;
newwriteOpts += extraOpt.str();
#endif
result = mbImpl->write_file(outFile.c_str(), NULL, newwriteOpts.c_str(), partSets);
PRINT_LAST_ERROR;
std::cout << "Wrote " << outFile << std::endl;
std::cout << "mbcoupler_test complete." << std::endl;
}
#ifdef USE_MPI
for (unsigned int i = 0; i < meshFiles.size(); i++) {
delete pcs[i];
}
#endif
delete mbImpl;
//may be leaking iMeshInst, don't care since it's end of program. Remove above deletes?
#ifdef USE_MPI
err = MPI_Finalize();
#endif
return 0;
}
| void print_usage | ( | char ** | argv | ) |
Definition at line 285 of file datacoupler_test.cpp.
{
std::cerr << "Usage: ";
std::cerr << argv[0] << " -meshes <source_mesh> <target_mesh> -itag <interp_tag> [-gnorm <gnorm_tag>] [-ssnorm <ssnorm_tag> <ssnorm_selection>] [-ropts <roptions>] [-outfile <out_file> [-wopts <woptions>]] [-dbgout [<dbg_file>]]" << std::endl;
std::cerr << " -meshes" << std::endl;
std::cerr << " Read in mesh files <source_mesh> and <target_mesh>." << std::endl;
std::cerr << " -itag" << std::endl;
std::cerr << " Interpolate tag <interp_tag> from source mesh to target mesh." << std::endl;
std::cerr << " -gnorm" << std::endl;
std::cerr << " Normalize the value of tag <gnorm_tag> over then entire mesh and save to" << std::endl;
std::cerr << " tag \"<gnorm_tag>_normf\" on the mesh set. Do this for all meshes." << std::endl;
std::cerr << " -ssnorm" << std::endl;
std::cerr << " Normalize the value of tag <ssnorm_tag> over subsets of a mesh and save to" << std::endl;
std::cerr << " tag \"<ssnorm_tag>_normf\" on the Entity Set for each subset. Subsets are selected" << std::endl;
std::cerr << " using criteria in <ssnorm_selection>. Do this for all meshes." << std::endl;
std::cerr << " -ropts" << std::endl;
std::cerr << " Read in the mesh files using options in <roptions>." << std::endl;
std::cerr << " -outfile" << std::endl;
std::cerr << " Write out target mesh to <out_file>." << std::endl;
std::cerr << " -wopts" << std::endl;
std::cerr << " Write out mesh files using options in <woptions>." << std::endl;
std::cerr << " -dbgout" << std::endl;
std::cerr << " Write stdout and stderr streams to the file \'<dbg_file>.txt\'." << std::endl;
std::cerr << " -eps" << std::endl;
std::cerr << " epsilon" << std::endl;
std::cerr << " -meth <method> (0=CONSTANT, 1=LINEAR_FE, 2=QUADRATIC_FE, 3=SPECTRAL)" << std::endl;
}
| void reduceMax | ( | double & | v | ) |
Definition at line 99 of file datacoupler_test.cpp.
{
double buf;
MPI_Barrier(MPI_COMM_WORLD);
MPI_Allreduce( &v, &buf, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
v=buf;
}
| ErrorCode test_interpolation | ( | Interface * | mbImpl, |
| DataCoupler::Method | method, | ||
| std::string & | interpTag, | ||
| std::string & | gNormTag, | ||
| std::string & | ssNormTag, | ||
| std::vector< const char * > & | ssTagNames, | ||
| std::vector< const char * > & | ssTagValues, | ||
| iBase_EntitySetHandle * | roots, | ||
| std::vector< ParallelComm * > & | pcs, | ||
| double & | instant_time, | ||
| double & | pointloc_time, | ||
| double & | interp_time, | ||
| double & | gnorm_time, | ||
| double & | ssnorm_time, | ||
| double & | toler | ||
| ) |
Definition at line 546 of file datacoupler_test.cpp.
{
assert(method >= DataCoupler::CONSTANT && method <= DataCoupler::SPECTRAL);
// source is 1st mesh, target is 2nd
Range src_elems, targ_elems, targ_verts;
ErrorCode result = pcs[0]->get_part_entities(src_elems, 3);
PRINT_LAST_ERROR;
CpuTimer timer;
// instantiate a coupler, which also initializes the tree
DataCoupler dc(mbImpl, src_elems, 0, pcs[0]);
// initialize spectral elements, if they exist
// bool specSou=false, specTar = false;
// result = mbc.initialize_spectral_elements((EntityHandle)roots[0], (EntityHandle)roots[1], specSou, specTar);
instant_time = timer.time_since_birth();
// get points from the target mesh to interpolate
// we have to treat differently the case when the target is a spectral mesh
// in that case, the points of interest are the GL points, not the vertex nodes
std::vector<double> vpos; // this will have the positions we are interested in
int numPointsOfInterest = 0;
#ifdef USE_MPI
result = pcs[1]->get_part_entities(targ_elems, 3);
#else
result = mbImpl->get_entities_by_dimension((EntityHandle)roots[1], 3);
#endif
PRINT_LAST_ERROR;
// first get all vertices adj to partition entities in target mesh
if (DataCoupler::CONSTANT == method)
targ_verts = targ_elems;
else
result = mbImpl->get_adjacencies(targ_elems, 0, false, targ_verts,
Interface::UNION);
PRINT_LAST_ERROR;
#ifdef USE_MPI
// then get non-owned verts and subtract
Range tmp_verts;
result = pcs[1]->get_pstatus_entities(0, PSTATUS_NOT_OWNED, tmp_verts);
PRINT_LAST_ERROR;
targ_verts = subtract( targ_verts, tmp_verts);
#endif
// get position of these entities; these are the target points
numPointsOfInterest = (int)targ_verts.size();
vpos.resize(3*targ_verts.size());
result = mbImpl->get_coords(targ_verts, &vpos[0]);
PRINT_LAST_ERROR;
// locate those points in the source mesh
#ifdef USE_MPI
std::cout << "rank "<< pcs[0]->proc_config().proc_rank();
#endif
std::cout << " points of interest: " << numPointsOfInterest << "\n";
result = dc.locate_points(&vpos[0], numPointsOfInterest, toler);
PRINT_LAST_ERROR;
pointloc_time = timer.time_elapsed();
// now interpolate tag onto target points
std::vector<double> field(numPointsOfInterest);
result = dc.interpolate(method, interpTag, &field[0]);
PRINT_LAST_ERROR;
interp_time = timer.time_elapsed();
/*
// do global normalization if specified
if (!gNormTag.empty()) {
int err;
// Normalize the source mesh
err = dc.normalize_mesh(roots[0], gNormTag.c_str(), DataCoupler::VOLUME, 4);
PRINT_LAST_ERR;
// Normalize the target mesh
err = dc.normalize_mesh(roots[1], gNormTag.c_str(), DataCoupler::VOLUME, 4);
PRINT_LAST_ERR;
}
gnorm_time = timer.time_elapsed();
// do subset normalization if specified
if (!ssNormTag.empty()) {
int err;
err = dc.normalize_subset(roots[0],
ssNormTag.c_str(),
&ssTagNames[0],
ssTagNames.size(),
&ssTagValues[0],
DataCoupler::VOLUME,
4);
PRINT_LAST_ERR;
err = dc.normalize_subset(roots[1],
ssNormTag.c_str(),
&ssTagNames[0],
ssTagNames.size(),
&ssTagValues[0],
DataCoupler::VOLUME,
4);
PRINT_LAST_ERR;
}
ssnorm_time = timer.time_elapsed();
*/
// set field values as tag on target vertices
// use original tag
Tag tag;
result = mbImpl->tag_get_handle(interpTag.c_str(), 1, MB_TYPE_DOUBLE, tag); PRINT_LAST_ERROR;
result = mbImpl->tag_set_data(tag, targ_verts, &field[0]);
PRINT_LAST_ERROR;
// done
return MB_SUCCESS;
}
| std::string TestDir | ( | "." | ) |
1.7.6.1