From a110369f81d7d8a0723893507342f9abc7524195 Mon Sep 17 00:00:00 2001 From: Aiden Woodruff Date: Sat, 16 Dec 2023 00:30:10 -0500 Subject: [PATCH] Recombine capCube,capCyl,capWing into one capVol Add force volume generation and verification (in -g case). Add command line flags and make it a util_exe_func. Add a fourth test field to capStoneSizeFields.h. Signed-off-by: Aiden Woodruff --- capstone_clis/capStoneSizeFields.h | 78 +++++++- test/CMakeLists.txt | 11 +- test/capCube.cc | 166 ----------------- test/capCyl.cc | 165 ----------------- test/capVol.cc | 281 +++++++++++++++++++++++++++++ test/capWing.cc | 166 ----------------- test/testing.cmake | 11 +- 7 files changed, 368 insertions(+), 510 deletions(-) delete mode 100644 test/capCube.cc delete mode 100644 test/capCyl.cc create mode 100644 test/capVol.cc delete mode 100644 test/capWing.cc diff --git a/capstone_clis/capStoneSizeFields.h b/capstone_clis/capStoneSizeFields.h index 66bd48571..eeb5b13f1 100644 --- a/capstone_clis/capStoneSizeFields.h +++ b/capstone_clis/capStoneSizeFields.h @@ -1,4 +1,5 @@ - +#include +#include // analytic size field classes // class UniformAniso : public ma::AnisotropicFunction @@ -104,6 +105,81 @@ class Shock : public ma::AnisotropicFunction ma::Mesh* mesh; }; +// Shock along a cylinder like the one in +// https://www.scorec.rpi.edu/~cwsmith/SCOREC//SimModSuite_2023.1-230428dev/MeshSimAdapt/group__MSAEx__GenAnisoMesh.html +class CylBoundaryLayer : public ma::AnisotropicFunction { + ma::Mesh* m; + enum Axis : int {X = 0, Y = 1, Z = 2} axis; + double radius, length; + +public: + CylBoundaryLayer(ma::Mesh* mesh): m(mesh) { + ma::Vector bmin(HUGE_VAL, HUGE_VAL, HUGE_VAL), + bmax(-HUGE_VAL, -HUGE_VAL, -HUGE_VAL); + ma::Iterator *it = m->begin(0); + ma::Entity *e; + while ((e = m->iterate(it))) { + ma::Vector p = ma::getPosition(m, e); + + if (p.x() < bmin.x()) bmin.x() = p.x(); + if (p.x() > bmax.x()) bmax.x() = p.x(); + if (p.y() < bmin.y()) bmin.y() = p.y(); + if (p.y() > bmax.y()) bmax.y() = p.y(); + if (p.z() < bmin.z()) bmin.z() = p.z(); + if (p.z() > bmax.z()) bmax.z() = p.z(); + } + m->end(it); + + ma::Vector len = bmax - bmin; + if (len.x() > len.y() && len.x() > len.z()) { + lion_oprint(1, "Shock2 along X axis.\n"); + length = len.x(); + axis = Axis::X; + // Radial lengths should be the same, but we find average diameter and + // then halve that. + radius = (len.y() / 2 + len.z() / 2) / 2; + } else if (len.y() > len.x() && len.y() > len.z()) { + lion_oprint(1, "Shock2 along Y axis.\n"); + length = len.y(); + axis = Axis::Y; + radius = (len.x() / 2 + len.z() / 2) / 2; + } else { + lion_oprint(1, "Shock2 along Z axis.\n"); + length = len.z(); + axis = Axis::Z; + radius = (len.x() / 2 + len.y() / 2) / 2; + } + } + + void getValue(ma::Entity *v, ma::Matrix &R, ma::Vector &H) { + const double M = 2.5, aspect_ratio = 16; + ma::Vector p = ma::getPosition(m, v); + H = ma::Vector(1, 1, 1) * radius / M; + R[axis] = ma::Vector(0, 0, 0); + R[axis][axis] = 1; + + int cosDirs[] = {1, 2, 0}, sinDirs[] = {2, 0, 1}; // Rotated sin/cos axes + int cosDir = cosDirs[axis], sinDir = sinDirs[axis]; + + double distFromAxis = std::sqrt(p[cosDir] * p[cosDir] + p[sinDir] * p[sinDir]); + + if (std::fabs(radius * radius - p[cosDir] * p[cosDir] - + p[sinDir] * p[sinDir]) < radius*radius / 4) { + H[cosDir] /= aspect_ratio; + } + + double cosT = p[cosDir] / distFromAxis; + double sinT = p[sinDir] / distFromAxis; + R[cosDir][cosDir] = cosT; R[cosDir][sinDir] = -sinT; R[cosDir][axis] = 0; + R[sinDir][cosDir] = sinT; R[sinDir][sinDir] = cosT; R[sinDir][axis] = 0; + + // Normalize R rows. Skip axis which is already normal. + for (int i = 0; i < 3; ++i) { + if (i != axis) R[i] = R[i].normalize(); + } + } +}; + class Uniform : public ma::IsotropicFunction { public: diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index d3acbc227..95100aa5b 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -210,15 +210,8 @@ if(ENABLE_SIMMETRIX) endif() if(ENABLE_CAPSTONE) - test_exe_func(capCyl capCyl.cc) - target_link_libraries(capCyl apf_cap) - target_include_directories(capCyl PUBLIC "${PROJECT_SOURCE_DIR}/capstone_clis") - test_exe_func(capWing capWing.cc) - target_link_libraries(capWing apf_cap) - target_include_directories(capWing PUBLIC "${PROJECT_SOURCE_DIR}/capstone_clis") - test_exe_func(capCube capCube.cc) - target_link_libraries(capCube apf_cap) - target_include_directories(capCube PUBLIC "${PROJECT_SOURCE_DIR}/capstone_clis") + util_exe_func(capVol capVol.cc) + target_include_directories(capVol PUBLIC "${PROJECT_SOURCE_DIR}/capstone_clis") endif() # send all the newly added utility executable targets diff --git a/test/capCube.cc b/test/capCube.cc deleted file mode 100644 index 87e41da71..000000000 --- a/test/capCube.cc +++ /dev/null @@ -1,166 +0,0 @@ -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - - -#include "CapstoneModule.h" - -using namespace CreateMG; -using namespace CreateMG::Attribution; -using namespace CreateMG::Mesh; -using namespace CreateMG::Geometry; - -#include "capStoneSizeFields.h" - -size_t countElements(apf::Mesh2* m) { - int dim = m->getDimension(); - size_t numElem = 0; - for (int i = 0; i <= dim; ++i) { - numElem += m->count(i); - } - return numElem; -} - -int main(int argc, char** argv) -{ - /* INIT CALLS */ - MPI_Init(&argc, &argv); - PCU_Comm_Init(); - - const char* createFileName = argv[1]; - - gmi_cap_start(); - gmi_register_cap(); - - /* LOAD CAPSTONE MESH */ - // create an instance of the Capstone Module activating CREATE/CREATE/CREATE - // for the Geometry/Mesh/Attribution databases - const std::string gdbName("Geometry Database : SMLIB");// Switch Create with SMLIB for CAD - const std::string mdbName("Mesh Database : Create"); - const std::string adbName("Attribution Database : Create"); - - CapstoneModule cs("test", gdbName.c_str(), mdbName.c_str(), adbName.c_str()); - - GeometryDatabaseInterface *g = cs.get_geometry(); - MeshDatabaseInterface *m = cs.get_mesh(); - AppContext *c = cs.get_context(); - - PCU_ALWAYS_ASSERT(g); - PCU_ALWAYS_ASSERT(m); - PCU_ALWAYS_ASSERT(c); - - v_string filenames; - filenames.push_back(createFileName); - - M_GModel gmodel = cs.load_files(filenames); - - M_MModel mmodel; - // Pick the volume mesh model from associated mesh models to this geom model - std::vector mmodels; - MG_API_CALL(m, get_associated_mesh_models(gmodel, mmodels)); - PCU_ALWAYS_ASSERT(mmodels.size() == 1); - MG_API_CALL(m, set_current_model(mmodels[0])); - - /* SET THE ADJACENCIES */ - MG_API_CALL(m, set_adjacency_state(REGION2FACE| - REGION2EDGE| - REGION2VERTEX| - FACE2EDGE| - FACE2VERTEX)); - MG_API_CALL(m, set_reverse_states()); - MG_API_CALL(m, set_adjacency_scope(TOPO_EDGE, SCOPE_FULL)); - MG_API_CALL(m, set_adjacency_scope(TOPO_FACE, SCOPE_FULL)); - MG_API_CALL(m, compute_adjacency()); - - /* CONVERT THE MESH TO APF::MESH2 */ - apf::Mesh2* apfCapMesh = apf::createMesh(m, g); - - /* ADD A TEST FIELD TO THE MESH TO DEMONSTRATE SOLUTION TRANSFER */ - apf::Field* tf = apf::createFieldOn(apfCapMesh, "test_field", apf::VECTOR); - apf::MeshEntity* ent; - apf::MeshIterator* it = apfCapMesh->begin(0); - while( (ent = apfCapMesh->iterate(it)) ) { - apf::Vector3 p = ma::getPosition(apfCapMesh, ent); - double x = p[0]; - double y = p[1]; - double z = p[2]; - apf::Vector3 s(y, z, 2*x); - apf::setVector(tf, ent, 0, s); - } - apfCapMesh->end(it); - - std::cout << "Element count: " << countElements(apfCapMesh) << std::endl; - std::cout << "Vertex count: " << apfCapMesh->count(0) << std::endl; - - /* SETUP AND CALL ADAPT */ - ma::Input* in; - - // choose a size field here - ma::AnisotropicFunction* sf = new Shock(apfCapMesh); - - // make pumi fields that hold the "frames" and "scales" for anisotropic size fields - // here we are using user-defined size-fields. Usually, "frames" and "scales" come - // from a solution driven error estimation procedure - apf::Field* frameField = apf::createFieldOn(apfCapMesh, "adapt_frames", apf::MATRIX); - apf::Field* scaleField = apf::createFieldOn(apfCapMesh, "adapt_scales", apf::VECTOR); - - apf::MeshEntity* v; - it = apfCapMesh->begin(0); - while( (v = apfCapMesh->iterate(it)) ) { - apf::Vector3 s; - apf::Matrix3x3 f; - sf->getValue(v, f, s); - apf::setVector(scaleField, v, 0, s); - apf::setMatrix(frameField, v, 0, f); - } - apfCapMesh->end(it); - - in = ma::makeAdvanced(ma::configure(apfCapMesh, scaleField, frameField)); - in->shouldSnap = true; - in->shouldTransferParametric = true; - in->shouldFixShape = true; - in->shouldForceAdaptation = true; - if (apfCapMesh->getDimension() == 2) - in->goodQuality = 0.04; // this is mean-ratio squared - else // 3D meshes - in->goodQuality = 0.027; // this is the mean-ratio cubed - in->maximumIterations = 10; - - ma::adaptVerbose(in, false); - - /* apfCapMesh->verify(); */ - - /* PRINT ADAPTED MESH INFO */ - M_MModel mmdl; - m->get_current_model(mmdl); - std::string info; - m->print_info(mmdl, info); - std::cout << "Element count: " << countElements(apfCapMesh) << std::endl; - std::cout << "Vertex count: " << apfCapMesh->count(0) << std::endl; - - // Reference counts based on test meshes on CreateMG v1230 - long numElem = countElements(apfCapMesh), numElemRef = 5528; - PCU_ALWAYS_ASSERT(std::abs(numElem - numElemRef) < 0.05*numElemRef); - - long numVert = apfCapMesh->count(0), numVertRef = 923; - PCU_ALWAYS_ASSERT(std::abs(numVert - numVertRef) < 0.05*numVertRef); - - /* CLEAN UPS */ - delete sf; - apf::destroyMesh(apfCapMesh); - - /* EXIT CALLS */ - gmi_cap_stop(); - PCU_Comm_Free(); - MPI_Finalize(); -} diff --git a/test/capCyl.cc b/test/capCyl.cc deleted file mode 100644 index 334174b27..000000000 --- a/test/capCyl.cc +++ /dev/null @@ -1,165 +0,0 @@ -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - - -#include "CapstoneModule.h" - -using namespace CreateMG; -using namespace CreateMG::Attribution; -using namespace CreateMG::Mesh; -using namespace CreateMG::Geometry; - -#include "capStoneSizeFields.h" - -size_t countElements(apf::Mesh2* m) { - int dim = m->getDimension(); - size_t numElem = 0; - for (int i = 0; i <= dim; ++i) { - numElem += m->count(i); - } - return numElem; -} - -int main(int argc, char** argv) -{ - /* INIT CALLS */ - MPI_Init(&argc, &argv); - PCU_Comm_Init(); - - const char* createFileName = argv[1]; - - gmi_cap_start(); - gmi_register_cap(); - - /* LOAD CAPSTONE MESH */ - // create an instance of the Capstone Module activating CREATE/CREATE/CREATE - // for the Geometry/Mesh/Attribution databases - const std::string gdbName("Geometry Database : SMLIB");// Switch Create with SMLIB for CAD - const std::string mdbName("Mesh Database : Create"); - const std::string adbName("Attribution Database : Create"); - - CapstoneModule cs("test", gdbName.c_str(), mdbName.c_str(), adbName.c_str()); - - GeometryDatabaseInterface *g = cs.get_geometry(); - MeshDatabaseInterface *m = cs.get_mesh(); - AppContext *c = cs.get_context(); - - PCU_ALWAYS_ASSERT(g); - PCU_ALWAYS_ASSERT(m); - PCU_ALWAYS_ASSERT(c); - - v_string filenames; - filenames.push_back(createFileName); - - M_GModel gmodel = cs.load_files(filenames); - - M_MModel mmodel; - // Pick the volume mesh model from associated mesh models to this geom model - std::vector mmodels; - MG_API_CALL(m, get_associated_mesh_models(gmodel, mmodels)); - PCU_ALWAYS_ASSERT(mmodels.size() == 1); - MG_API_CALL(m, set_current_model(mmodels[0])); - - /* SET THE ADJACENCIES */ - MG_API_CALL(m, set_adjacency_state(REGION2FACE| - REGION2EDGE| - REGION2VERTEX| - FACE2EDGE| - FACE2VERTEX)); - MG_API_CALL(m, set_reverse_states()); - MG_API_CALL(m, set_adjacency_scope(TOPO_EDGE, SCOPE_FULL)); - MG_API_CALL(m, set_adjacency_scope(TOPO_FACE, SCOPE_FULL)); - MG_API_CALL(m, compute_adjacency()); - - /* CONVERT THE MESH TO APF::MESH2 */ - apf::Mesh2* apfCapMesh = apf::createMesh(m, g); - - /* ADD A TEST FIELD TO THE MESH TO DEMONSTRATE SOLUTION TRANSFER */ - apf::Field* tf = apf::createFieldOn(apfCapMesh, "test_field", apf::VECTOR); - apf::MeshEntity* ent; - apf::MeshIterator* it = apfCapMesh->begin(0); - while( (ent = apfCapMesh->iterate(it)) ) { - apf::Vector3 p = ma::getPosition(apfCapMesh, ent); - double x = p[0]; - double y = p[1]; - double z = p[2]; - apf::Vector3 s(y, z, 2*x); - apf::setVector(tf, ent, 0, s); - } - apfCapMesh->end(it); - - std::cout << "Element count: " << countElements(apfCapMesh) << std::endl; - std::cout << "Vertex count: " << apfCapMesh->count(0) << std::endl; - - /* SETUP AND CALL ADAPT */ - ma::Input* in; - - // choose a size field here - ma::AnisotropicFunction* sf = new UniformAniso(apfCapMesh); - - // make pumi fields that hold the "frames" and "scales" for anisotropic size fields - // here we are using user-defined size-fields. Usually, "frames" and "scales" come - // from a solution driven error estimation procedure - apf::Field* frameField = apf::createFieldOn(apfCapMesh, "adapt_frames", apf::MATRIX); - apf::Field* scaleField = apf::createFieldOn(apfCapMesh, "adapt_scales", apf::VECTOR); - - apf::MeshEntity* v; - it = apfCapMesh->begin(0); - while( (v = apfCapMesh->iterate(it)) ) { - apf::Vector3 s; - apf::Matrix3x3 f; - sf->getValue(v, f, s); - apf::setVector(scaleField, v, 0, s); - apf::setMatrix(frameField, v, 0, f); - } - apfCapMesh->end(it); - - in = ma::makeAdvanced(ma::configure(apfCapMesh, scaleField, frameField)); - in->shouldSnap = true; - in->shouldTransferParametric = true; - in->shouldFixShape = true; - in->shouldForceAdaptation = true; - if (apfCapMesh->getDimension() == 2) - in->goodQuality = 0.04; // this is mean-ratio squared - else // 3D meshes - in->goodQuality = 0.027; // this is the mean-ratio cubed - in->maximumIterations = 10; - - ma::adaptVerbose(in, false); - - /* apfCapMesh->verify(); */ - - /* PRINT ADAPTED MESH INFO */ - M_MModel mmdl; - m->get_current_model(mmdl); - std::string info; - m->print_info(mmdl, info); - std::cout << "Element count: " << countElements(apfCapMesh) << std::endl; - std::cout << "Vertex count: " << apfCapMesh->count(0) << std::endl; - - long numElem = countElements(apfCapMesh), numElemRef = 5276; - PCU_ALWAYS_ASSERT(std::abs(numElem - numElemRef) < 0.05*numElemRef); - - long numVert = apfCapMesh->count(0), numVertRef = 881; - PCU_ALWAYS_ASSERT(std::abs(numVert - numVertRef) < 0.05*numVertRef); - - /* CLEAN UPS */ - delete sf; - apf::destroyMesh(apfCapMesh); - - /* EXIT CALLS */ - gmi_cap_stop(); - PCU_Comm_Free(); - MPI_Finalize(); -} diff --git a/test/capVol.cc b/test/capVol.cc new file mode 100644 index 000000000..5eaa0b52c --- /dev/null +++ b/test/capVol.cc @@ -0,0 +1,281 @@ +#include +#include + +// Output +#include + +// Parallelism +#include +#include + +// Mesh interfaces +#include +#include + +// Geometry interfaces +#include +#include + +// Mesh adapt +#include + +using namespace CreateMG; +using namespace CreateMG::Attribution; +using namespace CreateMG::Mesh; +using namespace CreateMG::Geometry; + +#include "capStoneSizeFields.h" + +namespace { + +void myExit(int exit_code = EXIT_SUCCESS) { + gmi_cap_stop(); + PCU_Comm_Free(); + MPI_Finalize(); + exit(exit_code); +} + +void writeCre(CapstoneModule& cs, const std::string& filename) { + GeometryDatabaseInterface *gdbi = cs.get_geometry(); + MeshDatabaseInterface *mdbi = cs.get_mesh(); + AppContext *ctx = cs.get_context(); + + // Get the CRE writer. + Writer *creWriter = get_writer(ctx, "Create Native Writer"); + if (!creWriter) { + lion_eprint(1, "FATAL: Could not find the CRE writer.\n"); + myExit(EXIT_FAILURE); + } + + IdMapping idmapping; + std::vector mmodels; + M_GModel gmodel; + M_MModel mmodel; + gdbi->get_current_model(gmodel); + mdbi->get_current_model(mmodel); + mmodels.clear(); + mmodels.push_back(mmodel); + creWriter->write(ctx, gmodel, mmodels, filename.c_str(), idmapping); +} + +void printUsage(char *argv0) { + printf("USAGE: %s [-agwv] \n", argv0); + printf("Flags:\n" + "-a\tEvaluate size-field analytically.\n" + "-g\tForce mesh generation.\n" + "-v\tEnable verbose output.\n" + "-w\tWrite before.vtk, after.vtk, and after.cre.\n" + "SIZE-FIELDS:\n" + "%d, for uniform anisotropic size-field\n" + "%d, for wing-shock size-field\n" + "%d, for cube-shock size-field\n" + "%d, for cylinder boundary-layer size-field\n", 1, 2, 3, 4); +} + +} // namespace + +int main(int argc, char** argv) { + // Initialize parallelism. + MPI_Init(&argc, &argv); + PCU_Comm_Init(); + + // Initialize logging. + lion_set_stdout(stdout); + lion_set_stderr(stderr); + + // Check arguments or print usage. + if (argc < 3) { + if (PCU_Comm_Self() == 0) { + printUsage(argv[0]); + } + myExit(EXIT_FAILURE); + } + + // Parse arguments. + bool volume_flag = false, write_flag = false, analytic_flag = false, + verbose_flag = false; + for (int i = 1; i < argc - 2; ++i) { + bool flagDone = false; + if (*argv[i] == '-') { + for (int j = 1; argv[i][j] != '\0'; ++j) { + switch(argv[i][j]) { + case '-': + flagDone = true; + break; + case 'a': + analytic_flag = true; + break; + case 'g': + volume_flag = true; + break; + case 'v': + verbose_flag = true; + lion_set_verbosity(1); + break; + case 'w': + write_flag = true; + break; + default: + printf("Error: invalid flag.\n"); + printUsage(argv[0]); + myExit(EXIT_FAILURE); + } + } + } + } + + const char* createFileName = argv[argc - 1]; + int mode = atoi(argv[argc - 2]); + + // Initialize GMI. + gmi_cap_start(); + gmi_register_cap(); + // create an instance of the Capstone Module activating SMLIB/CREATE/CREATE + // for the Geometry/Mesh/Attribution databases + const std::string gdbName("Geometry Database : SMLIB"); + const std::string mdbName("Mesh Database : Create"); + const std::string adbName("Attribution Database : Create"); + + CapstoneModule cs("capTest", gdbName.c_str(), mdbName.c_str(), adbName.c_str()); + + GeometryDatabaseInterface *g = cs.get_geometry(); + MeshDatabaseInterface *m = cs.get_mesh(); + AppContext *c = cs.get_context(); + + PCU_ALWAYS_ASSERT(g); + PCU_ALWAYS_ASSERT(m); + PCU_ALWAYS_ASSERT(c); + + // Load Capstone mesh. + v_string filenames; + filenames.push_back(createFileName); + M_GModel gmodel = cs.load_files(filenames); + + if (volume_flag) { + M_MModel mmodel = cs.generate_mesh(); + if (mmodel.is_invalid()) { + lion_eprint(1, "FATAL: Failed to mesh the model.\n"); + myExit(EXIT_FAILURE); + } + MG_API_CALL(m, set_current_model(mmodel)); + } else { + // Use the first existing mesh model. + std::vector mmodels; + MG_API_CALL(m, get_associated_mesh_models(gmodel, mmodels)); + PCU_ALWAYS_ASSERT(mmodels.size() == 1); + MG_API_CALL(m, set_current_model(mmodels[0])); + } + + if (write_flag) { + writeCre(cs, "core_capVol_before.cre"); + } + + // Calculate adjacencies. + MG_API_CALL(m, set_adjacency_state(REGION2FACE|REGION2EDGE|REGION2VERTEX| + FACE2EDGE|FACE2VERTEX)); + MG_API_CALL(m, set_reverse_states()); + MG_API_CALL(m, compute_adjacency()); + + // Make APF adapter over Capstone mesh. + ma::Mesh* apfCapMesh = apf::createMesh(m, g); + + // Choose appropriate size-field. + ma::AnisotropicFunction* sf = nullptr; + switch (mode) { + case 1: + sf = new UniformAniso(apfCapMesh); + break; + case 2: + sf = new WingShock(apfCapMesh, 50); + break; + case 3: + sf = new Shock(apfCapMesh); + break; + case 4: + sf = new CylBoundaryLayer(apfCapMesh); + break; + default: + lion_eprint(1, "FATAL: Invalid size-field.\n"); + myExit(EXIT_FAILURE); + } + + // Make pumi fields for the frames and scales for anisotropic size-fields. + apf::Field* frameField = nullptr; + apf::Field* scaleField = nullptr; + ma::Input *in = nullptr; + if (!analytic_flag || write_flag) { + frameField = apf::createFieldOn(apfCapMesh, "adapt_frames", apf::MATRIX); + scaleField = apf::createFieldOn(apfCapMesh, "adapt_scales", apf::VECTOR); + + ma::Entity *v; + ma::Iterator* it = apfCapMesh->begin(0); + while( (v = apfCapMesh->iterate(it)) ) { + ma::Vector s; + ma::Matrix f; + sf->getValue(v, f, s); + apf::setVector(scaleField, v, 0, s); + apf::setMatrix(frameField, v, 0, f); + } + apfCapMesh->end(it); + + if (write_flag) { + apf::writeVtkFiles("core_capVol_before", apfCapMesh); + } + } + + if (!analytic_flag) { + // Pass the field data. + in = ma::makeAdvanced(ma::configure(apfCapMesh, scaleField, frameField)); + } else { + // Pass the function. + in = ma::makeAdvanced(ma::configure(apfCapMesh, sf)); + } + + in->shouldSnap = true; + in->shouldTransferParametric = true; + in->shouldFixShape = true; + in->shouldForceAdaptation = true; + if (apfCapMesh->getDimension() == 2) + in->goodQuality = 0.04; // this is mean-ratio squared + else // 3D meshes + in->goodQuality = 0.027; // this is the mean-ratio cubed + in->maximumIterations = 10; + + if (verbose_flag) { + // Adapt with verbose logging but without intermediate VTKs. + ma::adaptVerbose(in, false); + } else { + ma::adapt(in); + } + + if (volume_flag) { + // We can't verify surface meshes. + apfCapMesh->verify(); + } + + if (write_flag) { + apf::writeVtkFiles("core_capVol_after", apfCapMesh); + writeCre(cs, "core_capVol_after.cre"); + } + + /* PRINT ADAPTED MESH INFO */ + if (verbose_flag) { + M_MModel mmdl; + m->get_current_model(mmdl); + std::string info; + m->print_info(mmdl, info); + lion_oprint(1, "%s", info.c_str()); + } + + // Clean up. + if (frameField) apf::destroyField(frameField); + if (scaleField) apf::destroyField(scaleField); + apf::destroyMesh(apfCapMesh); + delete sf; + + // Exit calls. + gmi_cap_stop(); + PCU_Comm_Free(); + MPI_Finalize(); +} + diff --git a/test/capWing.cc b/test/capWing.cc deleted file mode 100644 index 5c9d6f00d..000000000 --- a/test/capWing.cc +++ /dev/null @@ -1,166 +0,0 @@ -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - - -#include "CapstoneModule.h" - -using namespace CreateMG; -using namespace CreateMG::Attribution; -using namespace CreateMG::Mesh; -using namespace CreateMG::Geometry; - -#include "capStoneSizeFields.h" - -size_t countElements(apf::Mesh2* m) { - int dim = m->getDimension(); - size_t numElem = 0; - for (int i = 0; i <= dim; ++i) { - numElem += m->count(i); - } - return numElem; -} - -int main(int argc, char** argv) -{ - /* INIT CALLS */ - MPI_Init(&argc, &argv); - PCU_Comm_Init(); - - const char* createFileName = argv[1]; - - gmi_cap_start(); - gmi_register_cap(); - - /* LOAD CAPSTONE MESH */ - // create an instance of the Capstone Module activating CREATE/CREATE/CREATE - // for the Geometry/Mesh/Attribution databases - const std::string gdbName("Geometry Database : SMLIB");// Switch Create with SMLIB for CAD - const std::string mdbName("Mesh Database : Create"); - const std::string adbName("Attribution Database : Create"); - - CapstoneModule cs("test", gdbName.c_str(), mdbName.c_str(), adbName.c_str()); - - GeometryDatabaseInterface *g = cs.get_geometry(); - MeshDatabaseInterface *m = cs.get_mesh(); - AppContext *c = cs.get_context(); - - PCU_ALWAYS_ASSERT(g); - PCU_ALWAYS_ASSERT(m); - PCU_ALWAYS_ASSERT(c); - - v_string filenames; - filenames.push_back(createFileName); - - M_GModel gmodel = cs.load_files(filenames); - - M_MModel mmodel; - // Pick the volume mesh model from associated mesh models to this geom model - std::vector mmodels; - MG_API_CALL(m, get_associated_mesh_models(gmodel, mmodels)); - PCU_ALWAYS_ASSERT(mmodels.size() == 1); - MG_API_CALL(m, set_current_model(mmodels[0])); - - /* SET THE ADJACENCIES */ - MG_API_CALL(m, set_adjacency_state(REGION2FACE| - REGION2EDGE| - REGION2VERTEX| - FACE2EDGE| - FACE2VERTEX)); - MG_API_CALL(m, set_reverse_states()); - MG_API_CALL(m, set_adjacency_scope(TOPO_EDGE, SCOPE_FULL)); - MG_API_CALL(m, set_adjacency_scope(TOPO_FACE, SCOPE_FULL)); - MG_API_CALL(m, compute_adjacency()); - - /* CONVERT THE MESH TO APF::MESH2 */ - apf::Mesh2* apfCapMesh = apf::createMesh(m, g); - - /* ADD A TEST FIELD TO THE MESH TO DEMONSTRATE SOLUTION TRANSFER */ - apf::Field* tf = apf::createFieldOn(apfCapMesh, "test_field", apf::VECTOR); - apf::MeshEntity* ent; - apf::MeshIterator* it = apfCapMesh->begin(0); - while( (ent = apfCapMesh->iterate(it)) ) { - apf::Vector3 p = ma::getPosition(apfCapMesh, ent); - double x = p[0]; - double y = p[1]; - double z = p[2]; - apf::Vector3 s(y, z, 2*x); - apf::setVector(tf, ent, 0, s); - } - apfCapMesh->end(it); - - std::cout << "Element count: " << countElements(apfCapMesh) << std::endl; - std::cout << "Vertex count: " << apfCapMesh->count(0) << std::endl; - - /* SETUP AND CALL ADAPT */ - ma::Input* in; - - // choose a size field here - ma::AnisotropicFunction* sf = new WingShock(apfCapMesh, 50); - - // make pumi fields that hold the "frames" and "scales" for anisotropic size fields - // here we are using user-defined size-fields. Usually, "frames" and "scales" come - // from a solution driven error estimation procedure - apf::Field* frameField = apf::createFieldOn(apfCapMesh, "adapt_frames", apf::MATRIX); - apf::Field* scaleField = apf::createFieldOn(apfCapMesh, "adapt_scales", apf::VECTOR); - - apf::MeshEntity* v; - it = apfCapMesh->begin(0); - while( (v = apfCapMesh->iterate(it)) ) { - apf::Vector3 s; - apf::Matrix3x3 f; - sf->getValue(v, f, s); - apf::setVector(scaleField, v, 0, s); - apf::setMatrix(frameField, v, 0, f); - } - apfCapMesh->end(it); - - in = ma::makeAdvanced(ma::configure(apfCapMesh, scaleField, frameField)); - in->shouldSnap = true; - in->shouldTransferParametric = true; - in->shouldFixShape = true; - in->shouldForceAdaptation = true; - if (apfCapMesh->getDimension() == 2) - in->goodQuality = 0.04; // this is mean-ratio squared - else // 3D meshes - in->goodQuality = 0.027; // this is the mean-ratio cubed - in->maximumIterations = 10; - - ma::adaptVerbose(in, false); - - /* apfCapMesh->verify(); */ - - /* PRINT ADAPTED MESH INFO */ - M_MModel mmdl; - m->get_current_model(mmdl); - std::string info; - m->print_info(mmdl, info); - std::cout << "Element count: " << countElements(apfCapMesh) << std::endl; - std::cout << "Vertex count: " << apfCapMesh->count(0) << std::endl; - - // Reference counts based on test meshes on CreateMG v1230 - long numElem = countElements(apfCapMesh), numElemRef = 13016; - PCU_ALWAYS_ASSERT(std::abs(numElem - numElemRef) < 0.05*numElemRef); - - long numVert = apfCapMesh->count(0), numVertRef = 2171; - PCU_ALWAYS_ASSERT(std::abs(numVert - numVertRef) < 0.05*numVertRef); - - /* CLEAN UPS */ - delete sf; - apf::destroyMesh(apfCapMesh); - - /* EXIT CALLS */ - gmi_cap_stop(); - PCU_Comm_Free(); - MPI_Finalize(); -} diff --git a/test/testing.cmake b/test/testing.cmake index 6aa29b8e0..41469508e 100644 --- a/test/testing.cmake +++ b/test/testing.cmake @@ -873,7 +873,12 @@ if (PCU_COMPRESS) endif() if(ENABLE_CAPSTONE) - mpi_test(capCyl 1 ./capCyl ${MESHES}/cap/cyl_surf_only.cre) - mpi_test(capWing 1 ./capWing ${MESHES}/cap/wing_surf_only.cre) - mpi_test(capCube 1 ./capCube ${MESHES}/cap/cube_surf_only.cre) + mpi_test(capCyl 1 ./capVol -v 1 ${MESHES}/cap/cyl_surf_only.cre) + mpi_test(capWing 1 ./capVol -v 2 ${MESHES}/cap/wing_surf_only.cre) + mpi_test(capCube 1 ./capVol -v 3 ${MESHES}/cap/cube_surf_only.cre) + mpi_test(capCyl2 1 ./capVol -v 4 ${MESHES}/cap/cyl_surf_only.cre) + mpi_test(capVolCyl 1 ./capVol -vg 1 ${MESHES}/cap/cyl_surf_only.cre) + mpi_test(capVolWing 1 ./capVol -vg 2 ${MESHES}/cap/wing_surf_only.cre) + mpi_test(capVolCube 1 ./capVol -vg 3 ${MESHES}/cap/cube_surf_only.cre) + mpi_test(capVolCyl2 1 ./capVol -vg 4 ${MESHES}/cap/cyl_surf_only.cre) endif()