From d8871cab1746cd7460812d493ae5deb43b7fbff8 Mon Sep 17 00:00:00 2001 From: iasonmanolas Date: Fri, 30 Apr 2021 12:47:50 +0300 Subject: [PATCH] Refactoring --- src/main.cpp | 119 ++++++++++++----------- src/reducedmodeloptimizer.cpp | 173 +++++++++++++++++++++++----------- src/reducedmodeloptimizer.hpp | 13 +-- 3 files changed, 190 insertions(+), 115 deletions(-) diff --git a/src/main.cpp b/src/main.cpp index b83a77a..81e9dde 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -30,13 +30,10 @@ int main(int argc, char *argv[]) { ////Full pattern const std::string filepath_fullPattern = argv[1]; PatternGeometry fullPattern(filepath_fullPattern); - fullPattern.setLabel( - std::filesystem::path(filepath_fullPattern).stem().string()); fullPattern.scale(0.03,interfaceNodeIndex); ////Reduced pattern const std::string filepath_reducedPattern = argv[2]; PatternGeometry reducedPattern(filepath_reducedPattern); - reducedPattern.setLabel(std::filesystem::path(filepath_reducedPattern).stem().string()); reducedPattern.scale(0.03, interfaceNodeIndex); // Set the optization settings @@ -56,70 +53,82 @@ int main(int argc, char *argv[]) { settings_optimization.normalizationStrategy = ReducedPatternOptimization::Settings::NormalizationStrategy::Epsilon; settings_optimization.normalizationParameter = 0.0003; - settings_optimization.solutionAccuracy = 0.001; + settings_optimization.solverAccuracy = 0.001; settings_optimization.objectiveWeights.translational = std::atof(argv[4]); settings_optimization.objectiveWeights.rotational = 2 - std::atof(argv[4]); // Optimize pair - const std::string pairName = - fullPattern.getLabel() + "@" + reducedPattern.getLabel(); + const std::string pairName = fullPattern.getLabel() + "@" + reducedPattern.getLabel(); + const std::string optimizationName = pairName + "(" + + std::to_string(settings_optimization.numberOfFunctionCalls) + + "_" + + to_string_with_precision( + settings_optimization.objectiveWeights.translational) + + ")"; + const bool input_resultDirectoryDefined = argc >= 6; + const std::string optimizationResultsDirectory = input_resultDirectoryDefined + ? argv[5] + : std::filesystem::current_path().append( + "OptimizationResults"); + std::string resultsOutputDir; + bool optimizationResultFolderExists = false; + const std::filesystem::path crashedJobsDirPath(std::filesystem::path(optimizationResultsDirectory) + .append("CrashedJobs") + .append(optimizationName)); + optimizationResultFolderExists |= std::filesystem::exists(crashedJobsDirPath); + const std::filesystem::path convergedJobsDirPath( + std::filesystem::path(optimizationResultsDirectory) + .append("ConvergedJobs") + .append(optimizationName)); + optimizationResultFolderExists |= std::filesystem::exists(convergedJobsDirPath); + ReducedPatternOptimization::Results optimizationResults; + // bool optimizationAlreadyComputed = false; + // if (optimizationResultFolderExists) { + // if (optimizationResults.settings == settings_optimization + // && optimizationResults.load(resultsOutputDir)) { + // optimizationAlreadyComputed = true; + // } + // } + + // if (!optimizationAlreadyComputed) { const std::vector numberOfNodesPerSlot{1, 0, 0, 2, 1, 2, 1}; - assert(interfaceNodeIndex==numberOfNodesPerSlot[0]+numberOfNodesPerSlot[3]); + assert(interfaceNodeIndex == numberOfNodesPerSlot[0] + numberOfNodesPerSlot[3]); ReducedModelOptimizer optimizer(numberOfNodesPerSlot); optimizer.initializePatterns(fullPattern, reducedPattern, settings_optimization.xRanges.size()); - ReducedPatternOptimization::Results optimizationResults - // = optimizer.optimize( - // settings_optimization); - = optimizer.optimize(settings_optimization); + optimizer.optimize(settings_optimization, optimizationResults); + optimizationResults.label = optimizationName; + optimizationResults.baseTriangleFullPattern.copy(fullPattern); + optimizationResults.settings = settings_optimization; // Export results - const bool input_resultDirectoryDefined = argc >= 6; - std::string optimizationResultsDirectory = input_resultDirectoryDefined - ? argv[5] - : std::filesystem::current_path().append( - "OptimizationResults"); - std::string resultsOutputDir; if (optimizationResults.numberOfSimulationCrashes != 0) { - const auto crashedJobsDirPath = std::filesystem::path( - optimizationResultsDirectory.append("CrashedJobs") - .append( - pairName + "(" + std::to_string(settings_optimization.numberOfFunctionCalls) + "_" - + to_string_with_precision(settings_optimization.objectiveWeights.translational) - + ")")); - std::filesystem::create_directories(crashedJobsDirPath); resultsOutputDir = crashedJobsDirPath.string(); - } else { - std::filesystem::path convergedJobsDirPath( - std::filesystem::path(optimizationResultsDirectory) - .append("ConvergedJobs") - .append( - pairName + "(" + std::to_string(settings_optimization.numberOfFunctionCalls) + "_" - + to_string_with_precision(settings_optimization.objectiveWeights.translational) - + ")")); - std::filesystem::create_directories(convergedJobsDirPath); - resultsOutputDir = convergedJobsDirPath.string(); - } - optimizationResults.save(resultsOutputDir); - // Write results in csv - csvFile csv_results({}, false); - // csvFile csv_results(std::filesystem::path(dirPath_thisOptimization) - // .append("results.csv") - // .string(), - // false); - csv_results << "Name"; - optimizationResults.writeHeaderTo(csv_results); - settings_optimization.writeHeaderTo(csv_results); - csv_results << endrow; - csv_results << pairName; - optimizationResults.writeResultsTo(settings_optimization, csv_results); - settings_optimization.writeSettingsTo(csv_results); - csv_results << endrow; + } else { + resultsOutputDir = convergedJobsDirPath.string(); + } + optimizationResults.save(resultsOutputDir); + // } -#ifdef POLYSCOPE_DEFINED - optimizationResults.saveMeshFiles(); - optimizationResults.draw(); -#endif + // Write results in csv + csvFile csv_results({}, false); + // csvFile csv_results(std::filesystem::path(dirPath_thisOptimization) + // .append("results.csv") + // .string(), + // false); + csv_results << "Name"; + optimizationResults.writeHeaderTo(csv_results); + settings_optimization.writeHeaderTo(csv_results); + csv_results << endrow; + csv_results << pairName; + optimizationResults.writeResultsTo(settings_optimization, csv_results); + settings_optimization.writeSettingsTo(csv_results); + csv_results << endrow; - return 0; + //#ifdef POLYSCOPE_DEFINED + // optimizationResults.saveMeshFiles(); + // optimizationResults.draw(); + //#endif + + return 0; } diff --git a/src/reducedmodeloptimizer.cpp b/src/reducedmodeloptimizer.cpp index 31294c7..720978d 100644 --- a/src/reducedmodeloptimizer.cpp +++ b/src/reducedmodeloptimizer.cpp @@ -37,6 +37,7 @@ struct GlobalOptimizationVariables { constructScenarioFunction; FullPatternVertexIndex interfaceViForComputingScenarioError; double desiredMaxDisplacementValue; + double desiredMaxRotationAngle; } global; double ReducedModelOptimizer::computeDisplacementError( @@ -686,11 +687,10 @@ void ReducedModelOptimizer::computeDesiredReducedModelDisplacements( } } -ReducedPatternOptimization::Results ReducedModelOptimizer::getResults( - const dlib::function_evaluation &optimizationResult_dlib, const Settings &settings) +void ReducedModelOptimizer::getResults(const dlib::function_evaluation &optimizationResult_dlib, + const Settings &settings, + ReducedPatternOptimization::Results &results) { - ReducedPatternOptimization::Results results; - results.baseTriangle = global.baseTriangle; //Number of crashes results.numberOfSimulationCrashes = global.numOfSimulationCrashes; @@ -714,10 +714,14 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults( assert(global.minX[xVariableIndex] == optimizationResult_dlib.x(xVariableIndex)); optimalX[xVariableIndex] = optimizationResult_dlib.x(xVariableIndex); +#ifdef POLYSCOPE_DEFINED std::cout << results.optimalXNameValuePairs[xVariableIndex].first << ":" << optimalX[xVariableIndex] << " "; +#endif } +#ifdef POLYSCOPE_DEFINED std::cout << std::endl; +#endif // Compute obj value per simulation scenario and the raw objective value updateMesh(optimalX.size(), optimalX.data()); @@ -786,6 +790,8 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults( // global.rotationalDisplacementNormalizationValues[simulationScenarioIndex]); // assert(test_normalizedTranslationError == normalizedTranslationalError); // assert(test_normalizedRotationalError == normalizedRotationalError); + results.objectiveValue.totalRaw += rawTranslationalError + rawRotationalError; +#ifdef POLYSCOPE_DEFINED std::cout << "Simulation scenario:" << global.reducedPatternSimulationJobs[simulationScenarioIndex]->getLabel() << std::endl; @@ -803,8 +809,8 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults( // = normalizedTranslationalError + normalizedRotationalError; std::cout << "Total Error value:" << results.objectiveValue.perSimulationScenario_total[i] << std::endl; - results.objectiveValue.totalRaw += rawTranslationalError + rawRotationalError; std::cout << std::endl; +#endif } const bool printDebugInfo = false; @@ -828,11 +834,10 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults( // global.reducedPatternSimulationJobs[simulationScenarioIndex]->pMesh->registerForDrawing(); // global.reducedPatternSimulationJobs[simulationScenarioIndex]->pMesh->setLabel(temp); } - - return results; } -ReducedPatternOptimization::Results ReducedModelOptimizer::runOptimization(const Settings &settings) +void ReducedModelOptimizer::runOptimization(const Settings &settings, + ReducedPatternOptimization::Results &results) { global.objectiveValueHistory.clear(); dlib::matrix xMin(global.numberOfOptimizationParameters); @@ -850,13 +855,12 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::runOptimization(const xMax, dlib::max_function_calls(settings.numberOfFunctionCalls), std::chrono::hours(24 * 365 * 290), - settings.solutionAccuracy); + settings.solverAccuracy); auto end = std::chrono::system_clock::now(); auto elapsed = std::chrono::duration_cast(end - start); - ReducedPatternOptimization::Results results = getResults(result_dlib, settings); + getResults(result_dlib, settings, results); results.time = elapsed.count() / 1000.0; - return results; } void ReducedModelOptimizer::constructAxialSimulationScenario( @@ -911,6 +915,8 @@ void ReducedModelOptimizer::constructBendingSimulationScenario( } } +/*NOTE: From the results it seems as if the forced displacements are different in the linear and in the drm + * */ void ReducedModelOptimizer::constructDomeSimulationScenario( const double &forceMagnitude, const std::vector> @@ -948,10 +954,10 @@ void ReducedModelOptimizer::constructDomeSimulationScenario( } else { job.nodalExternalForces[viPair.first] = Vector6d({0, 0, 0, momentAxis[0], momentAxis[1], momentAxis[2]}) * forceMagnitude - / 5; + / 3; job.nodalExternalForces[viPair.second] = Vector6d({0, 0, 0, -momentAxis[0], -momentAxis[1], -momentAxis[2]}) - * forceMagnitude / 5; + * forceMagnitude / 3; job.constrainedVertices[viPair.first] = std::unordered_set{2}; job.constrainedVertices[viPair.second] = std::unordered_set{2}; } @@ -987,9 +993,8 @@ void ReducedModelOptimizer::constructSaddleSimulationScenario( } } -double fullPatternMaxSimulationForceObjective(const double &forceMagnitude) +double fullPatternMaxSimulationForceRotationalObjective(const double &forceMagnitude) { - const double &desiredDisplacementValue = global.desiredMaxDisplacementValue; SimulationJob job; job.pMesh = global.pFullPatternSimulationMesh; global.constructScenarioFunction(forceMagnitude, global.fullPatternInterfaceViPairs, job); @@ -1004,6 +1009,36 @@ double fullPatternMaxSimulationForceObjective(const double &forceMagnitude) // settings.totalTranslationalKineticEnergyThreshold = 1e-7; SimulationResults results = simulator.executeSimulation(std::make_shared(job), settings); + const double &desiredRotationAngle = global.desiredMaxRotationAngle; + const double error = std::abs( + // results.displacements[global.fullPatternInterfaceViPairs[1].first].getTranslation().norm() + results.rotationalDisplacementQuaternion[global.interfaceViForComputingScenarioError] + .angularDistance(Eigen::Quaterniond::Identity()) + - desiredRotationAngle); + + //#ifdef POLYSCOPE_DEFINED + // std::cout << "Force:" << forceMagnitude << " Error is:" << vcg::math::ToDeg(error) << std::endl; + //#endif + return error; +} + +double fullPatternMaxSimulationForceTranslationalObjective(const double &forceMagnitude) +{ + SimulationJob job; + job.pMesh = global.pFullPatternSimulationMesh; + global.constructScenarioFunction(forceMagnitude, global.fullPatternInterfaceViPairs, job); + // ReducedModelOptimizer::axialSimulationScenario(forceMagnitude, + // global.fullPatternInterfaceViPairs, + // job); + + DRMSimulationModel simulator; + DRMSimulationModel::Settings settings; + settings.totalExternalForcesNormPercentageTermination = 1e-2; + // settings.shouldUseTranslationalKineticEnergyThreshold = true; + // settings.totalTranslationalKineticEnergyThreshold = 1e-7; + SimulationResults results = simulator.executeSimulation(std::make_shared(job), + settings); + const double &desiredDisplacementValue = global.desiredMaxDisplacementValue; const double error = std::abs( // results.displacements[global.fullPatternInterfaceViPairs[1].first].getTranslation().norm() results.displacements[global.interfaceViForComputingScenarioError].getTranslation().norm() @@ -1020,30 +1055,45 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat (global.baseTriangle.cP(1) + global.baseTriangle.cP(2)) / 2); - const double optimizationEpsilon = global.desiredMaxDisplacementValue * 0.5; + const double optimizationEpsilon = 1e-2; switch (scenario) { case Axial: + global.desiredMaxDisplacementValue = 0.042 + * vcg::Distance(global.baseTriangle.cP(0), + (global.baseTriangle.cP(1) + + global.baseTriangle.cP(2)) + / 2); global.constructScenarioFunction = &ReducedModelOptimizer::constructAxialSimulationScenario; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first; - dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, + dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective, forceMagnitude, 1e-8, 1e8, optimizationEpsilon); break; case Shear: + global.desiredMaxDisplacementValue = 0.1 + * vcg::Distance(global.baseTriangle.cP(0), + (global.baseTriangle.cP(1) + + global.baseTriangle.cP(2)) + / 2); global.constructScenarioFunction = &ReducedModelOptimizer::constructShearSimulationScenario; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first; - dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, + dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective, forceMagnitude, 1e-8, 1e8, optimizationEpsilon); break; case Bending: + global.desiredMaxDisplacementValue = 0.2 + * vcg::Distance(global.baseTriangle.cP(0), + (global.baseTriangle.cP(1) + + global.baseTriangle.cP(2)) + / 2); global.constructScenarioFunction = &ReducedModelOptimizer::constructBendingSimulationScenario; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first; - dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, + dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective, forceMagnitude, 1e-8, 1e8, @@ -1051,21 +1101,30 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat 200); break; case Dome: - global.desiredMaxDisplacementValue *= 2; + global.desiredMaxRotationAngle = vcg::math::ToRad(40.0); global.constructScenarioFunction = &ReducedModelOptimizer::constructDomeSimulationScenario; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first; - dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, - forceMagnitude, - 1e-8, - 1e8, - optimizationEpsilon, - 200); + dlib::find_min_single_variable( + &fullPatternMaxSimulationForceRotationalObjective, + forceMagnitude, + 1e-8, + 1e8, + vcg::math::ToRad(0.1), + // global.desiredMaxRotationAngle * 0.5, + // optimizationEpsilon, + 500); break; case Saddle: - global.desiredMaxDisplacementValue *= 2; + // global.desiredMaxDisplacementValue *= 2; + global.desiredMaxDisplacementValue = 0.2 + * vcg::Distance(global.baseTriangle.cP(0), + (global.baseTriangle.cP(1) + + global.baseTriangle.cP(2)) + / 2); + // std::cout << "Saddle des disp:" << global.desiredMaxDisplacementValue << std::endl; global.constructScenarioFunction = &ReducedModelOptimizer::constructSaddleSimulationScenario; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first; - dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, + dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective, forceMagnitude, 1e-8, 1e8, @@ -1150,14 +1209,11 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( //// Bending const double maxForceMagnitude_bending = getFullPatternMaxSimulationForce( BaseSimulationScenario::Bending); - const double minForceMagnitude_bending = -maxForceMagnitude_bending; + const double minForceMagnitude_bending = 0; const int numberOfSimulationScenarios_bending = simulationScenariosResolution[BaseSimulationScenario::Bending]; - const double forceMagnitudeStep_bending = numberOfSimulationScenarios_bending == 1 - ? maxForceMagnitude_bending - : (maxForceMagnitude_bending - - minForceMagnitude_bending) - / (numberOfSimulationScenarios_bending - 1); + const double forceMagnitudeStep_bending = (maxForceMagnitude_bending - minForceMagnitude_bending) + / (numberOfSimulationScenarios_bending); const int baseSimulationScenarioIndexOffset_bending = std::accumulate(simulationScenariosResolution.begin(), simulationScenariosResolution.begin() + BaseSimulationScenario::Bending, @@ -1169,7 +1225,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( job.constrainedVertices.clear(); job.label = baseSimulationScenarioNames[BaseSimulationScenario::Bending] + "_" + std::to_string(bendingSimulationScenarioIndex); - const double forceMagnitude = (forceMagnitudeStep_bending * bendingSimulationScenarioIndex + const double forceMagnitude = (forceMagnitudeStep_bending + * (bendingSimulationScenarioIndex + 1) + minForceMagnitude_bending); constructBendingSimulationScenario(forceMagnitude, m_fullPatternOppositeInterfaceViPairs, @@ -1181,13 +1238,11 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( //// Dome const double maxForceMagnitude_dome = getFullPatternMaxSimulationForce( BaseSimulationScenario::Dome); - const double minForceMagnitude_dome = -maxForceMagnitude_dome; + const double minForceMagnitude_dome = 0; const int numberOfSimulationScenarios_dome = simulationScenariosResolution[BaseSimulationScenario::Dome]; - const double forceMagnitudeStep_dome = numberOfSimulationScenarios_dome == 1 - ? maxForceMagnitude_dome - : (maxForceMagnitude_dome - minForceMagnitude_dome) - / (numberOfSimulationScenarios_dome - 1); + const double forceMagnitudeStep_dome = (maxForceMagnitude_dome - minForceMagnitude_dome) + / (numberOfSimulationScenarios_dome); const int baseSimulationScenarioIndexOffset_dome = std::accumulate(simulationScenariosResolution.begin(), simulationScenariosResolution.begin() + BaseSimulationScenario::Dome, @@ -1200,7 +1255,7 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( job.nodalForcedDisplacements.clear(); job.label = baseSimulationScenarioNames[BaseSimulationScenario::Dome] + "_" + std::to_string(domeSimulationScenarioIndex); - const double forceMagnitude = (forceMagnitudeStep_dome * domeSimulationScenarioIndex + const double forceMagnitude = (forceMagnitudeStep_dome * (domeSimulationScenarioIndex + 1) + minForceMagnitude_dome); constructDomeSimulationScenario(forceMagnitude, m_fullPatternOppositeInterfaceViPairs, job); scenarios[baseSimulationScenarioIndexOffset_dome + domeSimulationScenarioIndex] @@ -1210,14 +1265,11 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( //// Saddle const double maxForceMagnitude_saddle = getFullPatternMaxSimulationForce( BaseSimulationScenario::Saddle); - const double minForceMagnitude_saddle = -maxForceMagnitude_saddle; + const double minForceMagnitude_saddle = 0; const int numberOfSimulationScenarios_saddle = simulationScenariosResolution[BaseSimulationScenario::Saddle]; - const double forceMagnitudeStep_saddle = numberOfSimulationScenarios_saddle == 1 - ? maxForceMagnitude_saddle - : (maxForceMagnitude_saddle - - minForceMagnitude_saddle) - / (numberOfSimulationScenarios_saddle - 1); + const double forceMagnitudeStep_saddle = (maxForceMagnitude_saddle - minForceMagnitude_saddle) + / numberOfSimulationScenarios_saddle; const int baseSimulationScenarioIndexOffset_saddle = std::accumulate(simulationScenariosResolution.begin(), simulationScenariosResolution.begin() + BaseSimulationScenario::Saddle, @@ -1230,7 +1282,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( job.nodalForcedDisplacements.clear(); job.label = baseSimulationScenarioNames[BaseSimulationScenario::Saddle] + "_" + std::to_string(saddleSimulationScenarioIndex); - const double forceMagnitude = (forceMagnitudeStep_saddle * saddleSimulationScenarioIndex + const double forceMagnitude = (forceMagnitudeStep_saddle + * (saddleSimulationScenarioIndex + 1) + minForceMagnitude_saddle); constructSaddleSimulationScenario(forceMagnitude, m_fullPatternOppositeInterfaceViPairs, @@ -1239,6 +1292,9 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios( = std::make_shared(job); } +#ifdef POLYSCOPE_DEFINED + std::cout << "Computed full pattern scenario magnitudes" << std::endl; +#endif return scenarios; } @@ -1313,8 +1369,9 @@ void ReducedModelOptimizer::computeObjectiveValueNormalizationFactors() { } } -Results ReducedModelOptimizer::optimize( +void ReducedModelOptimizer::optimize( const Settings &optimizationSettings, + ReducedPatternOptimization::Results &results, const std::vector &desiredBaseSimulationScenarioIndices) { for (int baseSimulationScenarioIndex : desiredBaseSimulationScenarioIndices) { @@ -1353,17 +1410,21 @@ Results ReducedModelOptimizer::optimize( DRMSimulationModel::Settings simulationSettings; simulationSettings.shouldDraw = false; - const bool drawFullPatternSimulationResults = false; +#ifdef POLYSCOPE_DEFINED + const bool drawFullPatternSimulationResults = true; + ; if (drawFullPatternSimulationResults) { global.fullPatternSimulationJobs[0]->pMesh->registerForDrawing( ReducedPatternOptimization::Colors::fullInitial); } +#endif // LinearSimulationModel linearSimulator; for (int simulationScenarioIndex : global.simulationScenarioIndices) { const std::shared_ptr &pFullPatternSimulationJob = global.fullPatternSimulationJobs[simulationScenarioIndex]; SimulationResults fullPatternResults = simulator.executeSimulation(pFullPatternSimulationJob, simulationSettings); +#ifdef POLYSCOPE_DEFINED if (drawFullPatternSimulationResults) { // SimulationResults fullPatternResults_linear = linearSimulator.executeSimulation( // pFullPatternSimulationJob); @@ -1378,6 +1439,7 @@ Results ReducedModelOptimizer::optimize( fullPatternResults.unregister(); // fullPatternResults_linear.unregister(); } +#endif global.fullPatternResults[simulationScenarioIndex] = fullPatternResults; SimulationJob reducedPatternSimulationJob; reducedPatternSimulationJob.pMesh = m_pReducedPatternSimulationMesh; @@ -1388,15 +1450,18 @@ Results ReducedModelOptimizer::optimize( = std::make_shared(reducedPatternSimulationJob); // std::cout << "Ran sim scenario:" << simulationScenarioIndex << std::endl; } + +#ifdef POLYSCOPE_DEFINED if (drawFullPatternSimulationResults) { global.fullPatternSimulationJobs[0]->pMesh->unregister(); } +#endif // if (global.optimizationSettings.normalizationStrategy // != Settings::NormalizationStrategy::NonNormalized) { computeObjectiveValueNormalizationFactors(); // } - PolyscopeInterface::deinitPolyscope(); - Results optResults = runOptimization(optimizationSettings); - - return optResults; +#ifdef POLYSCOPE_DEFINED + std::cout << "Running reduced model optimization" << std::endl; +#endif + runOptimization(optimizationSettings, results); } diff --git a/src/reducedmodeloptimizer.hpp b/src/reducedmodeloptimizer.hpp index 9b1f60d..41c1b0f 100644 --- a/src/reducedmodeloptimizer.hpp +++ b/src/reducedmodeloptimizer.hpp @@ -38,8 +38,9 @@ public: inline static int fanSize{6}; inline static double initialHexagonSize{0.3}; inline static VectorType patternPlaneNormal{0, 0, 1}; - ReducedPatternOptimization::Results optimize( + void optimize( const ReducedPatternOptimization::Settings &xRanges, + ReducedPatternOptimization::Results &results, const std::vector &simulationScenarios = std::vector()); double operator()(const Eigen::VectorXd &x, Eigen::VectorXd &) const; @@ -168,8 +169,8 @@ private: const SimulationResults &fullModelResults, const std::unordered_map &displacementsReducedToFullMap, Eigen::MatrixX3d &optimalDisplacementsOfReducedModel); - static ReducedPatternOptimization::Results runOptimization( - const ReducedPatternOptimization::Settings &settings); + static void runOptimization(const ReducedPatternOptimization::Settings &settings, + ReducedPatternOptimization::Results &results); std::vector> createFullPatternSimulationScenarios( const std::shared_ptr &pMesh); void computeMaps(PatternGeometry &fullModel, PatternGeometry &reducedPattern); @@ -180,9 +181,9 @@ private: static double objective(long n, const double *x); DRMSimulationModel simulator; void computeObjectiveValueNormalizationFactors(); - static ReducedPatternOptimization::Results getResults( - const dlib::function_evaluation &optimizationResult_dlib, - const ReducedPatternOptimization::Settings &settings); + static void getResults(const dlib::function_evaluation &optimizationResult_dlib, + const ReducedPatternOptimization::Settings &settings, + ReducedPatternOptimization::Results &results); std::vector getFullPatternMaxSimulationForces(); double getFullPatternMaxSimulationForce( const ReducedPatternOptimization::BaseSimulationScenario &scenario);