Normalizing per simulation scenario using the max(sum(fullPatternDisplacementNorm,d_e))

src/main.cpp

@@ -72,15 +72,15 @@ int main(int argc, char *argv[]) {
     FlatPattern fullPattern(filepathString);
       fullPattern.setLabel(filepath.stem().string());
       fullPattern.scale(0.03);
-    for (int reducedPatternIndex = 0;
+    //for (int reducedPatternIndex = 0;
-         reducedPatternIndex < reducedModels.size(); reducedPatternIndex++) {
+    //     reducedPatternIndex < reducedModels.size(); reducedPatternIndex++) {
       FlatPattern *pFullPattern = new FlatPattern();
       pFullPattern->copy(fullPattern);
       FlatPattern *pReducedPattern = new FlatPattern();
-      pReducedPattern->copy(*reducedModels[reducedPatternIndex]);
+      //pReducedPattern->copy(*reducedModels[reducedPatternIndex]);
-            //pReducedPattern->copy(*reducedModels[0]);
+            pReducedPattern->copy(*reducedModels[0]);
       patternPairs.push_back(std::make_pair(pFullPattern, pReducedPattern));
-    }
+    //}
   }
 
   //  for (double rangeOffset = 0.15; rangeOffset <= 0.95; rangeOffset += 0.05)
