manolas
/
ReducedModelOptimization
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

Refactoring

This commit is contained in:
iasonmanolas 2021-04-30 12:47:50 +03:00
parent 62ce79d959
commit d8871cab17
3 changed files with 190 additions and 115 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

119
src/main.cpp
View File

@ -30,13 +30,10 @@ int main(int argc, char *argv[]) {
////Full pattern ////Full pattern
const std::string filepath_fullPattern = argv[1]; const std::string filepath_fullPattern = argv[1];
PatternGeometry fullPattern(filepath_fullPattern); PatternGeometry fullPattern(filepath_fullPattern);
fullPattern.setLabel(
std::filesystem::path(filepath_fullPattern).stem().string());
fullPattern.scale(0.03,interfaceNodeIndex); fullPattern.scale(0.03,interfaceNodeIndex);
////Reduced pattern ////Reduced pattern
const std::string filepath_reducedPattern = argv[2]; const std::string filepath_reducedPattern = argv[2];
PatternGeometry reducedPattern(filepath_reducedPattern); PatternGeometry reducedPattern(filepath_reducedPattern);
reducedPattern.setLabel(std::filesystem::path(filepath_reducedPattern).stem().string());
reducedPattern.scale(0.03, interfaceNodeIndex); reducedPattern.scale(0.03, interfaceNodeIndex);
// Set the optization settings // Set the optization settings
@ -56,70 +53,82 @@ int main(int argc, char *argv[]) {
settings_optimization.normalizationStrategy settings_optimization.normalizationStrategy
= ReducedPatternOptimization::Settings::NormalizationStrategy::Epsilon; = ReducedPatternOptimization::Settings::NormalizationStrategy::Epsilon;
settings_optimization.normalizationParameter = 0.0003; 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.translational = std::atof(argv[4]);
settings_optimization.objectiveWeights.rotational = 2 - std::atof(argv[4]); settings_optimization.objectiveWeights.rotational = 2 - std::atof(argv[4]);
// Optimize pair // Optimize pair
const std::string pairName = const std::string pairName = fullPattern.getLabel() + "@" + reducedPattern.getLabel();
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<size_t> numberOfNodesPerSlot{1, 0, 0, 2, 1, 2, 1}; const std::vector<size_t> numberOfNodesPerSlot{1, 0, 0, 2, 1, 2, 1};
assert(interfaceNodeIndex==numberOfNodesPerSlot[0]+numberOfNodesPerSlot[3]); assert(interfaceNodeIndex == numberOfNodesPerSlot[0] + numberOfNodesPerSlot[3]);
ReducedModelOptimizer optimizer(numberOfNodesPerSlot); ReducedModelOptimizer optimizer(numberOfNodesPerSlot);
optimizer.initializePatterns(fullPattern, reducedPattern, settings_optimization.xRanges.size()); optimizer.initializePatterns(fullPattern, reducedPattern, settings_optimization.xRanges.size());
ReducedPatternOptimization::Results optimizationResults optimizer.optimize(settings_optimization, optimizationResults);
// = optimizer.optimize( optimizationResults.label = optimizationName;
// settings_optimization); optimizationResults.baseTriangleFullPattern.copy(fullPattern);
= optimizer.optimize(settings_optimization); optimizationResults.settings = settings_optimization;
// Export results // 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) { 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(); resultsOutputDir = crashedJobsDirPath.string();
} else { } else {
std::filesystem::path convergedJobsDirPath( resultsOutputDir = convergedJobsDirPath.string();
std::filesystem::path(optimizationResultsDirectory) }
.append("ConvergedJobs") optimizationResults.save(resultsOutputDir);
.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;
#ifdef POLYSCOPE_DEFINED // Write results in csv
optimizationResults.saveMeshFiles(); csvFile csv_results({}, false);
optimizationResults.draw(); // csvFile csv_results(std::filesystem::path(dirPath_thisOptimization)
#endif // .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;
} }

173
src/reducedmodeloptimizer.cpp
View File

