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

Added max iterations on all drm simulations.using global optimizer for dome magnitude.Added bool variable for success in force magnitude search.

This commit is contained in:
iasonmanolas 2021-06-30 10:24:54 +03:00
parent a9c5ca3fc6
commit 0863db658d
2 changed files with 127 additions and 84 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

208
src/reducedmodeloptimizer.cpp
View File

@ -692,9 +692,8 @@ void ReducedModelOptimizer::getResults(const dlib::function_evaluation &optimiza
const Settings &settings,
ReducedPatternOptimization::Results &results)
{
results.baseTriangle = global.baseTriangle;
//Number of crashes
results.numberOfSimulationCrashes = global.numOfSimulationCrashes;
// results.numberOfSimulationCrashes = global.numOfSimulationCrashes;
//Value of optimal objective Y
results.objectiveValue.total = optimizationResult_dlib.y;
//Optimal X values
@ -932,7 +931,7 @@ void ReducedModelOptimizer::constructDomeSimulationScenario(
= Eigen::Vector3d(-interfaceVector[0],
-interfaceVector[1],
0)
* 0.001
* 0.01
* std::abs(
forceMagnitude); //NOTE:Should the forced displacement change relatively to the magnitude?
// * std::abs(forceMagnitude / maxForceMagnitude_dome);
@ -1000,21 +999,31 @@ double fullPatternMaxSimulationForceRotationalObjective(const double &forceMagni
DRMSimulationModel simulator;
DRMSimulationModel::Settings settings;
settings.totalExternalForcesNormPercentageTermination = 1e-2;
// settings.shouldUseTranslationalKineticEnergyThreshold = true;
// settings.totalTranslationalKineticEnergyThreshold = 1e-7;
// settings.totalResidualForcesNormThreshold = 1e-3;
// settings.totalTranslationalKineticEnergyThreshold = 1e-6;
// settings.shouldUseTranslationalKineticEnergyThreshold = true;
settings.shouldDraw = true;
// settings.isDebugMode = true;
// settings.drawingStep = 500;
// settings.beVerbose = true;
// settings.debugModeStep = 200000;
// settings.shouldCreatePlots = true;
settings.maxDRMIterations = 100000;
SimulationResults results = simulator.executeSimulation(std::make_shared<SimulationJob>(job),
settings);
const double &desiredRotationAngle = global.desiredMaxRotationAngle;
if (!results.converged) {
return std::numeric_limits<double>::max();
}
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
#ifdef POLYSCOPE_DEFINED
std::cout << "Force:" << forceMagnitude << " Error is:" << vcg::math::ToDeg(error) << std::endl;
#endif
return error;
}
@ -1029,36 +1038,47 @@ double fullPatternMaxSimulationForceTranslationalObjective(const double &forceMa
DRMSimulationModel simulator;
DRMSimulationModel::Settings settings;
settings.totalExternalForcesNormPercentageTermination = 1e-2;
settings.shouldUseTranslationalKineticEnergyThreshold = true;
settings.totalTranslationalKineticEnergyThreshold = 1e-14;
settings.shouldUseTranslationalKineticEnergyThreshold = true;
settings.totalResidualForcesNormThreshold = 1e-3;
// settings.totalTranslationalKineticEnergyThreshold = 1e-10;
// settings.shouldUseTranslationalKineticEnergyThreshold = true;
// settings.shouldDraw = true;
// settings.isDebugMode = true;
// settings.drawingStep = 200000;
// settings.beVerbose = true;
// settings.debugModeStep = 200000;
settings.maxDRMIterations = 50000;
#ifdef POLYSCOPE_DEFINED
// settings.shouldDraw = true;
// settings.isDebugMode = true;
// settings.drawingStep = 100000;
// settings.debugModeStep = 10000;
// settings.shouldCreatePlots = true;
#endif
SimulationResults results = simulator.executeSimulation(std::make_shared<SimulationJob>(job),
settings);
SimulationResults results = simulator.executeSimulation(std::make_shared<SimulationJob>(job));
const double &desiredDisplacementValue = global.desiredMaxDisplacementValue;
if (!results.converged) {
return std::numeric_limits<double>::max();
}
const double error = std::abs(
// results.displacements[global.fullPatternInterfaceViPairs[1].first].getTranslation().norm()
results.displacements[global.interfaceViForComputingScenarioError].getTranslation().norm()
- desiredDisplacementValue);
#ifdef POLYSCOPE_DEFINED
std::cout << "Force:" << forceMagnitude << " Error is:" << error << std::endl;
#endif
return error;
}
double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulationScenario &scenario)
double ReducedModelOptimizer::getFullPatternMaxSimulationForce(
const BaseSimulationScenario &scenario, bool &wasSuccessful)
{
double forceMagnitude = 1;
global.desiredMaxDisplacementValue = 0.1
* vcg::Distance(global.baseTriangle.cP(0),
(global.baseTriangle.cP(1)
+ global.baseTriangle.cP(2))
/ 2);
const double optimizationEpsilon = 1e-1;
const double forceMagnitudeEpsilon = 1e-2;
double minimumError;
double translationalOptimizationEpsilon;
dlib::function_evaluation result;
switch (scenario) {
case Axial:
global.desiredMaxDisplacementValue = 0.03
@ -1068,11 +1088,14 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat
/ 2);
global.constructScenarioFunction = &ReducedModelOptimizer::constructAxialSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-8,
1e8,
optimizationEpsilon);
minimumError
= dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-2,
1e2,
forceMagnitudeEpsilon);
translationalOptimizationEpsilon = 1e-2 * global.desiredMaxDisplacementValue;
wasSuccessful = minimumError < translationalOptimizationEpsilon;
break;
case Shear:
global.desiredMaxDisplacementValue = 0.04
@ -1082,11 +1105,14 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat
/ 2);
global.constructScenarioFunction = &ReducedModelOptimizer::constructShearSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-8,
1e8,
optimizationEpsilon);
minimumError
= dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-2,
1e2,
forceMagnitudeEpsilon);
translationalOptimizationEpsilon = 1e-2 * global.desiredMaxDisplacementValue;
wasSuccessful = minimumError < translationalOptimizationEpsilon;
break;
case Bending:
global.desiredMaxDisplacementValue = 0.1
@ -1096,28 +1122,33 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat
/ 2);
global.constructScenarioFunction = &ReducedModelOptimizer::constructBendingSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-8,
1e8,
optimizationEpsilon);
minimumError
= dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-2,
1e2,
forceMagnitudeEpsilon);
translationalOptimizationEpsilon = 1e-2 * global.desiredMaxDisplacementValue;
wasSuccessful = minimumError < translationalOptimizationEpsilon;
break;
case Dome:
global.desiredMaxRotationAngle = vcg::math::ToRad(35.0);
global.desiredMaxRotationAngle = vcg::math::ToRad(25.0);
global.constructScenarioFunction = &ReducedModelOptimizer::constructDomeSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[1].first;
dlib::find_min_single_variable(
&fullPatternMaxSimulationForceRotationalObjective,
forceMagnitude,
1e-8,
1e8,
vcg::math::ToRad(1.0)
// global.desiredMaxRotationAngle * 0.5,
// optimizationEpsilon,
);
// minimumError
// = dlib::find_min_single_variable(&fullPatternMaxSimulationForceRotationalObjective,
// forceMagnitude,
// 1e-2,
// 1e2,
// forceMagnitudeEpsilon);
// wasSuccessful = minimumError < vcg::math::ToRad(3.0);
result = dlib::find_min_global(&fullPatternMaxSimulationForceRotationalObjective,
1e-2,
1,
dlib::max_function_calls(50));
wasSuccessful = result.y < vcg::math::ToRad(3.0);
break;
case Saddle:
// global.desiredMaxDisplacementValue *= 2;
global.desiredMaxDisplacementValue = 0.1
* vcg::Distance(global.baseTriangle.cP(0),
(global.baseTriangle.cP(1)
@ -1125,14 +1156,30 @@ double ReducedModelOptimizer::getFullPatternMaxSimulationForce(const BaseSimulat
/ 2);
global.constructScenarioFunction = &ReducedModelOptimizer::constructSaddleSimulationScenario;
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first;
dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-8,
1e8,
1e-2);
minimumError
= dlib::find_min_single_variable(&fullPatternMaxSimulationForceTranslationalObjective,
forceMagnitude,
1e-2,
1e2,
forceMagnitudeEpsilon);
translationalOptimizationEpsilon = 1e-2 * global.desiredMaxDisplacementValue;
wasSuccessful = minimumError < translationalOptimizationEpsilon;
break;
}
#ifdef POLYSCOPE_DEFINED
std::cout << "Max " << ReducedPatternOptimization::baseSimulationScenarioNames[scenario]
<< " magnitude:" << forceMagnitude << std::endl;
if (!wasSuccessful)
std::cout << "Was not successfull" << std::endl;
#endif
if (!wasSuccessful) {
const std::string debugMessage
= ReducedPatternOptimization::baseSimulationScenarioNames[scenario]
+ " max scenario magnitude was not succefully determined.";
optimizationNotes.append(debugMessage);
}
return forceMagnitude;
}
@ -1148,11 +1195,9 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
job.pMesh = pMesh;
//// Axial
#ifdef POLYSCOPE_DEFINED
std::cout << "Computing Axial scenarios.." << std::endl;
#endif
const double maxForceMagnitude_axial = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Axial);
bool wasSuccessful;
const double maxForceMagnitude_axial
= getFullPatternMaxSimulationForce(BaseSimulationScenario::Axial, wasSuccessful);
const double minForceMagnitude_axial = -maxForceMagnitude_axial;
const int numberOfSimulationScenarios_axial
= simulationScenariosResolution[BaseSimulationScenario::Axial];
@ -1167,9 +1212,6 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
for (int axialSimulationScenarioIndex = 0;
axialSimulationScenarioIndex < numberOfSimulationScenarios_axial;
axialSimulationScenarioIndex++) {
#ifdef POLYSCOPE_DEFINED
std::cout << "Computing Axial scenario:" << axialSimulationScenarioIndex << std::endl;
#endif
job.nodalExternalForces.clear();
job.constrainedVertices.clear();
job.nodalForcedDisplacements.clear();
@ -1185,8 +1227,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
}
//// Shear
const double maxForceMagnitude_shear = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Shear);
const double maxForceMagnitude_shear
= getFullPatternMaxSimulationForce(BaseSimulationScenario::Shear, wasSuccessful);
const double minForceMagnitude_shear = -maxForceMagnitude_shear;
const int numberOfSimulationScenarios_shear
= simulationScenariosResolution[BaseSimulationScenario::Shear];
@ -1201,9 +1243,6 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
for (int shearSimulationScenarioIndex = 0;
shearSimulationScenarioIndex < numberOfSimulationScenarios_shear;
shearSimulationScenarioIndex++) {
#ifdef POLYSCOPE_DEFINED
std::cout << "Computing shear scenario:" << shearSimulationScenarioIndex << std::endl;
#endif
job.constrainedVertices.clear();
job.nodalExternalForces.clear();
job.label = baseSimulationScenarioNames[BaseSimulationScenario::Shear] + "_"
@ -1216,8 +1255,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
}
//// Bending
const double maxForceMagnitude_bending = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Bending);
const double maxForceMagnitude_bending
= getFullPatternMaxSimulationForce(BaseSimulationScenario::Bending, wasSuccessful);
const double minForceMagnitude_bending = 0;
const int numberOfSimulationScenarios_bending
= simulationScenariosResolution[BaseSimulationScenario::Bending];
@ -1230,9 +1269,6 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
for (int bendingSimulationScenarioIndex = 0;
bendingSimulationScenarioIndex < numberOfSimulationScenarios_bending;
bendingSimulationScenarioIndex++) {
#ifdef POLYSCOPE_DEFINED
std::cout << "Computing bending scenario:" << bendingSimulationScenarioIndex << std::endl;
#endif
job.nodalExternalForces.clear();
job.constrainedVertices.clear();
job.label = baseSimulationScenarioNames[BaseSimulationScenario::Bending] + "_"
@ -1248,8 +1284,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
}
//// Dome
const double maxForceMagnitude_dome = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Dome);
const double maxForceMagnitude_dome
= getFullPatternMaxSimulationForce(BaseSimulationScenario::Dome, wasSuccessful);
const double minForceMagnitude_dome = 0;
const int numberOfSimulationScenarios_dome
= simulationScenariosResolution[BaseSimulationScenario::Dome];
@ -1262,9 +1298,6 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
for (int domeSimulationScenarioIndex = 0;
domeSimulationScenarioIndex < numberOfSimulationScenarios_dome;
domeSimulationScenarioIndex++) {
#ifdef POLYSCOPE_DEFINED
std::cout << "Computing dome scenario:" << domeSimulationScenarioIndex << std::endl;
#endif
job.constrainedVertices.clear();
job.nodalExternalForces.clear();
job.nodalForcedDisplacements.clear();
@ -1278,8 +1311,8 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
}
//// Saddle
const double maxForceMagnitude_saddle = getFullPatternMaxSimulationForce(
BaseSimulationScenario::Saddle);
const double maxForceMagnitude_saddle
= getFullPatternMaxSimulationForce(BaseSimulationScenario::Saddle, wasSuccessful);
const double minForceMagnitude_saddle = 0;
const int numberOfSimulationScenarios_saddle
= simulationScenariosResolution[BaseSimulationScenario::Saddle];
@ -1292,9 +1325,6 @@ ReducedModelOptimizer::createFullPatternSimulationScenarios(
for (int saddleSimulationScenarioIndex = 0;
saddleSimulationScenarioIndex < numberOfSimulationScenarios_saddle;
saddleSimulationScenarioIndex++) {
#ifdef POLYSCOPE_DEFINED
std::cout << "Computing saddle scenario:" << saddleSimulationScenarioIndex << std::endl;
#endif
job.constrainedVertices.clear();
job.nodalExternalForces.clear();
job.nodalForcedDisplacements.clear();
@ -1426,10 +1456,15 @@ void ReducedModelOptimizer::optimize(
m_pFullPatternSimulationMesh);
// polyscope::removeAllStructures();
results.baseTriangle = global.baseTriangle;
DRMSimulationModel::Settings simulationSettings;
simulationSettings.shouldDraw = false;
simulationSettings.maxDRMIterations = 200000;
// simulationSettings.shouldDraw = true;
// simulationSettings.isDebugMode = true;
// simulationSettings.debugModeStep = 100000;
#ifdef POLYSCOPE_DEFINED
const bool drawFullPatternSimulationResults = true;
const bool drawFullPatternSimulationResults = false;
if (drawFullPatternSimulationResults) {
global.fullPatternSimulationJobs[0]->pMesh->registerForDrawing(
ReducedPatternOptimization::Colors::fullInitial);
@ -1439,8 +1474,14 @@ void ReducedModelOptimizer::optimize(
for (int simulationScenarioIndex : global.simulationScenarioIndices) {
const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob
= global.fullPatternSimulationJobs[simulationScenarioIndex];
SimulationResults fullPatternResults = simulator.executeSimulation(pFullPatternSimulationJob,
simulationSettings);
if (!fullPatternResults.converged) {
results.wasSuccessful = false;
return;
}
results.wasSuccessful = true;
#ifdef POLYSCOPE_DEFINED
if (drawFullPatternSimulationResults) {
// SimulationResults fullPatternResults_linear = linearSimulator.executeSimulation(
@ -1481,4 +1522,5 @@ void ReducedModelOptimizer::optimize(
std::cout << "Running reduced model optimization" << std::endl;
#endif
runOptimization(optimizationSettings, results);
results.notes = optimizationNotes;
}

3
src/reducedmodeloptimizer.hpp
View File

@ -28,6 +28,7 @@ class ReducedModelOptimizer
m_fullPatternOppositeInterfaceViPairs;
std::unordered_map<size_t, size_t> nodeToSlot;
std::unordered_map<size_t, std::unordered_set<size_t>> slotToNode;
std::string optimizationNotes;
public:
constexpr static std::array<int, 5> simulationScenariosResolution = {10, 10, 20, 20, 20};
@ -186,7 +187,7 @@ private:
ReducedPatternOptimization::Results &results);
std::vector<double> getFullPatternMaxSimulationForces();
double getFullPatternMaxSimulationForce(
const ReducedPatternOptimization::BaseSimulationScenario &scenario);
const ReducedPatternOptimization::BaseSimulationScenario &scenario, bool &wasSuccessfull);
};
void updateMesh(long n, const double *x);
#endif // REDUCEDMODELOPTIMIZER_HPP