Refactored how the results of the optimization are saved.

src/main.cpp

@@ -3,9 +3,6 @@
 #include "edgemesh.hpp"
 #include "externvariables.hpp"
 #include "flatpattern.hpp"
-#include "polyscope/curve_network.h"
-#include "polyscope/point_cloud.h"
-#include "polyscope/polyscope.h"
 #include "reducedmodeloptimizer.hpp"
 #include "simulationhistoryplotter.hpp"
 #include "trianglepattterntopology.hpp"
@@ -57,33 +54,60 @@ int main(int argc, char *argv[]) {
   const bool input_numberOfFunctionCallsDefined = argc >= 5;
   settings_optimization.numberOfFunctionCalls =
       input_numberOfFunctionCallsDefined ? std::atoi(argv[4]) : 100;
+  settings_optimization.normalizeObjectiveValue = true;
+
   // Optimize pair
   const std::string pairName =
       fullPattern.getLabel() + "@" + reducedPattern.getLabel();
   if (printDebugInfo) {
     std::cout << "Optimizing " << pairName << std::endl;
   }
-
   const std::vector<size_t> numberOfNodesPerSlot{1, 0, 0, 2, 1, 2, 1};
   ReducedModelOptimizer optimizer(numberOfNodesPerSlot);
   optimizer.initializePatterns(fullPattern, reducedPattern, {});
   ReducedModelOptimizer::Results optimizationResults =
       optimizer.optimize(settings_optimization);
+
   // Export results
   const bool input_resultDirectoryDefined = argc >= 6;
-  std::string optimiziationResultsDirectory =
+  std::string optimizationResultsDirectory =
       input_resultDirectoryDefined ? argv[5] : "OptimizationResults";
-  std::filesystem::path dirPath_thisOptimization(
+  //// Get current date for creating the results folder
-      std::filesystem::path(optimiziationResultsDirectory).append(pairName));
+  std::time_t now = time(0);
-  std::filesystem::create_directories(dirPath_thisOptimization);
+  std::tm *ltm = std::localtime(&now);
+  std::string currentDate = std::to_string(ltm->tm_mday) + "_" +
+                            std::to_string(1 + ltm->tm_mon) + "_" +
+                            std::to_string(1900 + ltm->tm_year);
+  std::filesystem::path optimizationResultsDirectoryDatePath(
+      std::filesystem::path(optimizationResultsDirectory).append(currentDate));
+  if (optimizationResults.numberOfSimulationCrashes != 0) {
+    const auto crashedJobsDirPath =
+        std::filesystem::path(optimizationResultsDirectoryDatePath)
+            .append("CrashedJobs")
+            .append(pairName);
+    std::filesystem::create_directories(crashedJobsDirPath);
+    optimizationResults.save(crashedJobsDirPath.string());
+  } else {
+    std::filesystem::path dirPath_thisOptimization(
+        std::filesystem::path(optimizationResultsDirectoryDatePath)
+            .append("ConvergedJobs")
+            .append(pairName));
+    std::filesystem::create_directories(dirPath_thisOptimization);
+    optimizationResults.save(dirPath_thisOptimization.string());
+  }
   csvFile csv_results({}, false);
   //    csvFile csv_results(std::filesystem::path(dirPath_thisOptimization)
   //                            .append("results.csv")
   //                            .string(),
   //                        false);
-  settings_optimization.writeTo(csv_results);
+  csv_results << "Name";
-  optimizationResults.writeTo(settings_optimization, csv_results);
+  optimizationResults.writeHeaderTo(settings_optimization, csv_results);
-  optimizationResults.save(dirPath_thisOptimization.string());
+  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;
 
   //  optimizationResults.draw();
 

src/reducedmodeloptimizer.cpp

@@ -33,16 +33,7 @@ struct GlobalOptimizationVariables {
   int numberOfFunctionCalls{0};
   int numberOfOptimizationParameters{3};
   ReducedModelOptimizer::Settings optimizationSettings;
-};
+} global;
-
-// static GlobalOptimizationVariables global;
-
-const static int MAX_THREAD = 64;
-#if defined(_MSC_VER)
-__declspec(align(64)) GlobalOptimizationVariables tls[MAX_THREAD];
-#elif defined(__GNUC__)
-GlobalOptimizationVariables tls[MAX_THREAD] __attribute__((aligned(64)));
-#endif
 //#pragma omp threadprivate(global)
 
 // struct OptimizationCallback {