@@ -92,7 +92,7 @@ int main(int argc, char *argv[]) {
   //  for (settings_optimization.numberOfFunctionCalls = 100;
   //       settings_optimization.numberOfFunctionCalls < 5000;
   //       settings_optimization.numberOfFunctionCalls += 100) {
-  settings_optimization.numberOfFunctionCalls = 10;
+  settings_optimization.numberOfFunctionCalls = 100;
   const std::string optimizationSettingsString =
       settings_optimization.toString();
   std::string optimiziationResultsDirectory = "../OptimizationResults";
@@ -116,7 +116,7 @@ int main(int argc, char *argv[]) {
       patternPairs.size());
   auto start = std::chrono::high_resolution_clock::now();
 
-  #pragma omp parallel for
+  //#pragma omp parallel for
   for (int patternPairIndex = 0; patternPairIndex < patternPairs.size();
        patternPairIndex++) {
     const std::vector<size_t> numberOfNodesPerSlot{1, 0, 0, 2, 1, 2, 1};
@@ -124,12 +124,13 @@ int main(int argc, char *argv[]) {
     optimizer.initializePatterns(*patternPairs[patternPairIndex].first,
                                  *patternPairs[patternPairIndex].second, {});
     ReducedModelOptimizer::Results optimizationResults =
-        optimizer.optimize(settings_optimization);
+        optimizer.optimize(settings_optimization, {ReducedModelOptimizer::SimulationScenario::Dome});
 
     totalError += optimizationResults.objectiveValue;
     optimizationResults_testSet[patternPairIndex] = optimizationResults;
     totalNumberOfSimulationCrashes +=
         optimizationResults.numberOfSimulationCrashes;
+    optimizationResults.draw();
   }
   auto end = std::chrono::high_resolution_clock::now();
   auto runtime_ms =
@@ -144,7 +145,6 @@ int main(int argc, char *argv[]) {
     std::filesystem::create_directories(std::filesystem::path(saveToPath));
 
     optimizationResults_testSet[patternPairIndex].save(saveToPath.string());
-    optimizationResults_testSet[patternPairIndex].draw();
   }
   csvFile statistics(std::filesystem::path(thisOptimizationDirectory)
                          .append("statistics.csv")

src/reducedmodeloptimizer.cpp

@@ -9,12 +9,12 @@
 
 struct GlobalOptimizationVariables {
   std::vector<Eigen::MatrixX3d> g_optimalReducedModelDisplacements;
+  std::vector<std::vector<Vector6d>> fullPatternDisplacements;
+  std::vector<double> fullPatternDisplacementNormSum;
   std::vector<SimulationJob> g_fullPatternSimulationJob;
   std::vector<std::shared_ptr<SimulationJob>> reducedPatternSimulationJobs;
   std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
       reducedToFullInterfaceViMap;
-  std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
-      g_reducedToFullViMap;
   matplot::line_handle gPlotHandle;
   std::vector<double> gObjectiveValueHistory;
   Eigen::Vector2d g_initialX;
@@ -122,39 +122,9 @@ GlobalOptimizationVariables tls[MAX_THREAD] __attribute__((aligned(64)));
 //  }
 //};
 
-double ReducedModelOptimizer::computeError(
-    const SimulationResults &reducedPatternResults,
-    const Eigen::MatrixX3d &optimalReducedPatternDisplacements) {
-  double error = 0;
-  auto &global = tls[omp_get_thread_num()];
-  for (const auto reducedFullViPair : global.reducedToFullInterfaceViMap) {
-    VertexIndex reducedModelVi = reducedFullViPair.first;
-    //    const auto pos =
-    //        g_reducedPatternSimulationJob.mesh->vert[reducedModelVi].cP();
-    //    std::cout << "Interface vi " << reducedModelVi << " is at position "
-    //              << pos[0] << " " << pos[1] << " " << pos[2] << std::endl;
-    Eigen::Vector3d vertexDisplacement(
-        reducedPatternResults.displacements[reducedModelVi][0],
-        reducedPatternResults.displacements[reducedModelVi][1],
-        reducedPatternResults.displacements[reducedModelVi][2]);
-    if (!std::isfinite(vertexDisplacement[0]) ||
-        !std::isfinite(vertexDisplacement[1]) ||
-        !std::isfinite(vertexDisplacement[2])) {
-      return std::numeric_limits<double>::max();
-    }
-    Eigen::Vector3d errorVector =
-        Eigen::Vector3d(
-            optimalReducedPatternDisplacements.row(reducedModelVi)) -
-        vertexDisplacement;
-    //    error += errorVector.squaredNorm();
-    error += errorVector.norm();
-  }
-  return error;
-}
-
 double ReducedModelOptimizer::computeError(
     const std::vector<Vector6d> &reducedPatternDisplacements,
-    const std::vector<Vector6d> &fullPatternDisplacements,
+    const std::vector<Vector6d> &fullPatternDisplacements,const double& interfaceDisplacementsNormSum,
     const std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
         &reducedToFullInterfaceViMap) {
   double error = 0;
@@ -182,7 +152,9 @@ double ReducedModelOptimizer::computeError(
     //    error += errorVector.squaredNorm();
     error += errorVector.norm();
   }
-  return error;
+
+  return error/std::max(interfaceDisplacementsNormSum,0.00003);
+  //return error;
 }
 
 void updateMesh(long n, const double *x) {
@@ -317,9 +289,8 @@ double ReducedModelOptimizer::objective(long n, const double *x) {
       //          simulationSettings);
       global.numOfSimulationCrashes++;
     } else {
-      error += computeError(
+        error += computeError(
-          reducedModelResults,
+            reducedModelResults.displacements, global.fullPatternDisplacements[simulationScenarioIndex],global.fullPatternDisplacementNormSum[simulationScenarioIndex], global.reducedToFullInterfaceViMap);
-          global.g_optimalReducedModelDisplacements[simulationScenarioIndex]);
       filename = "/home/iason/Coding/Projects/Approximating shapes with flat "
                  "patterns/RodModelOptimizationForPatterns/build/"
                  "ProblematicSimulationJobs/conv_dimensions.txt";
@@ -1081,8 +1052,14 @@ void ReducedModelOptimizer::visualizeResults(
         reducedPatternSimulationJobs[simulationScenarioIndex];
     SimulationResults reducedModelResults =
         simulator.executeSimulation(pReducedPatternSimulationJob);
+    double interfaceDisplacementNormSum = 0;
+    for (const auto& interfaceViPair : reducedToFullInterfaceViMap) {
+        const int fullPatternInterfaceIndex = interfaceViPair.second;
+        Eigen::Vector3d fullPatternDisplacementVector(fullModelResults.displacements[fullPatternInterfaceIndex][0], fullModelResults.displacements[fullPatternInterfaceIndex][1], fullModelResults.displacements[fullPatternInterfaceIndex][2]);
+        interfaceDisplacementNormSum += fullPatternDisplacementVector.norm();
+    }
     double error = computeError(reducedModelResults.displacements,
-                                fullModelResults.displacements,
+                                fullModelResults.displacements,interfaceDisplacementNormSum,
                                 reducedToFullInterfaceViMap);
     std::cout << "Error of simulation scenario "
               << simulationScenarioStrings[simulationScenarioIndex] << " is "
@@ -1124,6 +1101,8 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
       createScenarios(m_pFullPatternSimulationMesh);
   global.g_optimalReducedModelDisplacements.resize(6);
   global.reducedPatternSimulationJobs.resize(6);
+  global.fullPatternDisplacements.resize(6);
+  global.fullPatternDisplacementNormSum.resize(6);
   global.g_firstRoundIterationIndex = 0;
   global.minY = std::numeric_limits<double>::max();
   global.numOfSimulationCrashes = 0;
@@ -1137,11 +1116,20 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
         simulationJobs[simulationScenarioIndex];
     SimulationResults fullModelResults =
         simulator.executeSimulation(pFullPatternSimulationJob, settings);
-    global.g_optimalReducedModelDisplacements[simulationScenarioIndex].resize(
+    global.fullPatternDisplacements[simulationScenarioIndex] = fullModelResults.displacements;
-        m_pReducedPatternSimulationMesh->VN(), 3);
+    double interfaceDisplacementNormSum = 0;
-    computeDesiredReducedModelDisplacements(
+    for (const auto& interfaceViPair : global.reducedToFullInterfaceViMap) {
-        fullModelResults, global.reducedToFullInterfaceViMap,
+        const int fullPatternInterfaceIndex = interfaceViPair.second;
-        global.g_optimalReducedModelDisplacements[simulationScenarioIndex]);
+        Eigen::Vector3d fullPatternDisplacementVector(fullModelResults.displacements[fullPatternInterfaceIndex][0], fullModelResults.displacements[fullPatternInterfaceIndex][1], fullModelResults.displacements[fullPatternInterfaceIndex][2]);
+        interfaceDisplacementNormSum += fullPatternDisplacementVector.norm();
+    }
+    global.fullPatternDisplacementNormSum[simulationScenarioIndex] = interfaceDisplacementNormSum;
+    //global.g_optimalReducedModelDisplacements[simulationScenarioIndex].resize(
+    //    m_pReducedPatternSimulationMesh->VN(), 3);
+    //computeDesiredReducedModelDisplacements(
+    //    fullModelResults, global.reducedToFullInterfaceViMap,
+    //    global.g_optimalReducedModelDisplacements[simulationScenarioIndex]);
+
     SimulationJob reducedPatternSimulationJob;
     reducedPatternSimulationJob.pMesh = m_pReducedPatternSimulationMesh;
     computeReducedModelSimulationJob(*pFullPatternSimulationJob,
@@ -1161,7 +1149,7 @@ ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
   //  optResults.draw();
 
   //  visualizeResults(simulationJobs, global.simulationScenarioIndices);
-  //  visualizeResults(simulationJobs, global.g_reducedPatternSimulationJob,
+    visualizeResults(simulationJobs, global.reducedPatternSimulationJobs,
-  //                   global.g_simulationScenarioIndices);
+                     global.simulationScenarioIndices,global.reducedToFullInterfaceViMap);
   return optResults;
 }

src/reducedmodeloptimizer.hpp

@@ -46,7 +46,8 @@ public:
   struct Settings {
     std::vector<xRange> xRanges;
     int numberOfFunctionCalls{100};
-    double solutionAccuracy{1e-5};
+    double solutionAccuracy{1e-2};
+    bool normalizeObjectiveValue{ true };
 
     std::string toString() const {
       std::string settingsString;
@@ -146,7 +147,7 @@ public:
                        &reducedToFullInterfaceViMap);
 
   static double computeError(const std::vector<Vector6d> &reducedPatternResults,
-                             const std::vector<Vector6d> &fullPatternResults,
+                             const std::vector<Vector6d> &fullPatternResults,const double& interfaceDisplacementNormSum,
                              const std::unordered_map<ReducedPatternVertexIndex,
                                                       FullPatternVertexIndex>
                                  &reducedToFullInterfaceViMap);