@ -37,6 +37,7 @@ struct GlobalOptimizationVariables {
constructScenarioFunction; constructScenarioFunction;
FullPatternVertexIndex interfaceViForComputingScenarioError; FullPatternVertexIndex interfaceViForComputingScenarioError;
double desiredMaxDisplacementValue; double desiredMaxDisplacementValue;
double desiredMaxRotationAngle;
} global; } global;
double ReducedModelOptimizer::computeDisplacementError( double ReducedModelOptimizer::computeDisplacementError(
@ -686,11 +687,10 @@ void ReducedModelOptimizer::computeDesiredReducedModelDisplacements(
} }
} }
ReducedPatternOptimization::Results ReducedModelOptimizer::getResults( void ReducedModelOptimizer::getResults(const dlib::function_evaluation &optimizationResult_dlib,
const dlib::function_evaluation &optimizationResult_dlib, const Settings &settings) const Settings &settings,
ReducedPatternOptimization::Results &results)
{ {
ReducedPatternOptimization::Results results;
results.baseTriangle = global.baseTriangle; results.baseTriangle = global.baseTriangle;
//Number of crashes //Number of crashes
results.numberOfSimulationCrashes = global.numOfSimulationCrashes; results.numberOfSimulationCrashes = global.numOfSimulationCrashes;
@ -714,10 +714,14 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults(
assert(global.minX[xVariableIndex] == optimizationResult_dlib.x(xVariableIndex)); assert(global.minX[xVariableIndex] == optimizationResult_dlib.x(xVariableIndex));
optimalX[xVariableIndex] = optimizationResult_dlib.x(xVariableIndex); optimalX[xVariableIndex] = optimizationResult_dlib.x(xVariableIndex);
#ifdef POLYSCOPE_DEFINED
std::cout << results.optimalXNameValuePairs[xVariableIndex].first << ":" std::cout << results.optimalXNameValuePairs[xVariableIndex].first << ":"
<< optimalX[xVariableIndex] << " "; << optimalX[xVariableIndex] << " ";
#endif
} }
#ifdef POLYSCOPE_DEFINED
std::cout << std::endl; std::cout << std::endl;
#endif
// Compute obj value per simulation scenario and the raw objective value // Compute obj value per simulation scenario and the raw objective value
updateMesh(optimalX.size(), optimalX.data()); updateMesh(optimalX.size(), optimalX.data());
@ -786,6 +790,8 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults(
// global.rotationalDisplacementNormalizationValues[simulationScenarioIndex]); // global.rotationalDisplacementNormalizationValues[simulationScenarioIndex]);
// assert(test_normalizedTranslationError == normalizedTranslationalError); // assert(test_normalizedTranslationError == normalizedTranslationalError);
// assert(test_normalizedRotationalError == normalizedRotationalError); // assert(test_normalizedRotationalError == normalizedRotationalError);
results.objectiveValue.totalRaw += rawTranslationalError + rawRotationalError;
#ifdef POLYSCOPE_DEFINED
std::cout << "Simulation scenario:" std::cout << "Simulation scenario:"
<< global.reducedPatternSimulationJobs[simulationScenarioIndex]->getLabel() << global.reducedPatternSimulationJobs[simulationScenarioIndex]->getLabel()
<< std::endl; << std::endl;
@ -803,8 +809,8 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults(
// = normalizedTranslationalError + normalizedRotationalError; // = normalizedTranslationalError + normalizedRotationalError;
std::cout << "Total Error value:" << results.objectiveValue.perSimulationScenario_total[i] std::cout << "Total Error value:" << results.objectiveValue.perSimulationScenario_total[i]
<< std::endl; << std::endl;
results.objectiveValue.totalRaw += rawTranslationalError + rawRotationalError;
std::cout << std::endl; std::cout << std::endl;
#endif
} }
const bool printDebugInfo = false; const bool printDebugInfo = false;
@ -828,11 +834,10 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::getResults(
// global.reducedPatternSimulationJobs[simulationScenarioIndex]->pMesh->registerForDrawing(); // global.reducedPatternSimulationJobs[simulationScenarioIndex]->pMesh->registerForDrawing();
// global.reducedPatternSimulationJobs[simulationScenarioIndex]->pMesh->setLabel(temp); // 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(); global.objectiveValueHistory.clear();
dlib::matrix<double, 0, 1> xMin(global.numberOfOptimizationParameters); dlib::matrix<double, 0, 1> xMin(global.numberOfOptimizationParameters);
@ -850,13 +855,12 @@ ReducedPatternOptimization::Results ReducedModelOptimizer::runOptimization(const
xMax, xMax,
dlib::max_function_calls(settings.numberOfFunctionCalls), dlib::max_function_calls(settings.numberOfFunctionCalls),
std::chrono::hours(24 * 365 * 290), std::chrono::hours(24 * 365 * 290),
settings.solutionAccuracy); settings.solverAccuracy);
auto end = std::chrono::system_clock::now(); auto end = std::chrono::system_clock::now();
auto elapsed = auto elapsed =
std::chrono::duration_cast<std::chrono::milliseconds>(end - start); std::chrono::duration_cast<std::chrono::milliseconds>(end - start);
ReducedPatternOptimization::Results results = getResults(result_dlib, settings); getResults(result_dlib, settings, results);
results.time = elapsed.count() / 1000.0; results.time = elapsed.count() / 1000.0;
return results;
} }
void ReducedModelOptimizer::constructAxialSimulationScenario( 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( void ReducedModelOptimizer::constructDomeSimulationScenario(
const double &forceMagnitude, const double &forceMagnitude,
const std::vector<std::pair<FullPatternVertexIndex, FullPatternVertexIndex>> const std::vector<std::pair<FullPatternVertexIndex, FullPatternVertexIndex>>
@ -948,10 +954,10 @@ void ReducedModelOptimizer::constructDomeSimulationScenario(
} else { } else {
job.nodalExternalForces[viPair.first] job.nodalExternalForces[viPair.first]
= Vector6d({0, 0, 0, momentAxis[0], momentAxis[1], momentAxis[2]}) * forceMagnitude = Vector6d({0, 0, 0, momentAxis[0], momentAxis[1], momentAxis[2]}) * forceMagnitude
/ 5; / 3;
job.nodalExternalForces[viPair.second] job.nodalExternalForces[viPair.second]
= Vector6d({0, 0, 0, -momentAxis[0], -momentAxis[1], -momentAxis[2]}) = Vector6d({0, 0, 0, -momentAxis[0], -momentAxis[1], -momentAxis[2]})
* forceMagnitude / 5; * forceMagnitude / 3;
job.constrainedVertices[viPair.first] = std::unordered_set<DoFType>{2}; job.constrainedVertices[viPair.first] = std::unordered_set<DoFType>{2};
job.constrainedVertices[viPair.second] = std::unordered_set<DoFType>{2}; job.constrainedVertices[viPair.second] = std::unordered_set<DoFType>{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; SimulationJob job;
job.pMesh = global.pFullPatternSimulationMesh; job.pMesh = global.pFullPatternSimulationMesh;
global.constructScenarioFunction(forceMagnitude, global.fullPatternInterfaceViPairs, job); global.constructScenarioFunction(forceMagnitude, global.fullPatternInterfaceViPairs, job);
@ -1004,6 +1009,36 @@ double fullPatternMaxSimulationForceObjective(const double &forceMagnitude)
// settings.totalTranslationalKineticEnergyThreshold = 1e-7; // settings.totalTranslationalKineticEnergyThreshold = 1e-7;
SimulationResults results = simulator.executeSimulation(std::make_shared<SimulationJob>(job), SimulationResults results = simulator.executeSimulation(std::make_shared<SimulationJob>(job),
settings); 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<SimulationJob>(job),
settings);
const double &desiredDisplacementValue = global.desiredMaxDisplacementValue;
const double error = std::abs( const double error = std::abs(
// results.displacements[global.fullPatternInterfaceViPairs[1].first].getTranslation().norm() // results.displacements[global.fullPatternInterfaceViPairs[1].first].getTranslation().norm()
results.displacements[global.interfaceViForComputingScenarioError].getTranslation().norm() results.displacements[global.interfaceViForComputingScenarioError].getTranslation().norm()
@ -1020,30 +1055,45 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat
(global.baseTriangle.cP(1) (global.baseTriangle.cP(1)
+ global.baseTriangle.cP(2)) + global.baseTriangle.cP(2))
/ 2); / 2);
const double optimizationEpsilon = global.desiredMaxDisplacementValue * 0.5; const double optimizationEpsilon = 1e-2;
switch (scenario) { switch (scenario) {
case Axial: 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.constructScenarioFunction = &ReducedModelOptimizer::constructAxialSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude, forceMagnitude,
1e-8, 1e-8,
1e8, 1e8,
optimizationEpsilon); optimizationEpsilon);
break; break;
case Shear: 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.constructScenarioFunction = &ReducedModelOptimizer::constructShearSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude, forceMagnitude,
1e-8, 1e-8,
1e8, 1e8,
optimizationEpsilon); optimizationEpsilon);
break; break;
case Bending: 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.constructScenarioFunction = &ReducedModelOptimizer::constructBendingSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude, forceMagnitude,
1e-8, 1e-8,
1e8, 1e8,
@ -1051,21 +1101,30 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat
200); 200);
break; break;
case Dome: case Dome:
global.desiredMaxDisplacementValue *= 2; global.desiredMaxRotationAngle = vcg::math::ToRad(40.0);
global.constructScenarioFunction = &ReducedModelOptimizer::constructDomeSimulationScenario; global.constructScenarioFunction = &ReducedModelOptimizer::constructDomeSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, dlib::find_min_single_variable(
forceMagnitude, &fullPatternMaxSimulationForceRotationalObjective,
1e-8, forceMagnitude,
1e8, 1e-8,
optimizationEpsilon, 1e8,
200); vcg::math::ToRad(0.1),
// global.desiredMaxRotationAngle * 0.5,
// optimizationEpsilon,
500);
break; break;
case Saddle: 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.constructScenarioFunction = &ReducedModelOptimizer::constructSaddleSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first; global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceObjective, dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude, forceMagnitude,
1e-8, 1e-8,
1e8, 1e8,
@ -1150,14 +1209,11 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
//// Bending //// Bending
const double maxForceMagnitude_bending = getFullPatternMaxSimulationForce( const double maxForceMagnitude_bending = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Bending); BaseSimulationScenario::Bending);
const double minForceMagnitude_bending = -maxForceMagnitude_bending; const double minForceMagnitude_bending = 0;
const int numberOfSimulationScenarios_bending const int numberOfSimulationScenarios_bending
= simulationScenariosResolution[BaseSimulationScenario::Bending]; = simulationScenariosResolution[BaseSimulationScenario::Bending];
const double forceMagnitudeStep_bending = numberOfSimulationScenarios_bending == 1 const double forceMagnitudeStep_bending = (maxForceMagnitude_bending - minForceMagnitude_bending)
? maxForceMagnitude_bending / (numberOfSimulationScenarios_bending);
: (maxForceMagnitude_bending
- minForceMagnitude_bending)
/ (numberOfSimulationScenarios_bending - 1);
const int baseSimulationScenarioIndexOffset_bending const int baseSimulationScenarioIndexOffset_bending
= std::accumulate(simulationScenariosResolution.begin(), = std::accumulate(simulationScenariosResolution.begin(),
simulationScenariosResolution.begin() + BaseSimulationScenario::Bending, simulationScenariosResolution.begin() + BaseSimulationScenario::Bending,
@ -1169,7 +1225,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
job.constrainedVertices.clear(); job.constrainedVertices.clear();
job.label = baseSimulationScenarioNames[BaseSimulationScenario::Bending] + "_" job.label = baseSimulationScenarioNames[BaseSimulationScenario::Bending] + "_"
+ std::to_string(bendingSimulationScenarioIndex); + std::to_string(bendingSimulationScenarioIndex);
const double forceMagnitude = (forceMagnitudeStep_bending * bendingSimulationScenarioIndex const double forceMagnitude = (forceMagnitudeStep_bending
* (bendingSimulationScenarioIndex + 1)
+ minForceMagnitude_bending); + minForceMagnitude_bending);
constructBendingSimulationScenario(forceMagnitude, constructBendingSimulationScenario(forceMagnitude,
m_fullPatternOppositeInterfaceViPairs, m_fullPatternOppositeInterfaceViPairs,
@ -1181,13 +1238,11 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
//// Dome //// Dome
const double maxForceMagnitude_dome = getFullPatternMaxSimulationForce( const double maxForceMagnitude_dome = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Dome); BaseSimulationScenario::Dome);
const double minForceMagnitude_dome = -maxForceMagnitude_dome; const double minForceMagnitude_dome = 0;
const int numberOfSimulationScenarios_dome const int numberOfSimulationScenarios_dome
= simulationScenariosResolution[BaseSimulationScenario::Dome]; = simulationScenariosResolution[BaseSimulationScenario::Dome];
const double forceMagnitudeStep_dome = numberOfSimulationScenarios_dome == 1 const double forceMagnitudeStep_dome = (maxForceMagnitude_dome - minForceMagnitude_dome)
? maxForceMagnitude_dome / (numberOfSimulationScenarios_dome);
: (maxForceMagnitude_dome - minForceMagnitude_dome)
/ (numberOfSimulationScenarios_dome - 1);
const int baseSimulationScenarioIndexOffset_dome const int baseSimulationScenarioIndexOffset_dome
= std::accumulate(simulationScenariosResolution.begin(), = std::accumulate(simulationScenariosResolution.begin(),
simulationScenariosResolution.begin() + BaseSimulationScenario::Dome, simulationScenariosResolution.begin() + BaseSimulationScenario::Dome,
@ -1200,7 +1255,7 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
job.nodalForcedDisplacements.clear(); job.nodalForcedDisplacements.clear();
job.label = baseSimulationScenarioNames[BaseSimulationScenario::Dome] + "_" job.label = baseSimulationScenarioNames[BaseSimulationScenario::Dome] + "_"
+ std::to_string(domeSimulationScenarioIndex); + std::to_string(domeSimulationScenarioIndex);
const double forceMagnitude = (forceMagnitudeStep_dome * domeSimulationScenarioIndex const double forceMagnitude = (forceMagnitudeStep_dome * (domeSimulationScenarioIndex + 1)
+ minForceMagnitude_dome); + minForceMagnitude_dome);
constructDomeSimulationScenario(forceMagnitude, m_fullPatternOppositeInterfaceViPairs, job); constructDomeSimulationScenario(forceMagnitude, m_fullPatternOppositeInterfaceViPairs, job);
scenarios[baseSimulationScenarioIndexOffset_dome + domeSimulationScenarioIndex] scenarios[baseSimulationScenarioIndexOffset_dome + domeSimulationScenarioIndex]
@ -1210,14 +1265,11 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
//// Saddle //// Saddle
const double maxForceMagnitude_saddle = getFullPatternMaxSimulationForce( const double maxForceMagnitude_saddle = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Saddle); BaseSimulationScenario::Saddle);
const double minForceMagnitude_saddle = -maxForceMagnitude_saddle; const double minForceMagnitude_saddle = 0;
const int numberOfSimulationScenarios_saddle const int numberOfSimulationScenarios_saddle
= simulationScenariosResolution[BaseSimulationScenario::Saddle]; = simulationScenariosResolution[BaseSimulationScenario::Saddle];
const double forceMagnitudeStep_saddle = numberOfSimulationScenarios_saddle == 1 const double forceMagnitudeStep_saddle = (maxForceMagnitude_saddle - minForceMagnitude_saddle)
? maxForceMagnitude_saddle / numberOfSimulationScenarios_saddle;
: (maxForceMagnitude_saddle
- minForceMagnitude_saddle)
/ (numberOfSimulationScenarios_saddle - 1);
const int baseSimulationScenarioIndexOffset_saddle const int baseSimulationScenarioIndexOffset_saddle
= std::accumulate(simulationScenariosResolution.begin(), = std::accumulate(simulationScenariosResolution.begin(),
simulationScenariosResolution.begin() + BaseSimulationScenario::Saddle, simulationScenariosResolution.begin() + BaseSimulationScenario::Saddle,
@ -1230,7 +1282,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
job.nodalForcedDisplacements.clear(); job.nodalForcedDisplacements.clear();
job.label = baseSimulationScenarioNames[BaseSimulationScenario::Saddle] + "_" job.label = baseSimulationScenarioNames[BaseSimulationScenario::Saddle] + "_"
+ std::to_string(saddleSimulationScenarioIndex); + std::to_string(saddleSimulationScenarioIndex);
const double forceMagnitude = (forceMagnitudeStep_saddle * saddleSimulationScenarioIndex const double forceMagnitude = (forceMagnitudeStep_saddle
* (saddleSimulationScenarioIndex + 1)
+ minForceMagnitude_saddle); + minForceMagnitude_saddle);
constructSaddleSimulationScenario(forceMagnitude, constructSaddleSimulationScenario(forceMagnitude,
m_fullPatternOppositeInterfaceViPairs, m_fullPatternOppositeInterfaceViPairs,
@ -1239,6 +1292,9 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
= std::make_shared<SimulationJob>(job); = std::make_shared<SimulationJob>(job);
} }
#ifdef POLYSCOPE_DEFINED
std::cout << "Computed full pattern scenario magnitudes" << std::endl;
#endif
return scenarios; return scenarios;
} }
@ -1313,8 +1369,9 @@ void ReducedModelOptimizer::computeObjectiveValueNormalizationFactors() {
} }
} }
Results ReducedModelOptimizer::optimize( void ReducedModelOptimizer::optimize(
const Settings &optimizationSettings, const Settings &optimizationSettings,
ReducedPatternOptimization::Results &results,
const std::vector<BaseSimulationScenario> &desiredBaseSimulationScenarioIndices) const std::vector<BaseSimulationScenario> &desiredBaseSimulationScenarioIndices)
{ {
for (int baseSimulationScenarioIndex : desiredBaseSimulationScenarioIndices) { for (int baseSimulationScenarioIndex : desiredBaseSimulationScenarioIndices) {
@ -1353,17 +1410,21 @@ Results ReducedModelOptimizer::optimize(
DRMSimulationModel::Settings simulationSettings; DRMSimulationModel::Settings simulationSettings;
simulationSettings.shouldDraw = false; simulationSettings.shouldDraw = false;
const bool drawFullPatternSimulationResults = false; #ifdef POLYSCOPE_DEFINED
const bool drawFullPatternSimulationResults = true;
;
if (drawFullPatternSimulationResults) { if (drawFullPatternSimulationResults) {
global.fullPatternSimulationJobs[0]->pMesh->registerForDrawing( global.fullPatternSimulationJobs[0]->pMesh->registerForDrawing(
ReducedPatternOptimization::Colors::fullInitial); ReducedPatternOptimization::Colors::fullInitial);
} }
#endif
// LinearSimulationModel linearSimulator; // LinearSimulationModel linearSimulator;
for (int simulationScenarioIndex : global.simulationScenarioIndices) { for (int simulationScenarioIndex : global.simulationScenarioIndices) {
const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob
= global.fullPatternSimulationJobs[simulationScenarioIndex]; = global.fullPatternSimulationJobs[simulationScenarioIndex];
SimulationResults fullPatternResults = simulator.executeSimulation(pFullPatternSimulationJob, SimulationResults fullPatternResults = simulator.executeSimulation(pFullPatternSimulationJob,
simulationSettings); simulationSettings);
#ifdef POLYSCOPE_DEFINED
if (drawFullPatternSimulationResults) { if (drawFullPatternSimulationResults) {
// SimulationResults fullPatternResults_linear = linearSimulator.executeSimulation( // SimulationResults fullPatternResults_linear = linearSimulator.executeSimulation(
// pFullPatternSimulationJob); // pFullPatternSimulationJob);
@ -1378,6 +1439,7 @@ Results ReducedModelOptimizer::optimize(
fullPatternResults.unregister(); fullPatternResults.unregister();
// fullPatternResults_linear.unregister(); // fullPatternResults_linear.unregister();
} }
#endif
global.fullPatternResults[simulationScenarioIndex] = fullPatternResults; global.fullPatternResults[simulationScenarioIndex] = fullPatternResults;
SimulationJob reducedPatternSimulationJob; SimulationJob reducedPatternSimulationJob;
reducedPatternSimulationJob.pMesh = m_pReducedPatternSimulationMesh; reducedPatternSimulationJob.pMesh = m_pReducedPatternSimulationMesh;
@ -1388,15 +1450,18 @@ Results ReducedModelOptimizer::optimize(
= std::make_shared<SimulationJob>(reducedPatternSimulationJob); = std::make_shared<SimulationJob>(reducedPatternSimulationJob);
// std::cout << "Ran sim scenario:" << simulationScenarioIndex << std::endl; // std::cout << "Ran sim scenario:" << simulationScenarioIndex << std::endl;
} }
#ifdef POLYSCOPE_DEFINED
if (drawFullPatternSimulationResults) { if (drawFullPatternSimulationResults) {
global.fullPatternSimulationJobs[0]->pMesh->unregister(); global.fullPatternSimulationJobs[0]->pMesh->unregister();
} }
#endif
// if (global.optimizationSettings.normalizationStrategy // if (global.optimizationSettings.normalizationStrategy
// != Settings::NormalizationStrategy::NonNormalized) { // != Settings::NormalizationStrategy::NonNormalized) {
computeObjectiveValueNormalizationFactors(); computeObjectiveValueNormalizationFactors();
// } // }
PolyscopeInterface::deinitPolyscope(); #ifdef POLYSCOPE_DEFINED
Results optResults = runOptimization(optimizationSettings); std::cout << "Running reduced model optimization" << std::endl;
#endif
return optResults; runOptimization(optimizationSettings, results);
} }

13
src/reducedmodeloptimizer.hpp
View File

@ -38,8 +38,9 @@ public:
inline static int fanSize{6}; inline static int fanSize{6};
inline static double initialHexagonSize{0.3}; inline static double initialHexagonSize{0.3};
inline static VectorType patternPlaneNormal{0, 0, 1}; inline static VectorType patternPlaneNormal{0, 0, 1};
ReducedPatternOptimization::Results optimize( void optimize(
const ReducedPatternOptimization::Settings &xRanges, const ReducedPatternOptimization::Settings &xRanges,
ReducedPatternOptimization::Results &results,
const std::vector<ReducedPatternOptimization::BaseSimulationScenario> &simulationScenarios const std::vector<ReducedPatternOptimization::BaseSimulationScenario> &simulationScenarios
= std::vector<ReducedPatternOptimization::BaseSimulationScenario>()); = std::vector<ReducedPatternOptimization::BaseSimulationScenario>());
double operator()(const Eigen::VectorXd &x, Eigen::VectorXd &) const; double operator()(const Eigen::VectorXd &x, Eigen::VectorXd &) const;
@ -168,8 +169,8 @@ private:
const SimulationResults &fullModelResults, const SimulationResults &fullModelResults,
const std::unordered_map<size_t, size_t> &displacementsReducedToFullMap, const std::unordered_map<size_t, size_t> &displacementsReducedToFullMap,
Eigen::MatrixX3d &optimalDisplacementsOfReducedModel); Eigen::MatrixX3d &optimalDisplacementsOfReducedModel);
static ReducedPatternOptimization::Results runOptimization( static void runOptimization(const ReducedPatternOptimization::Settings &settings,
const ReducedPatternOptimization::Settings &settings); ReducedPatternOptimization::Results &results);
std::vector<std::shared_ptr<SimulationJob>> createFullPatternSimulationScenarios( std::vector<std::shared_ptr<SimulationJob>> createFullPatternSimulationScenarios(
const std::shared_ptr<SimulationMesh> &pMesh); const std::shared_ptr<SimulationMesh> &pMesh);
void computeMaps(PatternGeometry &fullModel, PatternGeometry &reducedPattern); void computeMaps(PatternGeometry &fullModel, PatternGeometry &reducedPattern);
@ -180,9 +181,9 @@ private:
static double objective(long n, const double *x); static double objective(long n, const double *x);
DRMSimulationModel simulator; DRMSimulationModel simulator;
void computeObjectiveValueNormalizationFactors(); void computeObjectiveValueNormalizationFactors();
static ReducedPatternOptimization::Results getResults( static void getResults(const dlib::function_evaluation &optimizationResult_dlib,
const dlib::function_evaluation &optimizationResult_dlib, const ReducedPatternOptimization::Settings &settings,
const ReducedPatternOptimization::Settings &settings); ReducedPatternOptimization::Results &results);
std::vector<double> getFullPatternMaxSimulationForces(); std::vector<double> getFullPatternMaxSimulationForces();
double getFullPatternMaxSimulationForce( double getFullPatternMaxSimulationForce(
const ReducedPatternOptimization::BaseSimulationScenario &scenario); const ReducedPatternOptimization::BaseSimulationScenario &scenario);