From 77869347ddaf4250e01afea70a3b68a08c61cd63 Mon Sep 17 00:00:00 2001 From: Iason Date: Fri, 5 Feb 2021 19:58:15 +0200 Subject: [PATCH] Refactoring --- beamformfinder.cpp | 55 ++++++++++++++++++++------------------------ beamformfinder.hpp | 1 - elementalmesh.cpp | 11 +++++++++ flatpattern.cpp | 6 ++--- flatpattern.hpp | 2 +- simulationresult.hpp | 52 +++++++++++++++++++++++++++++------------ utilities.hpp | 12 ++++++++++ 7 files changed, 89 insertions(+), 50 deletions(-) diff --git a/beamformfinder.cpp b/beamformfinder.cpp index 8a530f2..b7bb397 100644 --- a/beamformfinder.cpp +++ b/beamformfinder.cpp @@ -1814,7 +1814,8 @@ FormFinder::executeSimulation(const std::shared_ptr &pJob, shouldApplyInitialDistortion = true; } updateNodalExternalForces(pJob->nodalExternalForces, constrainedVertices); - while (maxDRMIterations == 0 || mCurrentSimulationStep < maxDRMIterations) { + while (settings.maxDRMIterations == 0 || + mCurrentSimulationStep < settings.maxDRMIterations) { // while (true) { updateNormalDerivatives(); updateT1Derivatives(); @@ -1850,33 +1851,25 @@ FormFinder::executeSimulation(const std::shared_ptr &pJob, } if (std::isnan(pMesh->currentTotalKineticEnergy)) { - if (pMesh->elements[0].dimensions.b / pMesh->elements[0].dimensions.h < - 10) { - std::time_t now = time(0); - std::tm *ltm = std::localtime(&now); - std::string dir = "../ProblematicSimulationJobs/" + - std::to_string(ltm->tm_mday) + "_" + - std::to_string(1 + ltm->tm_mon) + "_" + - std::to_string(1900 + ltm->tm_year); - const std::string subDir = - std::filesystem::path(dir) - .append(std::to_string(ltm->tm_hour) + ":" + - std::to_string(ltm->tm_min)) - .string(); - std::filesystem::create_directories(subDir); - pJob->save(subDir); - // std::cout - // << "Non terminating simulation found. Saved simulation job - // to : " - // << dir << std::endl; - // std::terminate(); - } - // std::cout << "Exiting.." << std::endl; - // FormFinder debug; - // FormFinder::Settings settings; - // settings.shouldDraw = true; - // settings.beVerbose = true; - // debug.executeSimulation(pJob, settings); + std::time_t now = time(0); + std::tm *ltm = std::localtime(&now); + std::string dir = "../ProblematicSimulationJobs/" + + std::to_string(ltm->tm_mday) + "_" + + std::to_string(1 + ltm->tm_mon) + "_" + + std::to_string(1900 + ltm->tm_year); + const std::string subDir = std::filesystem::path(dir) + .append(std::to_string(ltm->tm_hour) + + ":" + std::to_string(ltm->tm_min)) + .string(); + std::filesystem::create_directories(subDir); + pJob->save(subDir); + std::cout << "Infinite kinetic energy detected.Exiting.." << std::endl; + FormFinder debug; + FormFinder::Settings settings; + settings.shouldDraw = true; + settings.beVerbose = true; + debug.executeSimulation(pJob, settings); + std::terminate(); break; } @@ -1966,10 +1959,12 @@ mesh->currentTotalPotentialEnergykN*/ SimulationResults results = computeResults(pJob); - if (mCurrentSimulationStep == maxDRMIterations) { + if (mCurrentSimulationStep == settings.maxDRMIterations && + mCurrentSimulationStep != 0) { std::cout << "Did not reach equilibrium before reaching the maximum number " "of DRM steps (" - << maxDRMIterations << "). Breaking simulation" << std::endl; + << settings.maxDRMIterations << "). Breaking simulation" + << std::endl; results.converged = false; } else if (std::isnan(pMesh->currentTotalKineticEnergy)) { results.converged = false; diff --git a/beamformfinder.hpp b/beamformfinder.hpp index 952f7b4..f1b0b7c 100644 --- a/beamformfinder.hpp +++ b/beamformfinder.hpp @@ -204,7 +204,6 @@ public: SimulationResults executeSimulation(const std::shared_ptr &pJob, const Settings &settings = Settings()); - inline static const size_t maxDRMIterations{0}; static void runUnitTests(); void setTotalResidualForcesNormThreshold(double value); diff --git a/elementalmesh.cpp b/elementalmesh.cpp index 0bcd38f..edce87b 100644 --- a/elementalmesh.cpp +++ b/elementalmesh.cpp @@ -326,6 +326,17 @@ bool SimulationMesh::loadPly(const string &plyFilename) { elements[ei].updateRigidity(); } + bool normalsAreAbsent = vert[0].cN().Norm() < 0.000001; + if (normalsAreAbsent) { + CoordType normalVector(0, 0, 1); + std::cout << "Warning: Normals are missing from " << plyFilename + << ". Added normal vector:" << toString(normalVector) + << std::endl; + for (auto &v : vert) { + v.N() = normalVector; + } + } + return true; } diff --git a/flatpattern.cpp b/flatpattern.cpp index a942cf3..309e4e6 100644 --- a/flatpattern.cpp +++ b/flatpattern.cpp @@ -93,9 +93,9 @@ bool FlatPattern::createHoneycombAtom() { return true; } -void FlatPattern::copy(FlatPattern &pattern) { - VCGEdgeMesh::copy(pattern); - baseTriangleHeight = pattern.getBaseTriangleHeight(); +void FlatPattern::copy(FlatPattern ©From) { + VCGEdgeMesh::copy(copyFrom); + baseTriangleHeight = copyFrom.getBaseTriangleHeight(); } void FlatPattern::deleteDanglingEdges() { diff --git a/flatpattern.hpp b/flatpattern.hpp index 2541c4c..4b6a360 100644 --- a/flatpattern.hpp +++ b/flatpattern.hpp @@ -11,7 +11,7 @@ public: const std::vector &edges); bool createHoneycombAtom(); - void copy(FlatPattern &pattern); + void copy(FlatPattern ©From); void tilePattern(VCGEdgeMesh &pattern, VCGTriMesh &tileInto); void tilePattern(VCGEdgeMesh &pattern, VCGPolyMesh &tileInto, diff --git a/simulationresult.hpp b/simulationresult.hpp index 8c3d72a..ae43a04 100644 --- a/simulationresult.hpp +++ b/simulationresult.hpp @@ -54,9 +54,15 @@ template <> struct adl_serializer> { class SimulationJob { // const std::unordered_map nodalForcedNormals; // json labels - inline static std::string jsonLabel_meshFilename{"mesh filename"}; - inline static std::string jsonLabel_constrainedVertices{"fixed vertices"}; - inline static std::string jsonLabel_nodalForces{"forces"}; + struct JSONLabels { + inline static std::string meshFilename{"mesh filename"}; + inline static std::string constrainedVertices{"fixed vertices"}; + inline static std::string nodalForces{"forces"}; + inline static std::string label{"label"}; + inline static std::string meshLabel{ + "label"}; // TODO: should be in the savePly function of the simulation + // mesh class + } jsonLabels; public: std::shared_ptr pMesh; @@ -78,14 +84,14 @@ public: std::string toString() const { nlohmann::json json; if (!constrainedVertices.empty()) { - json[jsonLabel_constrainedVertices] = constrainedVertices; + json[jsonLabels.constrainedVertices] = constrainedVertices; } if (!nodalExternalForces.empty()) { std::unordered_map> arrForces; for (const auto &f : nodalExternalForces) { arrForces[f.first] = f.second; } - json[jsonLabel_nodalForces] = arrForces; + json[jsonLabels.nodalForces] = arrForces; } return json.dump(); @@ -137,6 +143,7 @@ public: } }); + auto cn = polyscope::getCurveNetwork(meshLabel); if (!nodeColors.empty()) { polyscope::getCurveNetwork(meshLabel) ->addNodeColorQuantity("Boundary conditions_" + label, nodeColors) @@ -180,23 +187,23 @@ public: json << ifs; pMesh = std::make_shared(); - if (json.contains(jsonLabel_meshFilename)) { - pMesh->loadPly(json[jsonLabel_meshFilename]); + if (json.contains(jsonLabels.meshFilename)) { + pMesh->loadPly(json[jsonLabels.meshFilename]); } - if (json.contains(jsonLabel_constrainedVertices)) { + if (json.contains(jsonLabels.constrainedVertices)) { constrainedVertices = // auto conV = std::unordered_map>( - json[jsonLabel_constrainedVertices]); + json[jsonLabels.constrainedVertices]); std::cout << "Loaded constrained vertices. Number of constrained " "vertices found:" << constrainedVertices.size() << std::endl; } - if (json.contains(jsonLabel_nodalForces)) { + if (json.contains(jsonLabels.nodalForces)) { auto f = std::unordered_map>( - json[jsonLabel_nodalForces]); + json[jsonLabels.nodalForces]); for (const auto &forces : f) { nodalExternalForces[forces.first] = Vector6d(forces.second); } @@ -204,6 +211,14 @@ public: << nodalExternalForces.size() << std::endl; } + if (json.contains(jsonLabels.label)) { + label = json[jsonLabels.label]; + } + + if (json.contains(jsonLabels.meshLabel)) { + pMesh->setLabel(json[jsonLabels.meshLabel]); + } + return true; } @@ -224,21 +239,28 @@ public: .append(pMesh->getLabel() + ".ply"); returnValue = pMesh->savePly(meshFilename); nlohmann::json json; - json[jsonLabel_meshFilename] = meshFilename; + json[jsonLabels.meshFilename] = meshFilename; if (!constrainedVertices.empty()) { - json[jsonLabel_constrainedVertices] = constrainedVertices; + json[jsonLabels.constrainedVertices] = constrainedVertices; } if (!nodalExternalForces.empty()) { std::unordered_map> arrForces; for (const auto &f : nodalExternalForces) { arrForces[f.first] = f.second; } - json[jsonLabel_nodalForces] = arrForces; + json[jsonLabels.nodalForces] = arrForces; + } + + if (!label.empty()) { + json[jsonLabels.label] = label; + } + if (!pMesh->getLabel().empty()) { + json[jsonLabels.meshLabel] = pMesh->getLabel(); } std::string jsonFilename( std::filesystem::path(pathFolderDirectory) - .append(pMesh->getLabel() + "_simScenario.json")); + .append(label + "_" + pMesh->getLabel() + "_simulationJob.json")); std::ofstream jsonFile(jsonFilename); jsonFile << json; std::cout << "Saved simulation job as:" << jsonFilename << std::endl; diff --git a/utilities.hpp b/utilities.hpp index 92494a8..7918a0d 100644 --- a/utilities.hpp +++ b/utilities.hpp @@ -162,6 +162,13 @@ inline Eigen::MatrixXd toEigenMatrix(const std::vector &v) { #include "polyscope/curve_network.h" #include "polyscope/polyscope.h" +inline void deinitPolyscope() { + if (!polyscope::state::initialized) { + return; + } + + polyscope::shutdown(); +} inline void initPolyscope() { if (polyscope::state::initialized) { return; @@ -207,4 +214,9 @@ void constructInverseMap(const T1 &map, T2 &oppositeMap) { } } +template std::string toString(const T &v) { + return "(" + std::to_string(v[0]) + "," + std::to_string(v[1]) + "," + + std::to_string(v[2]) + ")"; +} + #endif // UTILITIES_H