diff --git a/src/main.cpp b/src/main.cpp index d0fd56e..85e265f 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -153,60 +153,63 @@ int main(int argc, char *argv[]) csv_resultsLocalFile << endrow; } + ReducedModelEvaluator::evaluateReducedModel(optimizationResults); #ifdef POLYSCOPE_DEFINED - std::vector scenarioLabels( - optimizationResults.objectiveValue.perSimulationScenario_total.size()); - const double colorAxial = 1; - const double colorShear = 3; - const double colorBending = 5; - const double colorDome = 0.1; - const double colorSaddle = 0; - std::vector colors(optimizationResults.objectiveValue.perSimulationScenario_total.size()); - for (int scenarioIndex = 0; scenarioIndex < scenarioLabels.size(); scenarioIndex++) { - scenarioLabels[scenarioIndex] - = optimizationResults.reducedPatternSimulationJobs[scenarioIndex]->getLabel(); - if (scenarioLabels[scenarioIndex].rfind("Axial", 0) == 0) { - colors[scenarioIndex] = colorAxial; - } else if (scenarioLabels[scenarioIndex].rfind("Shear", 0) == 0) { - colors[scenarioIndex] = colorShear; + std::vector scenarioLabels( + optimizationResults.objectiveValue.perSimulationScenario_total.size()); + const double colorAxial = 1; + const double colorShear = 3; + const double colorBending = 5; + const double colorDome = 0.1; + const double colorSaddle = 0; + std::vector colors( + optimizationResults.objectiveValue.perSimulationScenario_total.size()); + for (int scenarioIndex = 0; scenarioIndex < scenarioLabels.size(); scenarioIndex++) { + scenarioLabels[scenarioIndex] + = optimizationResults.reducedPatternSimulationJobs[scenarioIndex]->getLabel(); + if (scenarioLabels[scenarioIndex].rfind("Axial", 0) == 0) { + colors[scenarioIndex] = colorAxial; - } else if (scenarioLabels[scenarioIndex].rfind("Bending", 0) == 0) { - colors[scenarioIndex] = colorBending; + } else if (scenarioLabels[scenarioIndex].rfind("Shear", 0) == 0) { + colors[scenarioIndex] = colorShear; - } else if (scenarioLabels[scenarioIndex].rfind("Dome", 0) == 0) { - colors[scenarioIndex] = colorDome; - } else if (scenarioLabels[scenarioIndex].rfind("Saddle", 0) == 0) { - colors[scenarioIndex] = colorSaddle; - } else { - std::cerr << "Label could not be identified" << std::endl; - } - } - std::vector y(optimizationResults.objectiveValue.perSimulationScenario_total.size()); - for (int scenarioIndex = 0; scenarioIndex < scenarioLabels.size(); scenarioIndex++) { - y[scenarioIndex] - // = optimizationResults.objectiveValue.perSimulationScenario_rawTranslational[scenarioIndex] - // + optimizationResults.objectiveValue.perSimulationScenario_rawRotational[scenarioIndex]; - = optimizationResults.objectiveValue.perSimulationScenario_total_unweighted[scenarioIndex]; - } - std::vector x = matplot::linspace(0, y.size() - 1, y.size()); - std::vector markerSizes(y.size(), 5); - SimulationResultsReporter::createPlot("scenario index", - "error", - x, - y, - markerSizes, - colors, - std::filesystem::path(resultsOutputDir) - .append("perScenarioObjectiveValues.png")); - // optimizationResults.saveMeshFiles(); - std::cout << "Saved results to:" << resultsOutputDir << std::endl; - // optimizationResults.draw(); + } else if (scenarioLabels[scenarioIndex].rfind("Bending", 0) == 0) { + colors[scenarioIndex] = colorBending; - // ReducedModelOptimization::Results optResults; - // optResults.load("/home/iason/Desktop/dlib_ensmallen_comparison/TestSets/" - // "singlePattern_dlib_firstSubmission/12@single_reduced(100000_1.20)"); - ReducedModelEvaluator::evaluateReducedModel(optimizationResults); + } else if (scenarioLabels[scenarioIndex].rfind("Dome", 0) == 0) { + colors[scenarioIndex] = colorDome; + } else if (scenarioLabels[scenarioIndex].rfind("Saddle", 0) == 0) { + colors[scenarioIndex] = colorSaddle; + } else { + std::cerr << "Label could not be identified" << std::endl; + } + } + std::vector y(optimizationResults.objectiveValue.perSimulationScenario_total.size()); + for (int scenarioIndex = 0; scenarioIndex < scenarioLabels.size(); scenarioIndex++) { + y[scenarioIndex] + // = optimizationResults.objectiveValue.perSimulationScenario_rawTranslational[scenarioIndex] + // + optimizationResults.objectiveValue.perSimulationScenario_rawRotational[scenarioIndex]; + = optimizationResults.objectiveValue + .perSimulationScenario_total_unweighted[scenarioIndex]; + } + std::vector x = matplot::linspace(0, y.size() - 1, y.size()); + std::vector markerSizes(y.size(), 5); + SimulationResultsReporter::createPlot("scenario index", + "error", + x, + y, + markerSizes, + colors, + std::filesystem::path(resultsOutputDir) + .append("perScenarioObjectiveValues.png")); + // optimizationResults.saveMeshFiles(); + std::cout << "Saved results to:" << resultsOutputDir << std::endl; + // optimizationResults.draw(); + + // ReducedModelOptimization::Results optResults; + // optResults.load("/home/iason/Desktop/dlib_ensmallen_comparison/TestSets/" + // "singlePattern_dlib_firstSubmission/12@single_reduced(100000_1.20)"); #endif return 0; } diff --git a/src/reducedmodeloptimizer.cpp b/src/reducedmodeloptimizer.cpp index 2f681dd..4791a5f 100644 --- a/src/reducedmodeloptimizer.cpp +++ b/src/reducedmodeloptimizer.cpp @@ -108,9 +108,9 @@ double ReducedModelOptimizer::computeRawRotationalError( { double rawRotationalError = 0; for (const auto &reducedToFullInterfaceViPair : reducedToFullInterfaceViMap) { - const double vertexRotationalError - = rotatedQuaternion_fullPattern[reducedToFullInterfaceViPair.second].angularDistance( - rotatedQuaternion_reducedPattern[reducedToFullInterfaceViPair.first]); + const double vertexRotationalError = std::abs( + rotatedQuaternion_fullPattern[reducedToFullInterfaceViPair.second].angularDistance( + rotatedQuaternion_reducedPattern[reducedToFullInterfaceViPair.first])); rawRotationalError += vertexRotationalError; } @@ -247,7 +247,7 @@ double ReducedModelOptimizer::objective(const std::vector &x) // FormFinder simulator; std::for_each( - std::execution::par_unseq, + // std::execution::par_unseq, global.simulationScenarioIndices.begin(), global.simulationScenarioIndices.end(), [&](const int &simulationScenarioIndex) { @@ -267,15 +267,11 @@ double ReducedModelOptimizer::objective(const std::vector &x) std::cout << "Failed simulation" << std::endl; #endif } else { -// double simulationScenarioError = 0; -// constexpr bool usePotentialEnergy = false; -// if (usePotentialEnergy) { -// simulationScenarioError += std::abs( -// reducedModelResults.internalPotentialEnergy -// - global.fullPatternResults[simulationScenarioIndex].internalPotentialEnergy -// / global.fullPatternResults[simulationScenarioIndex] -// .internalPotentialEnergy); -// } else { + const double simulationScenarioError_PE = std::abs( + (reducedModelResults.internalPotentialEnergy + - global.fullPatternResults[simulationScenarioIndex].internalPotentialEnergy) + / global.fullPatternResults[simulationScenarioIndex].internalPotentialEnergy); +// else { #ifdef POLYSCOPE_DEFINED #ifdef USE_SCENARIO_WEIGHTS const double scenarioWeight = global.scenarioWeights[simulationScenarioIndex]; @@ -285,7 +281,7 @@ double ReducedModelOptimizer::objective(const std::vector &x) #else const double scenarioWeight = 1; #endif - const double simulationScenarioError = computeError( + const double simulationScenarioError_displacements = computeError( global.fullPatternResults[simulationScenarioIndex], reducedModelResults, global.reducedToFullInterfaceViMap, @@ -294,7 +290,8 @@ double ReducedModelOptimizer::objective(const std::vector &x) scenarioWeight, global.objectiveWeights[simulationScenarioIndex]); - simulationErrorsPerScenario[simulationScenarioIndex] = simulationScenarioError; + simulationErrorsPerScenario[simulationScenarioIndex] + = simulationScenarioError_displacements + simulationScenarioError_PE; // } // #ifdef POLYSCOPE_DEFINED // reducedJob->pMesh->registerForDrawing(Colors::reducedInitial); @@ -978,7 +975,8 @@ ReducedModelOptimizer::getFullPatternMaxSimulationForces( .append(m_pFullPatternSimulationMesh->getLabel() + ".json")); const bool fullPatternScenarioMagnitudesExist = std::filesystem::exists( patternMaxForceMagnitudesFilePath); - if (fullPatternScenarioMagnitudesExist) { + constexpr bool recomputeMagnitudes = false; + if (fullPatternScenarioMagnitudesExist && !recomputeMagnitudes) { nlohmann::json json; std::ifstream ifs(patternMaxForceMagnitudesFilePath.string()); ifs >> json; @@ -1798,9 +1796,9 @@ void ReducedModelOptimizer::computeObjectiveValueNormalizationFactors() continue; } } - angularDistanceSum += global.fullPatternResults[simulationScenarioIndex] - .rotationalDisplacementQuaternion[fullPatternVi] - .angularDistance(Eigen::Quaterniond::Identity()); + angularDistanceSum += std::abs(global.fullPatternResults[simulationScenarioIndex] + .rotationalDisplacementQuaternion[fullPatternVi] + .angularDistance(Eigen::Quaterniond::Identity())); } fullPatternTranslationalDisplacementNormSum[simulationScenarioIndex] @@ -1916,7 +1914,7 @@ void ReducedModelOptimizer::optimize( .append(m_pFullPatternSimulationMesh->getLabel()) .append(pFullPatternSimulationJob->getLabel())); // .append(pFullPatternSimulationJob->getLabel() + ".json") - constexpr bool recomputeFullPatternResults = false; + constexpr bool recomputeFullPatternResults = true; SimulationResults fullPatternResults; if (!recomputeFullPatternResults && std::filesystem::exists(jobResultsDirectoryPath)) { fullPatternResults.load(std::filesystem::path(jobResultsDirectoryPath).append("Results"), diff --git a/src/reducedmodeloptimizer.hpp b/src/reducedmodeloptimizer.hpp index 7271256..d22d0bb 100644 --- a/src/reducedmodeloptimizer.hpp +++ b/src/reducedmodeloptimizer.hpp @@ -66,7 +66,7 @@ public: // inline constexpr static ParameterLabels parameterLabels(); inline static std::array - parameterLabels = {"R", "A", "I2", "I3", "J", "Theta", "R"}; + parameterLabels = {"E", "A", "I2", "I3", "J", "Theta", "R"}; constexpr static std::array simulationScenariosResolution = {11, 11, 20, 20, 20}; constexpr static std::array