Normalizing based on a globally defined epsilon(0.0003)

This commit is contained in:
Iason 2021-02-24 22:10:42 +02:00
parent 94eeb0a9df
commit 55731d3131
3 changed files with 82 additions and 131 deletions

View File

@ -52,7 +52,7 @@ int main(int argc, char *argv[]) {
settings_optimization.numberOfFunctionCalls =
input_numberOfFunctionCallsDefined ? std::atoi(argv[3]) : 100;
settings_optimization.shouldNormalizeObjectiveValue = true;
settings_optimization.solutionAccuracy = 0.1;
settings_optimization.solutionAccuracy = 0.01;
// Optimize pair
const std::string pairName =

View File

@ -9,7 +9,7 @@
struct GlobalOptimizationVariables {
std::vector<Eigen::MatrixX3d> g_optimalReducedModelDisplacements;
std::vector<std::vector<Vector6d>> fullPatternDisplacements;
std::vector<double> reducedPatternMaximumDisplacementNormSum;
std::vector<double> fullPatternDisplacementNormSum;
std::vector<std::shared_ptr<SimulationJob>> fullPatternSimulationJobs;
std::vector<std::shared_ptr<SimulationJob>> reducedPatternSimulationJobs;
std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
@ -47,13 +47,7 @@ double ReducedModelOptimizer::computeError(
computeRawError(reducedPatternDisplacements, fullPatternDisplacements,
reducedToFullInterfaceViMap);
if (global.optimizationSettings.shouldNormalizeObjectiveValue) {
assert(rawError <= normalizationFactor);
if (rawError > normalizationFactor) {
std::cout << "Normalized error is bigger than one." << std::endl;
}
return rawError /
normalizationFactor; // std::max(interfaceDisplacementsNormSum,
// 0.00003);
return rawError / std::max(normalizationFactor, 0.0003);
}
return rawError;
@ -199,36 +193,11 @@ double ReducedModelOptimizer::objective(long n, const double *x) {
std::cout << "Failed simulation" << std::endl;
}
} else {
const double thisSimulationScenarioError =
computeError(reducedModelResults.displacements,
global.fullPatternDisplacements[simulationScenarioIndex],
global.reducedToFullInterfaceViMap,
global.reducedPatternMaximumDisplacementNormSum
[simulationScenarioIndex]);
if (thisSimulationScenarioError > 1) {
std::cout << "Simulation scenario "
<< simulationScenarioStrings[simulationScenarioIndex]
<< " results in an error bigger than one." << std::endl;
for (size_t parameterIndex = 0; parameterIndex < n; parameterIndex++) {
std::cout << "x[" + std::to_string(parameterIndex) + "]="
<< x[parameterIndex] << std::endl;
}
#ifdef POLYSCOPE_DEFINED
ReducedModelOptimizer::visualizeResults(
global.fullPatternSimulationJobs[simulationScenarioIndex],
global.reducedPatternSimulationJobs[simulationScenarioIndex],
global.reducedToFullInterfaceViMap, false);
ReducedModelOptimizer::visualizeResults(
global.fullPatternSimulationJobs[simulationScenarioIndex],
std::make_shared<SimulationJob>(
reducedPatternMaximumDisplacementSimulationJobs
[simulationScenarioIndex]),
global.reducedToFullInterfaceViMap, true);
polyscope::removeAllStructures();
#endif // POLYSCOPE_DEFINED
}
const double thisSimulationScenarioError = computeError(
reducedModelResults.displacements,
global.fullPatternDisplacements[simulationScenarioIndex],
global.reducedToFullInterfaceViMap,
global.fullPatternDisplacementNormSum[simulationScenarioIndex]);
error += thisSimulationScenarioError;
}
// std::ofstream out(filename, std::ios_base::app);
@ -269,7 +238,7 @@ double ReducedModelOptimizer::objective(long n, const double *x) {
//}
// compute error and return it
global.gObjectiveValueHistory.push_back(error);
// global.gObjectiveValueHistory.push_back(error);
// auto xPlot = matplot::linspace(0, gObjectiveValueHistory.size(),
// gObjectiveValueHistory.size());
// std::vector<double> colors(gObjectiveValueHistory.size(), 2);
@ -577,7 +546,7 @@ void ReducedModelOptimizer::computeReducedModelSimulationJob(
}
#if POLYSCOPE_DEFINED
void ReducedModelOptimizer::visualizeResults(
void ReducedModelOptimizer::registerResultsForDrawing(
const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob,
const std::shared_ptr<SimulationJob> &pReducedPatternSimulationJob,
const std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
@ -609,19 +578,6 @@ void ReducedModelOptimizer::visualizeResults(
global.reducedPatternMaximumDisplacementNormSum[simulationScenarioIndex]*/);
std::cout << "Error is " << error << std::endl;
reducedModelResults.registerForDrawing();
if (draw) {
polyscope::show();
const std::string screenshotFilename =
"/home/iason/Coding/Projects/Approximating shapes with flat "
"patterns/RodModelOptimizationForPatterns/build/OptimizationResults/"
"Images/" +
pFullPatternSimulationMesh->getLabel() + "_" + "noScenarioName";
// simulationScenarioStrings[simulationScenarioIndex];
polyscope::screenshot(screenshotFilename, false);
fullModelResults.unregister();
reducedModelResults.unregister();
pFullPatternSimulationMesh->unregister();
}
}
#endif // POLYSCOPE_DEFINED
@ -720,16 +676,37 @@ ReducedModelOptimizer::runOptimization(const Settings &settings) {
results.numberOfSimulationCrashes = global.numOfSimulationCrashes;
results.x = global.minX;
results.objectiveValue = global.minY;
if (settings.shouldNormalizeObjectiveValue) {
global.optimizationSettings.shouldNormalizeObjectiveValue = false;
if (global.optimizeInnerHexagonSize) {
results.objectiveValue = objective(4, results.x.data());
// if (settings.shouldNormalizeObjectiveValue) {
// global.optimizationSettings.shouldNormalizeObjectiveValue = false;
// if (global.optimizeInnerHexagonSize) {
// results.objectiveValue = objective(4, results.x.data());
} else {
results.objectiveValue =
objective(results.x[0], results.x[1], results.x[2]);
}
global.optimizationSettings.shouldNormalizeObjectiveValue = true;
// } else {
// results.objectiveValue =
// objective(results.x[0], results.x[1], results.x[2]);
// }
// global.optimizationSettings.shouldNormalizeObjectiveValue = true;
// }
// Compute obj value per simulation scenario
updateMesh(results.x.size(), results.x.data());
results.objectiveValuePerSimulationScenario.resize(
global.simulationScenarioIndices.size());
FormFinder::Settings simulationSettings;
FormFinder simulator;
for (int simulationScenarioIndex = 0;
simulationScenarioIndex <
SimulationScenario::NumberOfSimulationScenarios;
simulationScenarioIndex++) {
SimulationResults reducedModelResults = simulator.executeSimulation(
global.reducedPatternSimulationJobs[simulationScenarioIndex],
simulationSettings);
const double error = computeError(
reducedModelResults.displacements,
global.fullPatternDisplacements[simulationScenarioIndex],
global.reducedToFullInterfaceViMap,
global.fullPatternDisplacementNormSum[simulationScenarioIndex]);
results.objectiveValuePerSimulationScenario[simulationScenarioIndex] =
error;
}
results.time = elapsed.count() / 1000.0;
const bool printDebugInfo = false;
@ -1044,7 +1021,6 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
global.g_optimalReducedModelDisplacements.resize(6);
global.reducedPatternSimulationJobs.resize(6);
global.fullPatternDisplacements.resize(6);
global.reducedPatternMaximumDisplacementNormSum.resize(6);
global.g_firstRoundIterationIndex = 0;
global.minY = std::numeric_limits<double>::max();
global.numOfSimulationCrashes = 0;
@ -1052,7 +1028,7 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
global.optimizationSettings = optimizationSettings;
global.fullPatternSimulationJobs =
createScenarios(m_pFullPatternSimulationMesh);
reducedPatternMaximumDisplacementSimulationJobs.resize(6);
global.fullPatternDisplacementNormSum.resize(6);
// polyscope::removeAllStructures();
FormFinder::Settings simulationSettings;
@ -1073,51 +1049,17 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
std::make_shared<SimulationJob>(reducedPatternSimulationJob);
if (optimizationSettings.shouldNormalizeObjectiveValue) {
// Compute the quantities for normalizing the obj value
const double minB = optimizationSettings.xRanges[0].min;
const double maxRatio = optimizationSettings.xRanges[1].max;
const double minE = optimizationSettings.xRanges[2].min;
const double minHS = 0.3;
std::vector<double> mostFlexibleOptimizationParameters{minB, maxRatio,
minE, minHS};
if (global.optimizeInnerHexagonSize) {
updateMesh(4, mostFlexibleOptimizationParameters.data());
} else {
updateMesh(3, mostFlexibleOptimizationParameters.data());
double displacementNormSum = 0;
for (auto interfaceViPair : global.reducedToFullInterfaceViMap) {
const int interfaceFullPatternVi = interfaceViPair.second;
const Eigen::Vector3d displacement(
fullModelResults.displacements[interfaceFullPatternVi][0],
fullModelResults.displacements[interfaceFullPatternVi][1],
fullModelResults.displacements[interfaceFullPatternVi][2]);
displacementNormSum += displacement.norm();
}
reducedPatternMaximumDisplacementSimulationJobs[simulationScenarioIndex] =
global.reducedPatternSimulationJobs[simulationScenarioIndex]
->getCopy();
reducedPatternMaximumDisplacementSimulationJobs[simulationScenarioIndex]
.pMesh->setLabel("reduced_maxDisplacement");
SimulationResults reducedModelResults = simulator.executeSimulation(
global.reducedPatternSimulationJobs[simulationScenarioIndex],
simulationSettings);
const double errorOfMaxDisplacedReduced = computeRawError(
reducedModelResults.displacements, fullModelResults.displacements,
global.reducedToFullInterfaceViMap);
const double errorOfNonDisplacedReduced = computeRawError(
std::vector<Vector6d>(reducedModelResults.displacements.size(),
Vector6d(0)),
fullModelResults.displacements, global.reducedToFullInterfaceViMap);
double displacementMultiplier = 2;
// if (global.optimizeInnerHexagonSize) {
// displacementMultiplier = 2;
// } else {
// displacementMultiplier = 2;
// }
global.reducedPatternMaximumDisplacementNormSum[simulationScenarioIndex] =
displacementMultiplier *
std::max(errorOfMaxDisplacedReduced, errorOfNonDisplacedReduced);
// if (errorOfMaxDisplacedReduced > errorOfNonDisplacedReduced) {
// std::cout << "Max disp results in a bigger error for scenario "
// << simulationScenarioStrings[simulationScenarioIndex]
// << std::endl;
// } else {
// std::cout << "Zero disp results in a bigger error for scenario "
// << simulationScenarioStrings[simulationScenarioIndex]
// << std::endl;
// }
global.fullPatternDisplacementNormSum[simulationScenarioIndex] =
displacementNormSum;
}
}
Results optResults = runOptimization(optimizationSettings);
@ -1132,9 +1074,7 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
// visualizeResults(simulationJobs, global.simulationScenarioIndices);
#ifdef POLYSCOPE_DEFINED
// visualizeResults(
// global.fullPatternSimulationJobs, global.reducedPatternSimulationJobs,
// global.simulationScenarioIndices, global.reducedToFullInterfaceViMap);
optResults.draw();
#endif // POLYSCOPE_DEFINED
return optResults;
}

View File

@ -95,6 +95,7 @@ public:
int numberOfSimulationCrashes{0};
std::vector<double> x;
double objectiveValue;
std::vector<double> objectiveValuePerSimulationScenario;
std::vector<std::shared_ptr<SimulationJob>> fullPatternSimulationJobs;
std::vector<std::shared_ptr<SimulationJob>> reducedPatternSimulationJobs;
@ -189,21 +190,22 @@ public:
reducedModelResults.registerForDrawing();
polyscope::show();
// Save a screensh
// const std::string screenshotFilename =
// "/home/iason/Coding/Projects/Approximating shapes with flat "
// "patterns/RodModelOptimizationForPatterns/build/OptimizationResults/"
// + m_pFullPatternSimulationMesh->getLabel() + "_" +
// simulationScenarioStrings[simulationScenarioIndex];
// polyscope::screenshot(screenshotFilename, false);
const std::string screenshotFilename =
"/home/iason/Coding/Projects/Approximating shapes with flat "
"patterns/RodModelOptimizationForPatterns/Results/Images/" +
pFullPatternSimulationJob->pMesh->getLabel() + "_" +
simulationScenarioStrings[simulationJobIndex];
polyscope::screenshot(screenshotFilename, false);
fullModelResults.unregister();
reducedModelResults.unregister();
// double error = computeError(
// reducedModelResults,
// global.g_optimalReducedModelDisplacements[simulationScenarioIndex]);
// std::cout << "Error of simulation scenario "
// << simulationScenarioStrings[simulationScenarioIndex] <<
// " is "
// << error << std::endl;
// double error = computeError(
// reducedModelResults.displacements,fullModelResults.displacements,
// );
// std::cout << "Error of simulation scenario "
// <<
// simulationScenarioStrings[simulationScenarioIndex]
// << " is "
// << error << std::endl;
}
}
#endif // POLYSCOPE_DEFINED
@ -239,7 +241,13 @@ public:
void
writeHeaderTo(const ReducedModelOptimizer::Settings &settings_optimization,
csvFile &os) {
os << "Obj value";
os << "Total Obj value";
for (int simulationScenarioIndex = 0;
simulationScenarioIndex <
SimulationScenario::NumberOfSimulationScenarios;
simulationScenarioIndex++) {
os << "Obj value " + simulationScenarioStrings[simulationScenarioIndex];
}
for (const ReducedModelOptimizer::xRange &range :
settings_optimization.xRanges) {
os << range.label;
@ -252,6 +260,9 @@ public:
writeResultsTo(const ReducedModelOptimizer::Settings &settings_optimization,
csvFile &os) const {
os << objectiveValue;
for (double scenarioObjValue : objectiveValuePerSimulationScenario) {
os << scenarioObjValue;
}
for (const double &optimalX : x) {
os << optimalX;
}
@ -319,7 +330,7 @@ public:
&fullToReducedInterfaceViMap,
std::unordered_map<FullPatternVertexIndex, ReducedPatternVertexIndex>
&fullPatternOppositeInterfaceViMap);
static void visualizeResults(
static void registerResultsForDrawing(
const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob,
const std::shared_ptr<SimulationJob> &pReducedPatternSimulationJob,
const std::unordered_map<ReducedPatternVertexIndex,
@ -343,10 +354,10 @@ public:
const double &normalizationFactor);
private:
void
visualizeResults(const std::vector<std::shared_ptr<SimulationJob>>
&fullPatternSimulationJobs,
const std::vector<SimulationScenario> &simulationScenarios);
void registerResultsForDrawing(
const std::vector<std::shared_ptr<SimulationJob>>
&fullPatternSimulationJobs,
const std::vector<SimulationScenario> &simulationScenarios);
static void computeDesiredReducedModelDisplacements(
const SimulationResults &fullModelResults,
const std::unordered_map<size_t, size_t> &displacementsReducedToFullMap,