From d7d1951be66ee3fd1ca381c78174f153b8139160 Mon Sep 17 00:00:00 2001 From: iasonmanolas Date: Tue, 30 Mar 2021 19:30:49 +0300 Subject: [PATCH] Refactoring --- src/main.cpp | 4 +- src/reducedmodeloptimizer.cpp | 171 ++++++++++++++++++++-------------- src/reducedmodeloptimizer.hpp | 87 ++++++++++------- 3 files changed, 157 insertions(+), 105 deletions(-) diff --git a/src/main.cpp b/src/main.cpp index f1ddb01..b20fca1 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -65,7 +65,7 @@ int main(int argc, char *argv[]) { const std::vector numberOfNodesPerSlot{1, 0, 0, 2, 1, 2, 1}; assert(interfaceNodeIndex==numberOfNodesPerSlot[0]+numberOfNodesPerSlot[3]); ReducedModelOptimizer optimizer(numberOfNodesPerSlot); - optimizer.initializePatterns(fullPattern, reducedPattern, {},settings_optimization.xRanges.size()); + optimizer.initializePatterns(fullPattern, reducedPattern, settings_optimization.xRanges.size()); ReducedModelOptimization::Results optimizationResults = optimizer.optimize(settings_optimization); @@ -106,7 +106,7 @@ int main(int argc, char *argv[]) { csv_results << endrow; #ifdef POLYSCOPE_DEFINED - optimizationResults.draw(); + optimizationResults.draw(); #endif return 0; } diff --git a/src/reducedmodeloptimizer.cpp b/src/reducedmodeloptimizer.cpp index 7ab6843..48114e4 100644 --- a/src/reducedmodeloptimizer.cpp +++ b/src/reducedmodeloptimizer.cpp @@ -39,14 +39,14 @@ struct GlobalOptimizationVariables { std::vector reducedPatternMaximumDisplacementSimulationJobs; double ReducedModelOptimizer::computeDisplacementError( - const std::vector &reducedPatternDisplacements, const std::vector &fullPatternDisplacements, + const std::vector &reducedPatternDisplacements, const std::unordered_map &reducedToFullInterfaceViMap, const double &normalizationFactor) { - const double rawError = computeRawDisplacementError(reducedPatternDisplacements, - fullPatternDisplacements, + const double rawError = computeRawDisplacementError(fullPatternDisplacements, + reducedPatternDisplacements, reducedToFullInterfaceViMap); if (global.optimizationSettings.normalizationStrategy != Settings::NormalizationStrategy::NonNormalized) { @@ -56,8 +56,8 @@ double ReducedModelOptimizer::computeDisplacementError( } double ReducedModelOptimizer::computeRawDisplacementError( - const std::vector &reducedPatternDisplacements, const std::vector &fullPatternDisplacements, + const std::vector &reducedPatternDisplacements, const std::unordered_map &reducedToFullInterfaceViMap) { @@ -81,6 +81,58 @@ double ReducedModelOptimizer::computeRawDisplacementError( return error; } + +double ReducedModelOptimizer::computeRawRotationalError( + const std::vector> &rotatedQuaternion_fullPattern, + const std::vector> &rotatedQuaternion_reducedPattern, + const std::unordered_map + &reducedToFullInterfaceViMap) +{ + double rawRotationalError = 0; + for (const auto &reducedToFullInterfaceViPair : reducedToFullInterfaceViMap) { + const double vertexRotationalError + = rotatedQuaternion_fullPattern[reducedToFullInterfaceViPair.second].angularDistance( + rotatedQuaternion_reducedPattern[reducedToFullInterfaceViPair.first]); + rawRotationalError += vertexRotationalError; + } + + return rawRotationalError; +} + +double ReducedModelOptimizer::computeRotationalError( + const std::vector> &rotatedQuaternion_fullPattern, + const std::vector> &rotatedQuaternion_reducedPattern, + const std::unordered_map + &reducedToFullInterfaceViMap, + const double &normalizationFactor) +{ + const double rawRotationalError = computeRawRotationalError(rotatedQuaternion_fullPattern, + rotatedQuaternion_reducedPattern, + reducedToFullInterfaceViMap); + return rawRotationalError / normalizationFactor; +} + +double ReducedModelOptimizer::computeError( + const SimulationResults &simulationResults_fullPattern, + const SimulationResults &simulationResults_reducedPattern, + const std::unordered_map + &reducedToFullInterfaceViMap, + const double &normalizationFactor_translationalDisplacement, + const double &normalizationFactor_rotationalDisplacement) +{ + const double translationalError + = computeDisplacementError(simulationResults_fullPattern.displacements, + simulationResults_reducedPattern.displacements, + reducedToFullInterfaceViMap, + normalizationFactor_translationalDisplacement); + const double rotationalError + = computeRotationalError(simulationResults_fullPattern.rotationalDisplacementQuaternion, + simulationResults_reducedPattern.rotationalDisplacementQuaternion, + reducedToFullInterfaceViMap, + normalizationFactor_rotationalDisplacement); + return translationalError + rotationalError; +} + double ReducedModelOptimizer::objective(double E,double A,double J,double I2,double I3, double innerHexagonSize, double innerHexagonRotationAngle) { @@ -129,24 +181,12 @@ double ReducedModelOptimizer::objective(long n, const double *x) { reducedModelResults.internalPotentialEnergy - global.fullPatternResults[simulationScenarioIndex].internalPotentialEnergy); } else { - const double translationalError = computeDisplacementError( - reducedModelResults.displacements, - global.fullPatternResults[simulationScenarioIndex].displacements, + simulationScenarioError += computeError( + global.fullPatternResults[simulationScenarioIndex], + reducedModelResults, global.reducedToFullInterfaceViMap, - global.translationalDisplacementNormalizationValues[simulationScenarioIndex]); - double rotationalError = 0; - for (const auto &reducedToFullInterfaceViPair : global.reducedToFullInterfaceViMap) { - const double vertexRotationalError - = global.fullPatternResults[simulationScenarioIndex] - .rotationalDisplacementQuaternion[reducedToFullInterfaceViPair.second] - .angularDistance(reducedModelResults.rotationalDisplacementQuaternion - [reducedToFullInterfaceViPair.first]); - rotationalError += vertexRotationalError; - } - simulationScenarioError += translationalError - + rotationalError - / global.rotationalDisplacementNormalizationValues - [simulationScenarioIndex]; + global.translationalDisplacementNormalizationValues[simulationScenarioIndex], + global.rotationalDisplacementNormalizationValues[simulationScenarioIndex]); } //#ifdef POLYSCOPE_DEFINED // std::cout << "sim error:" << simulationScenarioError << std::endl; @@ -246,7 +286,6 @@ void ReducedModelOptimizer::createSimulationMeshes( } void ReducedModelOptimizer::computeMaps( - const std::unordered_set &reducedModelExcludedEdges, const std::unordered_map> &slotToNode, PatternGeometry &fullPattern, PatternGeometry &reducedPattern, std::unordered_map @@ -375,13 +414,15 @@ void ReducedModelOptimizer::computeMaps( } } -void ReducedModelOptimizer::computeMaps( - PatternGeometry &fullPattern, PatternGeometry &reducedPattern, - const std::unordered_set &reducedModelExcludedEdges) { - ReducedModelOptimizer::computeMaps( - reducedModelExcludedEdges, slotToNode, fullPattern, reducedPattern, - global.reducedToFullInterfaceViMap, m_fullToReducedInterfaceViMap, - m_fullPatternOppositeInterfaceViMap); +void ReducedModelOptimizer::computeMaps(PatternGeometry &fullPattern, + PatternGeometry &reducedPattern) +{ + ReducedModelOptimizer::computeMaps(slotToNode, + fullPattern, + reducedPattern, + global.reducedToFullInterfaceViMap, + m_fullToReducedInterfaceViMap, + m_fullPatternOppositeInterfaceViMap); } ReducedModelOptimizer::ReducedModelOptimizer( @@ -390,11 +431,10 @@ ReducedModelOptimizer::ReducedModelOptimizer( FlatPatternTopology::constructSlotToNodeMap(nodeToSlot, slotToNode); } -void ReducedModelOptimizer::initializePatterns( - PatternGeometry &fullPattern, PatternGeometry &reducedPattern, - const std::unordered_set &reducedModelExcludedEdges,const int& optimizationParameters) { - // fullPattern.setLabel("full_pattern_" + fullPattern.getLabel()); - // reducedPattern.setLabel("reduced_pattern_" + reducedPattern.getLabel()); +void ReducedModelOptimizer::initializePatterns(PatternGeometry &fullPattern, + PatternGeometry &reducedPattern, + const int &optimizationParameters) +{ assert(fullPattern.VN() == reducedPattern.VN() && fullPattern.EN() >= reducedPattern.EN()); #if POLYSCOPE_DEFINED @@ -406,31 +446,8 @@ void ReducedModelOptimizer::initializePatterns( copyFullPattern.copy(fullPattern); copyReducedPattern.copy(reducedPattern); global.baseTriangle = copyReducedPattern.getBaseTriangle(); - /* - * Here we create the vector that connects the central node with the bottom - * left node of the base triangle. During the optimization the vi%2==0 nodes - * move on these vectors. - * */ - // const double h = copyReducedPattern.getBaseTriangleHeight(); - // double baseTriangle_bottomEdgeSize = 2 * h / tan(vcg::math::ToRad(60.0)); - // VectorType baseTriangle_leftBottomNode(-baseTriangle_bottomEdgeSize / 2, -h, - // 0); - // for (int rotationCounter = 0; rotationCounter < fanSize; rotationCounter++) { - // VectorType rotatedVector = - // vcg::RotationMatrix(patternPlaneNormal, - // vcg::math::ToRad(rotationCounter * 60.0)) * - // baseTriangle_leftBottomNode; - // global.g_innerHexagonVectors[rotationCounter] = rotatedVector; - // } - // const double innerHexagonInitialPos_x = - // copyReducedPattern.vert[0].cP()[0] / global.g_innerHexagonVectors[0][0]; - // const double innerHexagonInitialPos_y = - // copyReducedPattern.vert[0].cP()[1] / global.g_innerHexagonVectors[0][1]; - // global.innerHexagonInitialPos = innerHexagonInitialPos_x; - // global.innerHexagonInitialRotationAngle = - // 30; /* NOTE: Here I assume that the CW reduced pattern is given as input. - // This is not very generic */ - computeMaps(copyFullPattern, copyReducedPattern, reducedModelExcludedEdges); + + computeMaps(copyFullPattern, copyReducedPattern); createSimulationMeshes(copyFullPattern, copyReducedPattern); initializeOptimizationParameters(m_pReducedPatternSimulationMesh,optimizationParameters); } @@ -716,17 +733,31 @@ ReducedModelOptimizer::runOptimization(const Settings &settings) { for (int simulationScenarioIndex:global.simulationScenarioIndices) { SimulationResults reducedModelResults = simulator.executeSimulation( global.reducedPatternSimulationJobs[simulationScenarioIndex]); - const double error = computeDisplacementError( - reducedModelResults.displacements, + + const double rawTranslationalError = computeRawDisplacementError( global.fullPatternResults[simulationScenarioIndex].displacements, - global.reducedToFullInterfaceViMap, - global.translationalDisplacementNormalizationValues[simulationScenarioIndex]); - results.rawObjectiveValue += computeRawDisplacementError( reducedModelResults.displacements, - global.fullPatternResults[simulationScenarioIndex].displacements, global.reducedToFullInterfaceViMap); - results.objectiveValuePerSimulationScenario[simulationScenarioIndex] = - error; + const double rawRotationalError = computeRawRotationalError( + global.fullPatternResults[simulationScenarioIndex].rotationalDisplacementQuaternion, + reducedModelResults.rotationalDisplacementQuaternion, + global.reducedToFullInterfaceViMap); + + results.rawObjectiveValue += rawTranslationalError + rawRotationalError; + + const double normalizedTranslationalError + = rawTranslationalError + / global.translationalDisplacementNormalizationValues[simulationScenarioIndex]; + const double normalizedRotationalError + = rawRotationalError + / global.rotationalDisplacementNormalizationValues[simulationScenarioIndex]; + std::cout << "Simulation scenario:" << simulationScenarioStrings[simulationScenarioIndex] + << std::endl; + std::cout << "Translational error:" << normalizedTranslationalError << std::endl; + std::cout << "Rotational error:" << normalizedRotationalError << std::endl; + std::cout << std::endl; + results.objectiveValuePerSimulationScenario[simulationScenarioIndex] + = normalizedTranslationalError + normalizedRotationalError; } results.time = elapsed.count() / 1000.0; const bool printDebugInfo = false; @@ -880,8 +911,6 @@ ReducedModelOptimizer::createScenarios( } void ReducedModelOptimizer::computeObjectiveValueNormalizationFactors() { - if (global.optimizationSettings.normalizationStrategy == - Settings::NormalizationStrategy::Epsilon) { // Compute the sum of the displacement norms std::vector fullPatternTranslationalDisplacementNormSum(NumberOfSimulationScenarios); @@ -939,10 +968,12 @@ void ReducedModelOptimizer::computeObjectiveValueNormalizationFactors() { global.rotationalDisplacementNormalizationValues[simulationScenarioIndex] = std::max(fullPatternAngularDistance[simulationScenarioIndex], epsilon_rotationalDisplacement); + } else { + global.translationalDisplacementNormalizationValues[simulationScenarioIndex] = 1; + global.rotationalDisplacementNormalizationValues[simulationScenarioIndex] = 1; } } } -} Results ReducedModelOptimizer::optimize( const Settings &optimizationSettings, diff --git a/src/reducedmodeloptimizer.hpp b/src/reducedmodeloptimizer.hpp index 4ea8f4e..ca4561d 100644 --- a/src/reducedmodeloptimizer.hpp +++ b/src/reducedmodeloptimizer.hpp @@ -47,8 +47,9 @@ public: SimulationJob getReducedSimulationJob(const SimulationJob &fullModelSimulationJob); - void initializePatterns(PatternGeometry &fullPattern, PatternGeometry &reducedPatterm, - const std::unordered_set &reducedModelExcludedEges, const int &optimizationParameters); + void initializePatterns(PatternGeometry &fullPattern, + PatternGeometry &reducedPatterm, + const int &optimizationParameters); static void runSimulation(const std::string &filename, std::vector &x); @@ -66,16 +67,15 @@ public: PatternGeometry &fullModel, PatternGeometry &reducedModel, std::shared_ptr &pFullPatternSimulationMesh, std::shared_ptr &pReducedPatternSimulationMesh); - static void computeMaps( - const std::unordered_set &reducedModelExcludedEdges, - const std::unordered_map> &slotToNode, - PatternGeometry &fullPattern, PatternGeometry &reducedPattern, - std::unordered_map - &reducedToFullInterfaceViMap, - std::unordered_map - &fullToReducedInterfaceViMap, - std::unordered_map - &fullPatternOppositeInterfaceViMap); + static void computeMaps(const std::unordered_map> &slotToNode, + PatternGeometry &fullPattern, + PatternGeometry &reducedPattern, + std::unordered_map + &reducedToFullInterfaceViMap, + std::unordered_map + &fullToReducedInterfaceViMap, + std::unordered_map + &fullPatternOppositeInterfaceViMap); static void visualizeResults( const std::vector> &fullPatternSimulationJobs, const std::vector> &reducedPatternSimulationJobs, @@ -90,38 +90,59 @@ public: &reducedToFullInterfaceViMap); static double computeRawDisplacementError( - const std::vector &reducedPatternDisplacements, const std::vector &fullPatternDisplacements, + const std::vector &reducedPatternDisplacements, const std::unordered_map &reducedToFullInterfaceViMap); static double computeDisplacementError( - const std::vector &reducedPatternDisplacements, const std::vector &fullPatternDisplacements, + const std::vector &reducedPatternDisplacements, const std::unordered_map &reducedToFullInterfaceViMap, const double &normalizationFactor); - static double objective(double E, double A, double innerHexagonSize, double innerHexagonRotationAngle); - private: - static void computeDesiredReducedModelDisplacements( - const SimulationResults &fullModelResults, - const std::unordered_map &displacementsReducedToFullMap, - Eigen::MatrixX3d &optimalDisplacementsOfReducedModel); - static ReducedModelOptimization::Results runOptimization( - const ReducedModelOptimization::Settings &settings); - std::vector> createScenarios( - const std::shared_ptr &pMesh); - void computeMaps(PatternGeometry &fullModel, - PatternGeometry &reducedPattern, - const std::unordered_set &reducedModelExcludedEges); - void createSimulationMeshes(PatternGeometry &fullModel, PatternGeometry &reducedModel); - static void initializeOptimizationParameters(const std::shared_ptr &mesh, - const int &optimizationParamters); + static double objective(double E, + double A, + double innerHexagonSize, + double innerHexagonRotationAngle); + static double computeRawRotationalError( + const std::vector> &rotatedQuaternion_fullPattern, + const std::vector> &rotatedQuaternion_reducedPattern, + const std::unordered_map + &reducedToFullInterfaceViMap); - static double objective(long n, const double *x); - DRMSimulationModel simulator; - void computeObjectiveValueNormalizationFactors(); + static double computeRotationalError( + const std::vector> &rotatedQuaternion_fullPattern, + const std::vector> &rotatedQuaternion_reducedPattern, + const std::unordered_map + &reducedToFullInterfaceViMap, + const double &normalizationFactor); + static double computeError( + const SimulationResults &simulationResults_fullPattern, + const SimulationResults &simulationResults_reducedPattern, + const std::unordered_map + &reducedToFullInterfaceViMap, + const double &normalizationFactor_translationalDisplacement, + const double &normalizationFactor_rotationalDisplacement); + + private: + static void computeDesiredReducedModelDisplacements( + const SimulationResults &fullModelResults, + const std::unordered_map &displacementsReducedToFullMap, + Eigen::MatrixX3d &optimalDisplacementsOfReducedModel); + static ReducedModelOptimization::Results runOptimization( + const ReducedModelOptimization::Settings &settings); + std::vector> createScenarios( + const std::shared_ptr &pMesh); + void computeMaps(PatternGeometry &fullModel, PatternGeometry &reducedPattern); + void createSimulationMeshes(PatternGeometry &fullModel, PatternGeometry &reducedModel); + static void initializeOptimizationParameters(const std::shared_ptr &mesh, + const int &optimizationParamters); + + static double objective(long n, const double *x); + DRMSimulationModel simulator; + void computeObjectiveValueNormalizationFactors(); }; void updateMesh(long n, const double *x); #endif // REDUCEDMODELOPTIMIZER_HPP