@@ -128,7 +119,6 @@ double ReducedModelOptimizer::computeError(
     const double &interfaceDisplacementsNormSum,
     const std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
         &reducedToFullInterfaceViMap) {
-  auto &global = tls[omp_get_thread_num()];
   double error = 0;
   for (const auto reducedFullViPair : reducedToFullInterfaceViMap) {
     VertexIndex reducedModelVi = reducedFullViPair.first;
@@ -162,7 +152,6 @@ double ReducedModelOptimizer::computeError(
 }
 
 void updateMesh(long n, const double *x) {
-  auto &global = tls[omp_get_thread_num()];
   std::shared_ptr<SimulationMesh> &pReducedPatternSimulationMesh =
       global.reducedPatternSimulationJobs[global.simulationScenarioIndices[0]]
           ->pMesh;
@@ -235,7 +224,6 @@ double ReducedModelOptimizer::objective(double b, double h, double E,
 }
 
 double ReducedModelOptimizer::objective(long n, const double *x) {
-  auto &global = tls[omp_get_thread_num()];
   //  std::cout.precision(17);
 
   //  for (size_t parameterIndex = 0; parameterIndex < n; parameterIndex++) {
@@ -263,34 +251,8 @@ double ReducedModelOptimizer::objective(long n, const double *x) {
         global.reducedPatternSimulationJobs[simulationScenarioIndex],
         simulationSettings);
     std::string filename;
-    if (!reducedModelResults.converged /*&&
+    if (!reducedModelResults.converged) {
-        g_reducedPatternSimulationJob[g_simulationScenarioIndices[0]]
-                    ->pMesh->elements[0]
-                    .A > 1e-8 &
-            x[0] / x[1] < 60*/) {
-      std::cout << "Failed simulation" << std::endl;
       error += std::numeric_limits<double>::max();
-      filename = "/home/iason/Coding/Projects/Approximating shapes with flat "
-                 "patterns/RodModelOptimizationForPatterns/build/"
-                 "ProblematicSimulationJobs/nonConv_dimensions.txt";
-      //      if (failedSimulationsXRatio.empty()) {
-      //        failedSimulationsXRatio.resize(2);
-      //      }
-      //      failedSimulationsXRatio[0].push_back(std::log(x[0] / x[1]));
-      //      failedSimulationsXRatio[1].push_back(
-      //          std::log(g_reducedPatternSimulationJob[g_simulationScenarioIndices[0]]
-      //                       ->pMesh->elements[0]
-      //                       .A));
-
-      //      SimulationResultsReporter::createPlot(
-      //          "log(b/h)", "log(A)", failedSimulationsXRatio[0],
-      //          failedSimulationsXRatio[1], "ratioToAPlot.png");
-      //      std::cout << "Failed simulation" << std::endl;
-      //      simulationSettings.shouldDraw = true;
-      //      simulationSettings.debugMessages = true;
-      //      simulator.executeSimulation(
-      //          g_reducedPatternSimulationJob[simulationScenarioIndex],
-      //          simulationSettings);
       global.numOfSimulationCrashes++;
     } else {
       error += computeError(
@@ -298,9 +260,6 @@ double ReducedModelOptimizer::objective(long n, const double *x) {
           global.fullPatternDisplacements[simulationScenarioIndex],
           global.fullPatternDisplacementNormSum[simulationScenarioIndex],
           global.reducedToFullInterfaceViMap);
-      filename = "/home/iason/Coding/Projects/Approximating shapes with flat "
-                 "patterns/RodModelOptimizationForPatterns/build/"
-                 "ProblematicSimulationJobs/conv_dimensions.txt";
     }
     std::ofstream out(filename, std::ios_base::app);
     auto pMesh =
@@ -391,7 +350,6 @@ void ReducedModelOptimizer::computeMaps(
     std::unordered_map<FullPatternVertexIndex, ReducedPatternVertexIndex>
         &fullPatternOppositeInterfaceViMap) {
 
-  auto &global = tls[omp_get_thread_num()];
   // Compute the offset between the interface nodes
   const size_t interfaceSlotIndex = 4; // bottom edge
   assert(slotToNode.find(interfaceSlotIndex) != slotToNode.end() &&
@@ -456,65 +414,11 @@ void ReducedModelOptimizer::computeMaps(
     assert(vi0 < fullPattern.VN() && vi1 < fullPattern.VN());
     fullPatternOppositeInterfaceViMap[vi0] = vi1;
   }
-
-  const bool debugMapping = false;
-  if (debugMapping) {
-    reducedPattern.registerForDrawing();
-    std::vector<glm::vec3> colors_reducedPatternExcludedEdges(
-        reducedPattern.EN(), glm::vec3(0, 0, 0));
-    for (const size_t ei : global.reducedPatternExludedEdges) {
-      colors_reducedPatternExcludedEdges[ei] = glm::vec3(1, 0, 0);
-    }
-    const std::string label = reducedPattern.getLabel();
-    polyscope::getCurveNetwork(label)
-        ->addEdgeColorQuantity("Excluded edges",
-                               colors_reducedPatternExcludedEdges)
-        ->setEnabled(true);
-    polyscope::show();
-
-    std::vector<glm::vec3> nodeColorsOpposite(fullPattern.VN(),
-                                              glm::vec3(0, 0, 0));
-    for (const std::pair<size_t, size_t> oppositeVerts :
-         fullPatternOppositeInterfaceViMap) {
-      auto color = polyscope::getNextUniqueColor();
-      nodeColorsOpposite[oppositeVerts.first] = color;
-      nodeColorsOpposite[oppositeVerts.second] = color;
-    }
-    fullPattern.registerForDrawing();
-    polyscope::getCurveNetwork(fullPattern.getLabel())
-        ->addNodeColorQuantity("oppositeMap", nodeColorsOpposite)
-        ->setEnabled(true);
-    polyscope::show();
-
-    std::vector<glm::vec3> nodeColorsReducedToFull_reduced(reducedPattern.VN(),
-                                                           glm::vec3(0, 0, 0));
-    std::vector<glm::vec3> nodeColorsReducedToFull_full(fullPattern.VN(),
-                                                        glm::vec3(0, 0, 0));
-    for (size_t vi = 0; vi < reducedPattern.VN(); vi++) {
-      if (global.reducedToFullInterfaceViMap.contains(vi)) {
-
-        auto color = polyscope::getNextUniqueColor();
-        nodeColorsReducedToFull_reduced[vi] = color;
-        nodeColorsReducedToFull_full[global.reducedToFullInterfaceViMap[vi]] =
-            color;
-      }
-    }
-    polyscope::getCurveNetwork(reducedPattern.getLabel())
-        ->addNodeColorQuantity("reducedToFull_reduced",
-                               nodeColorsReducedToFull_reduced)
-        ->setEnabled(true);
-    polyscope::getCurveNetwork(fullPattern.getLabel())
-        ->addNodeColorQuantity("reducedToFull_full",
-                               nodeColorsReducedToFull_full)
-        ->setEnabled(true);
-    polyscope::show();
-  }
 }
 
 void ReducedModelOptimizer::computeMaps(
     FlatPattern &fullPattern, FlatPattern &reducedPattern,
     const std::unordered_set<size_t> &reducedModelExcludedEdges) {
-  auto &global = tls[omp_get_thread_num()];
   ReducedModelOptimizer::computeMaps(
       reducedModelExcludedEdges, slotToNode, fullPattern, reducedPattern,
       global.reducedToFullInterfaceViMap, m_fullToReducedInterfaceViMap,
@@ -534,13 +438,11 @@ void ReducedModelOptimizer::initializePatterns(
   // reducedPattern.setLabel("reduced_pattern_" + reducedPattern.getLabel());
   assert(fullPattern.VN() == reducedPattern.VN() &&
          fullPattern.EN() >= reducedPattern.EN());
-  polyscope::removeAllStructures();
   // Create copies of the input models
   FlatPattern copyFullPattern;
   FlatPattern copyReducedPattern;
   copyFullPattern.copy(fullPattern);
   copyReducedPattern.copy(reducedPattern);
-  auto &global = tls[omp_get_thread_num()];
   global.optimizeInnerHexagonSize = copyReducedPattern.EN() == 2;
   if (global.optimizeInnerHexagonSize) {
     const double h = copyReducedPattern.getBaseTriangleHeight();
@@ -571,7 +473,6 @@ void ReducedModelOptimizer::initializePatterns(
 
 void ReducedModelOptimizer::initializeOptimizationParameters(
     const std::shared_ptr<SimulationMesh> &mesh) {
-  auto &global = tls[omp_get_thread_num()];
   global.numberOfOptimizationParameters = 3;
   global.g_initialParameters.resize(
       global.optimizeInnerHexagonSize ? ++global.numberOfOptimizationParameters
@@ -666,7 +567,6 @@ void ReducedModelOptimizer::computeDesiredReducedModelDisplacements(
 
 ReducedModelOptimizer::Results
 ReducedModelOptimizer::runOptimization(const Settings &settings) {
-  auto &global = tls[omp_get_thread_num()];
 
   global.gObjectiveValueHistory.clear();
 
@@ -1035,71 +935,9 @@ ReducedModelOptimizer::createScenarios(
 //  reducedModelResults.unregister();
 //}
 
-void ReducedModelOptimizer::visualizeResults(
-    const std::vector<std::shared_ptr<SimulationJob>>
-        &fullPatternSimulationJobs,
-    const std::vector<std::shared_ptr<SimulationJob>>
-        &reducedPatternSimulationJobs,
-    const std::vector<SimulationScenario> &simulationScenarios,
-    const std::unordered_map<ReducedPatternVertexIndex, FullPatternVertexIndex>
-        &reducedToFullInterfaceViMap) {
-  FormFinder simulator;
-  std::shared_ptr<SimulationMesh> pFullPatternSimulationMesh =
-      fullPatternSimulationJobs[0]->pMesh;
-  pFullPatternSimulationMesh->registerForDrawing();
-  double totalError = 0;
-  for (const int simulationScenarioIndex : simulationScenarios) {
-    const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob =
-        fullPatternSimulationJobs[simulationScenarioIndex];
-    pFullPatternSimulationJob->registerForDrawing(
-        pFullPatternSimulationMesh->getLabel());
-    SimulationResults fullModelResults =
-        simulator.executeSimulation(pFullPatternSimulationJob);
-    fullModelResults.registerForDrawing();
-    // fullModelResults.saveDeformedModel();
-    const std::shared_ptr<SimulationJob> &pReducedPatternSimulationJob =
-        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,
-        interfaceDisplacementNormSum, reducedToFullInterfaceViMap);
-    std::cout << "Error of simulation scenario "
-              << simulationScenarioStrings[simulationScenarioIndex] << " is "
-              << error << std::endl;
-    totalError += error;
-    reducedModelResults.registerForDrawing();
-    //    firstOptimizationRoundResults[simulationScenarioIndex].registerForDrawing();
-    //    reducedModelResults.saveDeformedModel();
-    //    registerWorldAxes();
-    const std::string screenshotFilename =
-        "/home/iason/Coding/Projects/Approximating shapes with flat "
-        "patterns/RodModelOptimizationForPatterns/build/OptimizationResults/"
-        "Images/" +
-        pFullPatternSimulationMesh->getLabel() + "_" +
-        simulationScenarioStrings[simulationScenarioIndex];
-    polyscope::show();
-    polyscope::screenshot(screenshotFilename, false);
-    fullModelResults.unregister();
-    reducedModelResults.unregister();
-    //    firstOptimizationRoundResults[simulationScenarioIndex].unregister();
-  }
-  std::cout << "Total error:" << totalError << std::endl;
-}
-
 ReducedModelOptimizer::Results ReducedModelOptimizer::optimize(
     const Settings &optimizationSettings,
     const std::vector<SimulationScenario> &simulationScenarios) {
-  auto &global = tls[omp_get_thread_num()];
 
   global.simulationScenarioIndices = simulationScenarios;
   if (global.simulationScenarioIndices.empty()) {

src/reducedmodeloptimizer.hpp

@@ -65,21 +65,19 @@ public:
       return settingsString;
     }
 
-    void writeTo(csvFile &os, const bool writeHeader = true) const {
+    void writeHeaderTo(csvFile &os) const {
-      // Create settings csv header
+      if (!xRanges.empty()) {
-      if (writeHeader) {
+        for (const xRange &range : xRanges) {
-        if (!xRanges.empty()) {
+          os << range.label + " max";
-          for (const xRange &range : xRanges) {
+          os << range.label + " min";
-            os << range.label + " max";
-            os << range.label + " min";
-          }
         }
-        os << "Function Calls";
-        os << "Solution Accuracy";
-        //        os << std::endl;
-        os << endrow;
       }
+      os << "Function Calls";
+      os << "Solution Accuracy";
+      //        os << std::endl;
+    }
 
+    void writeSettingsTo(csvFile &os) const {
       if (!xRanges.empty()) {
         for (const xRange &range : xRanges) {
           os << range.max;
@@ -88,8 +86,6 @@ public:
       }
       os << numberOfFunctionCalls;
       os << solutionAccuracy;
-      //      os << std::endl;
-      os << endrow;
     }
   };
 
@@ -191,50 +187,6 @@ struct ReducedModelOptimizer::Results {
   std::vector<std::shared_ptr<SimulationJob>> fullPatternSimulationJobs;
   std::vector<std::shared_ptr<SimulationJob>> reducedPatternSimulationJobs;
 
-  void draw() const {
-    initPolyscope();
-    FormFinder simulator;
-    assert(fullPatternSimulationJobs.size() ==
-           reducedPatternSimulationJobs.size());
-    fullPatternSimulationJobs[0]->pMesh->registerForDrawing();
-    reducedPatternSimulationJobs[0]->pMesh->registerForDrawing();
-
-    const int numberOfSimulationJobs = fullPatternSimulationJobs.size();
-    for (int simulationJobIndex = 0;
-         simulationJobIndex < numberOfSimulationJobs; simulationJobIndex++) {
-      // Drawing of full pattern results
-      const std::shared_ptr<SimulationJob> &pFullPatternSimulationJob =
-          fullPatternSimulationJobs[simulationJobIndex];
-      pFullPatternSimulationJob->registerForDrawing(
-          fullPatternSimulationJobs[0]->pMesh->getLabel());
-      SimulationResults fullModelResults =
-          simulator.executeSimulation(pFullPatternSimulationJob);
-      fullModelResults.registerForDrawing();
-      // Drawing of reduced pattern results
-      const std::shared_ptr<SimulationJob> &pReducedPatternSimulationJob =
-          reducedPatternSimulationJobs[simulationJobIndex];
-      SimulationResults reducedModelResults =
-          simulator.executeSimulation(pReducedPatternSimulationJob);
-      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);
-      fullModelResults.unregister();
-      reducedModelResults.unregister();
-      //    double error = computeError(
-      //        reducedModelResults,
-      //        global.g_optimalReducedModelDisplacements[simulationScenarioIndex]);
-      //    std::cout << "Error of simulation scenario "
-      //              << simulationScenarioStrings[simulationScenarioIndex] << "
-      //              is "
-      //              << error << std::endl;
-    }
-  }
   void save(const string &saveToPath) const {
     assert(std::filesystem::is_directory(saveToPath));
 
@@ -325,20 +277,21 @@ struct ReducedModelOptimizer::Results {
     }
   }
 
-  void writeTo(const ReducedModelOptimizer::Settings &settings_optimization,
+  void
-               csvFile &os, const bool writeHeader = true) const {
+  writeHeaderTo(const ReducedModelOptimizer::Settings &settings_optimization,
-    if (writeHeader) {
+                csvFile &os) {
-      //// Write header to csv
+    os << "Obj value";
-      os << "Obj value";
+    for (const ReducedModelOptimizer::xRange &range :
-      for (const ReducedModelOptimizer::xRange &range :
+         settings_optimization.xRanges) {
-           settings_optimization.xRanges) {
+      os << range.label;
-        os << range.label;
-      }
-      os << "Time(s)";
-      os << "#Crashes";
-      //      os << std::endl;
-      os << endrow;
     }
+    os << "Time(s)";
+    os << "#Crashes";
+  }
+
+  void
+  writeResultsTo(const ReducedModelOptimizer::Settings &settings_optimization,
+                 csvFile &os) const {
     os << objectiveValue;
     for (const double &optimalX : x) {
       os << optimalX;
@@ -356,8 +309,6 @@ struct ReducedModelOptimizer::Results {
     } else {
       os << numberOfSimulationCrashes;
     }
-    //    os << std::endl;
-    os << endrow;
   }
 };