diff --git a/CMakeLists.txt b/CMakeLists.txt index 5b7c38d..d1452db 100755 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,7 +1,10 @@ cmake_minimum_required(VERSION 2.8) project(MySources) -set(CMAKE_CXX_STANDARD 20) -set(CMAKE_CXX_STANDARD_REQUIRED ON) + +file(GLOB MySourcesFiles ${CMAKE_CURRENT_LIST_DIR}/*.hpp ${CMAKE_CURRENT_LIST_DIR}/*.cpp) +add_library(${PROJECT_NAME} ${MySourcesFiles} ) +target_compile_features(${PROJECT_NAME} PUBLIC cxx_std_20) + #Download external dependencies NOTE: If the user has one of these libs it shouldn't be downloaded again. include(${CMAKE_CURRENT_LIST_DIR}/cmake/DownloadProject.cmake) if (CMAKE_VERSION VERSION_LESS 3.2) @@ -26,16 +29,18 @@ if(${USE_POLYSCOPE}) ) add_subdirectory(${POLYSCOPE_SOURCE_DIR} ${POLYSCOPE_BINARY_DIR}) add_compile_definitions(POLYSCOPE_DEFINED) + target_sources(${PROJECT_NAME} PUBLIC ${POLYSCOPE_SOURCE_DIR}/deps/imgui/imgui/misc/cpp/imgui_stdlib.cpp) endif() ##vcglib devel branch download_project(PROJ vcglib_devel - GIT_REPOSITORY https://github.com/IasonManolas/vcglib.git + GIT_REPOSITORY https://gitea-s2i2s.isti.cnr.it/manolas/vcglib.git GIT_TAG devel PREFIX ${EXTERNAL_DEPS_DIR} ${UPDATE_DISCONNECTED_IF_AVAILABLE} ) add_subdirectory(${vcglib_devel_SOURCE_DIR} ${vcglib_devel_BINARY_DIR}) +target_sources(${PROJECT_NAME} PUBLIC ${vcglib_devel_SOURCE_DIR}/wrap/ply/plylib.cpp) ##matplot++ lib set(MATPLOTPLUSPLUS_BINARY_DIR ${CMAKE_CURRENT_BINARY_DIR}/matplot) @@ -77,23 +82,26 @@ if(MSVC) add_compile_definitions(_HAS_STD_BYTE=0) endif(MSVC) #link_directories(${CMAKE_CURRENT_LIST_DIR}/boost_graph/libs) -file(GLOB MySourcesFiles ${CMAKE_CURRENT_LIST_DIR}/*.hpp ${CMAKE_CURRENT_LIST_DIR}/*.cpp) -add_library(${PROJECT_NAME} ${MySourcesFiles} ${vcglib_devel_SOURCE_DIR}/wrap/ply/plylib.cpp) -target_include_directories(${PROJECT_NAME} - PUBLIC ${CMAKE_CURRENT_LIST_DIR}/boost_graph - PUBLIC ${vcglib_devel_SOURCE_DIR} - PUBLIC ${threed-beam-fea_SOURCE_DIR} - ) + if(${MYSOURCES_STATIC_LINK}) message("Linking statically") - target_link_libraries(${PROJECT_NAME} -static Eigen3::Eigen matplot ThreedBeamFEA ${TBB_BINARY_DIR}/libtbb_static.a pthread gfortran quadmath) + target_link_libraries(${PROJECT_NAME} PUBLIC -static Eigen3::Eigen matplot ThreedBeamFEA ${TBB_BINARY_DIR}/libtbb_static.a #[[tbb_static]] pthread gfortran quadmath) else() - target_link_libraries(${PROJECT_NAME} Eigen3::Eigen matplot ThreedBeamFEA tbb pthread) +# set_target_properties(${PROJECT_NAME} PROPERTIES POSITION_INDEPENDENT_CODE TRUE) + target_link_libraries(${PROJECT_NAME} PUBLIC Eigen3::Eigen matplot ThreedBeamFEA tbb pthread) if(${USE_POLYSCOPE}) - target_link_libraries(${PROJECT_NAME} polyscope) + message("Using polyscope") + target_link_libraries(${PROJECT_NAME} PUBLIC polyscope) endif() endif() target_link_directories(MySources PUBLIC ${CMAKE_CURRENT_LIST_DIR}/boost_graph/libs) +target_include_directories(${PROJECT_NAME} + PUBLIC ${CMAKE_CURRENT_LIST_DIR}/boost_graph + PUBLIC ${vcglib_devel_SOURCE_DIR} +# PUBLIC ${threed-beam-fea_SOURCE_DIR} + PUBLIC ThreedBeamFEA + PUBLIC matplot + ) if(USE_ENSMALLEN) ##ENSMALLEN @@ -113,5 +121,5 @@ if(USE_ENSMALLEN) else() target_link_libraries(${PROJECT_NAME} PUBLIC armadillo ensmallen) endif() - target_include_directories(${PROJECT_NAME} PUBLIC ${ARMADILLO_SOURCE_DIR}/include) + target_include_directories(${PROJECT_NAME} PUBLIC ${ARMADILLO_SOURCE_DIR}/include ${ENSMALLEN_SOURCE_DIR}) endif() diff --git a/drmsimulationmodel.cpp b/drmsimulationmodel.cpp index 4326d6e..3b68b57 100755 --- a/drmsimulationmodel.cpp +++ b/drmsimulationmodel.cpp @@ -1889,7 +1889,7 @@ void DRMSimulationModel::printCurrentState() const void DRMSimulationModel::printDebugInfo() const { std::cout << pMesh->elements[0].rigidity.toString() << std::endl; - std::cout << "Number of dampings:" << numOfDampings << endl; + std::cout << "Number of dampings:" << numOfDampings << std::endl; printCurrentState(); } diff --git a/drmsimulationmodel.hpp b/drmsimulationmodel.hpp index c3e7c39..4d7e052 100755 --- a/drmsimulationmodel.hpp +++ b/drmsimulationmodel.hpp @@ -152,7 +152,7 @@ private: const std::unordered_map> &fixedVertices); #ifdef POLYSCOPE_DEFINED - void draw(const string &screenshotsFolder= {}); + void draw(const std::string &screenshotsFolder = {}); #endif void updateNodalInternalForce(Vector6d &nodalInternalForce, diff --git a/edgemesh.cpp b/edgemesh.cpp index aa040fa..3fa7665 100755 --- a/edgemesh.cpp +++ b/edgemesh.cpp @@ -25,7 +25,7 @@ Eigen::MatrixX3d VCGEdgeMesh::getEigenEdgeNormals() const { return eigenEdgeNormals; } -bool VCGEdgeMesh::save(const string &plyFilename) +bool VCGEdgeMesh::save(const std::string &plyFilename) { std::string filename = plyFilename; if (filename.empty()) { @@ -197,36 +197,36 @@ bool VCGEdgeMesh::createSpanGrid(const size_t desiredWidth, return true; } -bool VCGEdgeMesh::load(const string &plyFilename) { +bool VCGEdgeMesh::load(const std::string &plyFilename) +{ + std::string usedPath = plyFilename; + if (std::filesystem::path(plyFilename).is_relative()) { + usedPath = std::filesystem::absolute(plyFilename).string(); + } + assert(std::filesystem::exists(usedPath)); + Clear(); + // if (!loadUsingNanoply(usedPath)) { + // std::cerr << "Error: Unable to open " + usedPath << std::endl; + // return false; + // } + if (!loadUsingDefaultLoader(usedPath)) { + std::cerr << "Error: Unable to open " + usedPath << std::endl; + return false; + } - std::string usedPath = plyFilename; - if (std::filesystem::path(plyFilename).is_relative()) { - usedPath = std::filesystem::absolute(plyFilename).string(); - } - assert(std::filesystem::exists(usedPath)); - Clear(); - // if (!loadUsingNanoply(usedPath)) { - // std::cerr << "Error: Unable to open " + usedPath << std::endl; - // return false; - // } - if (!loadUsingDefaultLoader(usedPath)) { - std::cerr << "Error: Unable to open " + usedPath << std::endl; - return false; - } + getEdges(eigenEdges); + getVertices(eigenVertices); + vcg::tri::UpdateTopology::VertexEdge(*this); + label = std::filesystem::path(plyFilename).stem().string(); - getEdges(eigenEdges); - getVertices(eigenVertices); - vcg::tri::UpdateTopology::VertexEdge(*this); - label = std::filesystem::path(plyFilename).stem().string(); + const bool printDebugInfo = false; + if (printDebugInfo) { + std::cout << plyFilename << " was loaded successfuly." << std::endl; + std::cout << "Mesh has " << EN() << " edges." << std::endl; + } - const bool printDebugInfo = false; - if (printDebugInfo) { - std::cout << plyFilename << " was loaded successfuly." << std::endl; - std::cout << "Mesh has " << EN() << " edges." << std::endl; - } - - label=std::filesystem::path(plyFilename).stem().string(); - return true; + label = std::filesystem::path(plyFilename).stem().string(); + return true; } //bool VCGEdgeMesh::loadUsingNanoply(const std::string &plyFilename) { diff --git a/linearsimulationmodel.cpp b/linearsimulationmodel.cpp index 48a6f96..e1f19c6 100644 --- a/linearsimulationmodel.cpp +++ b/linearsimulationmodel.cpp @@ -108,8 +108,6 @@ SimulationResults LinearSimulationModel::getResults( Eigen::Quaternion q_nz; q_nz = Eigen::AngleAxis(nodalDisplacement[5], Eigen::Vector3d(0, 0, 1)); results.rotationalDisplacementQuaternion[vi] = q_nx * q_ny * q_nz; - // results.rotationalDisplacementQuaternion[vi] = q_nz * q_ny * q_nx; - // results.rotationalDisplacementQuaternion[vi] = q_nz * q_nx * q_ny; } results.setLabelPrefix("Linear"); @@ -214,7 +212,11 @@ SimulationResults LinearSimulationModel::executeSimulation( const std::shared_ptr &pSimulationJob) { assert(pSimulationJob->pMesh->VN() != 0); - if (!simulator.wasInitialized()) { + const bool wasInitializedWithDifferentStructure = pSimulationJob->pMesh->EN() + != simulator.structure.elems.size() + || pSimulationJob->pMesh->VN() + != simulator.structure.nodes.size(); + if (!simulator.wasInitialized() || wasInitializedWithDifferentStructure) { setStructure(pSimulationJob->pMesh); } // printInfo(job); diff --git a/reducedmodeloptimizer_structs.hpp b/reducedmodeloptimizer_structs.hpp index 4cc695a..e79f126 100644 --- a/reducedmodeloptimizer_structs.hpp +++ b/reducedmodeloptimizer_structs.hpp @@ -34,9 +34,9 @@ struct Colors struct xRange { - std::string label; - double min; - double max; + std::string label{}; + double min{0}; + double max{0}; inline bool operator<(const xRange &other) { return label < other.label; } std::string toString() const @@ -60,120 +60,321 @@ struct xRange { return label == xrange.label && min == xrange.min && max == xrange.max; } + + std::tuple toTuple() const + { + return std::make_tuple(label, min, max); + } + + void set(const std::tuple &inputTuple) + { + if (std::get<1>(inputTuple) > std::get<2>(inputTuple)) { + std::cerr + << "Invalid xRange tuple. Second argument(min) cant be smaller than the third(max)" + << std::endl; + std::terminate(); + // return; + } + std::tie(label, min, max) = inputTuple; + } }; -enum OptimizationParameterIndex { E, A, I2, I3, J, R, Theta, NumberOfOptimizationParameters }; +enum OptimizationParameterIndex { E, A, I2, I3, J, R, Theta, NumberOfOptimizationVariables }; + +inline int getParameterIndex(const std::string &s) +{ + if ("E" == s) { + return E; + } else if ("A" == s) { + return A; + } else if ("I2" == s) { + return I2; + } else if ("I3" == s) { + return I3; + } else if ("J" == s) { + return J; + } else if ("R" == s || "HexSize" == s) { + return R; + } else if ("Theta" == s || "HexAngle" == s) { + return Theta; + } else { + std::cerr << "Input is not recognized as a valid optimization variable index:" << s + << std::endl; + return -1; + } +} struct Settings { - std::filesystem::path intermediateResultsDirectoryPath; - std::vector> optimizationVariables; + inline static std::string defaultFilename{"OptimizationSettings.json"}; + std::vector> optimizationStrategy; + std::optional> + optimizationVariablesGroupsWeights; //TODO:should be removed in the future if not a splitting strategy is used for the optimization enum NormalizationStrategy { NonNormalized, Epsilon }; - std::vector parameterRanges; - inline static vector normalizationStrategyStrings{"NonNormalized", "Epsilon"}; - int numberOfFunctionCalls{100000}; - double solverAccuracy{1e-2}; + inline static std::vector normalizationStrategyStrings{"NonNormalized", "Epsilon"}; NormalizationStrategy normalizationStrategy{Epsilon}; - double translationNormalizationParameter{0.003}; - double rotationNormalizationParameter{3}; + std::array variablesRanges; + int numberOfFunctionCalls{100000}; + double solverAccuracy{1e-3}; + double translationNormalizationEpsilon{1e-4}; + double rotationNormalizationEpsilon{vcg::math::ToRad(3.0)}; struct ObjectiveWeights { double translational{1}; double rotational{1}; - } objectiveWeights; + bool operator==(const ObjectiveWeights &other) const + { + return this->translational == other.translational + && this->rotational == other.rotational; + } + }; + std::array perBaseScenarioObjectiveWeights; + std::array, NumberOfBaseSimulationScenarios> + convertObjectiveWeightsToPairs() const + { + std::array, NumberOfBaseSimulationScenarios> objectiveWeightsPairs; + for (int baseScenario = Axial; baseScenario != NumberOfBaseSimulationScenarios; + baseScenario++) { + objectiveWeightsPairs[baseScenario] + = std::make_pair(perBaseScenarioObjectiveWeights[baseScenario].translational, + perBaseScenarioObjectiveWeights[baseScenario].rotational); + } + return objectiveWeightsPairs; + } struct JsonKeys { - inline static std::string filename{"OptimizationSettings.json"}; + inline static std::string OptimizationStrategy{"OptimizationStrategy"}; + inline static std::string OptimizationStrategyGroupWeights{ + "OptimizationStrategyGroupWeights"}; inline static std::string OptimizationVariables{"OptimizationVariables"}; inline static std::string NumberOfFunctionCalls{"NumberOfFunctionCalls"}; inline static std::string SolverAccuracy{"SolverAccuracy"}; - inline static std::string TranslationalObjectiveWeight{"TransObjWeight"}; - inline static std::string OptimizationParameters{"OptimizationParameters"}; + inline static std::string ObjectiveWeights{"ObjectiveWeight"}; }; + void setDefault() + { + ReducedPatternOptimization::xRange beamE{"E", 0.001, 1000}; + ReducedPatternOptimization::xRange beamA{"A", 0.001, 1000}; + ReducedPatternOptimization::xRange beamI2{"I2", 0.001, 1000}; + ReducedPatternOptimization::xRange beamI3{"I3", 0.001, 1000}; + ReducedPatternOptimization::xRange beamJ{"J", 0.001, 1000}; + ReducedPatternOptimization::xRange innerHexagonSize{"R", 0.05, 0.95}; + ReducedPatternOptimization::xRange innerHexagonAngle{"Theta", -30.0, 30.0}; + variablesRanges = {beamE, beamA, beamI2, beamI3, beamJ, innerHexagonSize, innerHexagonAngle}; + numberOfFunctionCalls = 100000; + + enum OptimizationParameterComparisonScenarioIndex { + AllVar, + GeoYM, + MatGeo, + YMMat_Geo, + YM_MatGeo, + MatGeo_YM, + Geo_YM_Mat, + YM_Geo_Mat, + Geo_Mat, + YMGeo_Mat, + NumberOfScenarios + }; + const std::vector< + std::vector>> + optimizationParameters = [&]() { + std::vector< + std::vector>> + optimizationParameters(NumberOfScenarios); + using namespace ReducedPatternOptimization; + optimizationParameters[AllVar] = {{E, A, I2, I3, J, R, Theta}}; + optimizationParameters[GeoYM] = {{R, Theta, E}}; + optimizationParameters[MatGeo] = {{A, I2, I3, J, R, Theta}}; + optimizationParameters[YMMat_Geo] = {{E, A, I2, I3, J}, {R, Theta}}; + optimizationParameters[YM_MatGeo] = {{E}, {A, I2, I3, J, R, Theta}}; + optimizationParameters[MatGeo_YM] = {{A, I2, I3, J, R, Theta}, {E}}; + optimizationParameters[Geo_YM_Mat] = {{R, Theta}, {E}, {A, I2, I3, J}}; + optimizationParameters[YM_Geo_Mat] = {{E}, {R, Theta}, {A, I2, I3, J}}; + optimizationParameters[Geo_Mat] = {{R, Theta}, {A, I2, I3, J}}; + optimizationParameters[YMGeo_Mat] = {{E, R, Theta}, {A, I2, I3, J}}; + return optimizationParameters; + }(); + + constexpr OptimizationParameterComparisonScenarioIndex scenario = MatGeo; + optimizationStrategy = optimizationParameters[scenario]; + if (scenario == YMGeo_Mat) { + optimizationVariablesGroupsWeights = {0.15, 0.85}; + } + + normalizationStrategy = ReducedPatternOptimization::Settings::NormalizationStrategy::Epsilon; + //#ifdef POLYSCOPE_DEFINED + //#else + translationNormalizationEpsilon = 4e-3; + rotationNormalizationEpsilon = vcg::math::ToRad(3.0); + + // translationNormalizationParameter = 0; + // rotationNormalizationParameter = 0; + // translationNormalizationParameter = 1e-3; + // rotationNormalizationParameter = vcg::math::ToRad(3.0); + // translationNormalizationParameter = 5e-2; + // rotationNormalizationParameter = vcg::math::ToRad(3.0); + //#endif + solverAccuracy = 1e-2; + //TODO refactor that + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Axial].translational = 1.5; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Axial].rotational + = 2 - perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Axial].translational; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Shear].translational = 1.5; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Shear].rotational + = 2 - perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Shear].translational; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Bending].translational = 1.2; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Bending].rotational + = 2 + - perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Bending].translational; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Dome].translational = 1.95; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Dome].rotational + = 2 - perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Dome].translational; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Saddle].translational = 1.2; + perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Saddle].rotational + = 2 - perBaseScenarioObjectiveWeights[ReducedPatternOptimization::Saddle].translational; + } + void save(const std::filesystem::path &saveToPath) { assert(std::filesystem::is_directory(saveToPath)); nlohmann::json json; - //write x ranges - for (size_t xRangeIndex = 0; xRangeIndex < parameterRanges.size(); xRangeIndex++) { - const auto &xRange = parameterRanges[xRangeIndex]; - json[JsonKeys::OptimizationVariables + "_" + std::to_string(xRangeIndex)] - = xRange.toString(); + json[JsonKeys::OptimizationStrategy] = optimizationStrategy; + if (optimizationVariablesGroupsWeights.has_value()) { + json[JsonKeys::OptimizationStrategyGroupWeights] = optimizationVariablesGroupsWeights + .value(); } + //write x ranges + const std::array, NumberOfOptimizationVariables> + xRangesAsTuples = [=]() { + std::array, NumberOfOptimizationVariables> + xRangesAsTuples; + for (int optimizationParameterIndex = E; + optimizationParameterIndex != NumberOfOptimizationVariables; + optimizationParameterIndex++) { + xRangesAsTuples[optimizationParameterIndex] + = variablesRanges[optimizationParameterIndex].toTuple(); + } + return xRangesAsTuples; + }(); + json[JsonKeys::OptimizationVariables] = xRangesAsTuples; + // for (size_t xRangeIndex = 0; xRangeIndex < variablesRanges.size(); xRangeIndex++) { + // const auto &xRange = variablesRanges[xRangeIndex]; + // json[JsonKeys::OptimizationVariables + "_" + std::to_string(xRangeIndex)] + // = xRange.toString(); + // } json[JsonKeys::NumberOfFunctionCalls] = numberOfFunctionCalls; json[JsonKeys::SolverAccuracy] = solverAccuracy; - json[JsonKeys::TranslationalObjectiveWeight] = objectiveWeights.translational; - json[JsonKeys::OptimizationParameters] = optimizationVariables; + //Objective weights + std::array, NumberOfBaseSimulationScenarios> objectiveWeightsPairs; + std::transform(perBaseScenarioObjectiveWeights.begin(), + perBaseScenarioObjectiveWeights.end(), + objectiveWeightsPairs.begin(), + [](const ObjectiveWeights &objectiveWeights) { + return std::make_pair(objectiveWeights.translational, + objectiveWeights.rotational); + }); + json[JsonKeys::ObjectiveWeights] = objectiveWeightsPairs; std::filesystem::path jsonFilePath( - std::filesystem::path(saveToPath).append(JsonKeys::filename)); + std::filesystem::path(saveToPath).append(defaultFilename)); std::ofstream jsonFile(jsonFilePath.string()); jsonFile << json; } - bool load(const std::filesystem::path &loadFromPath) + bool load(const std::filesystem::path &jsonFilepath) { - assert(std::filesystem::is_directory(loadFromPath)); - //Load optimal X nlohmann::json json; - std::filesystem::path jsonFilepath( - std::filesystem::path(loadFromPath).append(JsonKeys::filename)); if (!std::filesystem::exists(jsonFilepath)) { + std::cerr << "Optimization settings could not be loaded because input filepath does " + "not exist:" + << jsonFilepath << std::endl; return false; } std::ifstream ifs(jsonFilepath); ifs >> json; + if (json.contains(JsonKeys::OptimizationStrategy)) { + optimizationStrategy + = std::vector>( + (json.at(JsonKeys::OptimizationStrategy))); + } + if (json.contains(JsonKeys::OptimizationStrategyGroupWeights)) { + optimizationVariablesGroupsWeights = std::vector( + json[JsonKeys::OptimizationStrategyGroupWeights]); + } //read x ranges - size_t xRangeIndex = 0; - while (true) { - const std::string jsonXRangeKey = JsonKeys::OptimizationVariables + "_" - + std::to_string(xRangeIndex); - if (!json.contains(jsonXRangeKey)) { - break; + if (json.contains(JsonKeys::OptimizationVariables)) { + const std::array, NumberOfOptimizationVariables> + xRangesAsTuples = json.at(JsonKeys::OptimizationVariables); + for (const auto &rangeTuple : xRangesAsTuples) { + variablesRanges[getParameterIndex(std::get<0>(rangeTuple))].set(rangeTuple); + } + } else { //NOTE:legacy compatibility + size_t xRangeIndex = 0; + while (true) { + const std::string jsonXRangeKey = JsonKeys::OptimizationVariables + "_" + + std::to_string(xRangeIndex++); + if (!json.contains(jsonXRangeKey)) { + break; + } + xRange x; + x.fromString(json.at(jsonXRangeKey)); + variablesRanges[getParameterIndex(x.label)] = x; } - xRange x; - x.fromString(json.at(jsonXRangeKey)); - parameterRanges.push_back(x); - xRangeIndex++; } - optimizationVariables - = std::vector>( - (json.at(JsonKeys::OptimizationParameters))); numberOfFunctionCalls = json.at(JsonKeys::NumberOfFunctionCalls); solverAccuracy = json.at(JsonKeys::SolverAccuracy); - objectiveWeights.translational = json.at(JsonKeys::TranslationalObjectiveWeight); - objectiveWeights.rotational = 2 - objectiveWeights.translational; + //Objective weights + if (json.contains(JsonKeys::ObjectiveWeights)) { + std::array, NumberOfBaseSimulationScenarios> + objectiveWeightsPairs = json.at(JsonKeys::ObjectiveWeights); + std::transform(objectiveWeightsPairs.begin(), + objectiveWeightsPairs.end(), + perBaseScenarioObjectiveWeights.begin(), + [](const std::pair &objectiveWeightsPair) { + return ObjectiveWeights{objectiveWeightsPair.first, + objectiveWeightsPair.second}; + }); + } + + // perBaseScenarioObjectiveWeights = json.at(JsonKeys::ObjectiveWeights); + // objectiveWeights.translational = json.at(JsonKeys::ObjectiveWeights); + // objectiveWeights.rotational = 2 - objectiveWeights.translational; return true; } std::string toString() const { std::string settingsString; - if (!parameterRanges.empty()) { + if (!variablesRanges.empty()) { std::string xRangesString; - for (const xRange &range : parameterRanges) { + for (const xRange &range : variablesRanges) { xRangesString += range.toString() + " "; } settingsString += xRangesString; } settingsString += "FuncCalls=" + std::to_string(numberOfFunctionCalls) - + " Accuracy=" + std::to_string(solverAccuracy) - + " TransWeight=" + std::to_string(objectiveWeights.translational) - + " RotWeight=" + std::to_string(objectiveWeights.rotational); + + " Accuracy=" + std::to_string(solverAccuracy); + for (int baseScenario = Axial; baseScenario != NumberOfBaseSimulationScenarios; + baseScenario++) { + +" W_t=" + std::to_string(perBaseScenarioObjectiveWeights[baseScenario].translational) + + " W_r=" + + std::to_string(perBaseScenarioObjectiveWeights[baseScenario].rotational); + } return settingsString; } void writeHeaderTo(csvFile &os) const { - if (!parameterRanges.empty()) { - for (const xRange &range : parameterRanges) { + if (!variablesRanges.empty()) { + for (const xRange &range : variablesRanges) { os << range.label + " max"; os << range.label + " min"; } @@ -181,16 +382,15 @@ struct Settings os << "Function Calls"; os << "Solution Accuracy"; os << "Normalization strategy"; - os << "Trans weight"; - os << "Rot weight"; + os << JsonKeys::ObjectiveWeights; os << "Optimization parameters"; // os << std::endl; } void writeSettingsTo(csvFile &os) const { - if (!parameterRanges.empty()) { - for (const xRange &range : parameterRanges) { + if (!variablesRanges.empty()) { + for (const xRange &range : variablesRanges) { os << range.max; os << range.min; } @@ -198,13 +398,22 @@ struct Settings os << numberOfFunctionCalls; os << solverAccuracy; os << normalizationStrategyStrings[normalizationStrategy] + "_" - + std::to_string(translationNormalizationParameter); - os << objectiveWeights.translational; - os << objectiveWeights.rotational; + + std::to_string(translationNormalizationEpsilon); + std::string objectiveWeightsString; + objectiveWeightsString += "{"; + for (int baseScenario = Axial; baseScenario != NumberOfBaseSimulationScenarios; + baseScenario++) { + objectiveWeightsString + += "{" + std::to_string(perBaseScenarioObjectiveWeights[baseScenario].translational) + + "," + std::to_string(perBaseScenarioObjectiveWeights[baseScenario].rotational) + + "}"; + } + objectiveWeightsString += "}"; + os << objectiveWeightsString; //export optimization parameters std::vector> vv; - for (const std::vector &v : optimizationVariables) { + for (const std::vector &v : optimizationStrategy) { std::vector vi; vi.reserve(v.size()); for (const OptimizationParameterIndex ¶meter : v) { @@ -218,12 +427,17 @@ struct Settings inline bool operator==(const Settings &settings1, const Settings &settings2) noexcept { + const bool haveTheSameObjectiveWeights + = std::mismatch(settings1.perBaseScenarioObjectiveWeights.begin(), + settings1.perBaseScenarioObjectiveWeights.end(), + settings2.perBaseScenarioObjectiveWeights.begin()) + .first + != settings1.perBaseScenarioObjectiveWeights.end(); return settings1.numberOfFunctionCalls == settings2.numberOfFunctionCalls - && settings1.parameterRanges == settings2.parameterRanges - && settings1.solverAccuracy == settings2.solverAccuracy - && settings1.objectiveWeights.translational == settings2.objectiveWeights.translational - && settings1.translationNormalizationParameter - == settings2.translationNormalizationParameter; + && settings1.variablesRanges == settings2.variablesRanges + && settings1.solverAccuracy == settings2.solverAccuracy && haveTheSameObjectiveWeights + && settings1.translationNormalizationEpsilon + == settings2.translationNormalizationEpsilon; } inline void updateMeshWithOptimalVariables(const std::vector &x, SimulationMesh &m) @@ -306,6 +520,7 @@ struct Settings std::vector perSimulationScenario_translational; std::vector perSimulationScenario_rotational; std::vector perSimulationScenario_total; + std::vector perSimulationScenario_total_unweighted; } objectiveValue; std::vector perScenario_fullPatternPotentialEnergy; std::vector objectiveValueHistory; @@ -332,7 +547,7 @@ struct Settings inline static std::string fullPatternYoungsModulus{"youngsModulus"}; }; - void save(const string &saveToPath, const bool shouldExportDebugFiles = false) + void save(const std::string &saveToPath, const bool shouldExportDebugFiles = false) { //clear directory if (std::filesystem::exists(saveToPath)) { @@ -405,7 +620,8 @@ struct Settings = fullPatternSimulationJobs[simulationScenarioIndex]; std::filesystem::path simulationJobFolderPath( std::filesystem::path(simulationJobsPath) - .append(pFullPatternSimulationJob->label)); + .append(std::to_string(simulationScenarioIndex) + "_" + + pFullPatternSimulationJob->label)); std::filesystem::create_directories(simulationJobFolderPath); const auto fullPatternDirectoryPath = std::filesystem::path(simulationJobFolderPath).append("Full"); @@ -435,7 +651,7 @@ struct Settings #endif } - bool load(const string &loadFromPath, const bool &shouldLoadDebugFiles = false) + bool load(const std::string &loadFromPath, const bool &shouldLoadDebugFiles = false) { assert(std::filesystem::is_directory(loadFromPath)); std::filesystem::path jsonFilepath( @@ -496,40 +712,36 @@ struct Settings if (shouldLoadDebugFiles) { const std::filesystem::path simulationJobsFolderPath( std::filesystem::path(loadFromPath).append("SimulationJobs")); - for (const auto &directoryEntry : - filesystem::directory_iterator(simulationJobsFolderPath)) { - const auto simulationScenarioPath = directoryEntry.path(); + + std::vector sortedByName; + for (auto &entry : std::filesystem::directory_iterator(simulationJobsFolderPath)) + sortedByName.push_back(entry.path()); + + std::sort(sortedByName.begin(), sortedByName.end(), &Utilities::compareNat); + + for (const auto &simulationScenarioPath : sortedByName) { if (!std::filesystem::is_directory(simulationScenarioPath)) { continue; } - // Load full pattern files - for (const auto &fileEntry : filesystem::directory_iterator( - std::filesystem::path(simulationScenarioPath).append("Full"))) { - const auto filepath = fileEntry.path(); - if (filepath.extension() == ".json") { - SimulationJob job; - job.load(filepath.string()); - job.pMesh->setBeamMaterial(0.3, fullPatternYoungsModulus); - fullPatternSimulationJobs.push_back(std::make_shared(job)); - } - } - // Load reduced pattern files and apply optimal parameters - for (const auto &fileEntry : filesystem::directory_iterator( - std::filesystem::path(simulationScenarioPath).append("Reduced"))) { - const auto filepath = fileEntry.path(); - if (filepath.extension() == ".json") { - SimulationJob job; - job.load(filepath.string()); - // applyOptimizationResults_innerHexagon(*this, baseTriangle, *job.pMesh); - applyOptimizationResults_elements(*this, job.pMesh); - reducedPatternSimulationJobs.push_back( - std::make_shared(job)); - } - } + const auto fullJobFilepath = Utilities::getFilepathWithExtension( + std::filesystem::path(simulationScenarioPath).append("Full"), ".json"); + SimulationJob fullJob; + fullJob.load(fullJobFilepath.string()); + fullJob.pMesh->setBeamMaterial(0.3, fullPatternYoungsModulus); + fullPatternSimulationJobs.push_back(std::make_shared(fullJob)); + + const auto reducedJobFilepath = Utilities::getFilepathWithExtension( + std::filesystem::path(simulationScenarioPath).append("Reduced"), ".json"); + SimulationJob reducedJob; + reducedJob.load(reducedJobFilepath.string()); + // applyOptimizationResults_innerHexagon(*this, baseTriangle, *job.pMesh); + applyOptimizationResults_elements(*this, reducedJob.pMesh); + reducedPatternSimulationJobs.push_back( + std::make_shared(reducedJob)); } } - settings.load(loadFromPath); + settings.load(std::filesystem::path(loadFromPath).append(Settings::defaultFilename)); return true; } @@ -543,9 +755,16 @@ struct Settings std::unordered_map optimalXVariables(reducedPattern_optimizationResults.optimalXNameValuePairs.begin(), reducedPattern_optimizationResults.optimalXNameValuePairs.end()); - assert(optimalXVariables.contains("HexSize") && optimalXVariables.contains("HexAngle")); - applyOptimizationResults_innerHexagon(optimalXVariables.at("HexSize"), - optimalXVariables.at("HexAngle"), + assert(optimalXVariables.contains("R") && optimalXVariables.contains("Theta")); + if (optimalXVariables.contains("HexSize")) { + applyOptimizationResults_innerHexagon(optimalXVariables.at("HexSize"), + optimalXVariables.at("HexAngle"), + patternBaseTriangle, + reducedPattern); + return; + } + applyOptimizationResults_innerHexagon(optimalXVariables.at("R"), + optimalXVariables.at("Theta"), patternBaseTriangle, reducedPattern); } @@ -589,7 +808,7 @@ struct Settings static void applyOptimizationResults_elements( const ReducedPatternOptimization::Results &reducedPattern_optimizationResults, - const shared_ptr &pTiledReducedPattern_simulationMesh) + const std::shared_ptr &pTiledReducedPattern_simulationMesh) { assert(CrossSectionType::name == RectangularBeamDimensions::name); // Set reduced parameters fron the optimization results diff --git a/simulation_structs.hpp b/simulation_structs.hpp index fcc2a4e..e30240b 100755 --- a/simulation_structs.hpp +++ b/simulation_structs.hpp @@ -133,32 +133,32 @@ public: return constrainedVertices.empty() && nodalExternalForces.empty() && nodalForcedDisplacements.empty() && pMesh == nullptr; } - void remap(const std::unordered_map &thisToOtherViMap, - SimulationJob &simulationJobMapped) const + void remap(const std::unordered_map &sourceToDestinationViMap, + SimulationJob &destination_simulationJob) const { - assert(simulationJobMapped.pMesh->VN() != 0); - std::unordered_map> reducedModelFixedVertices; - for (auto fullModelFixedVertex : this->constrainedVertices) { - reducedModelFixedVertices[thisToOtherViMap.at(fullModelFixedVertex.first)] - = fullModelFixedVertex.second; + std::unordered_map> destination_fixedVertices; + for (const auto &source_fixedVertex : this->constrainedVertices) { + destination_fixedVertices[sourceToDestinationViMap.at(source_fixedVertex.first)] + = source_fixedVertex.second; } - std::unordered_map reducedModelNodalForces; - for (auto fullModelNodalForce : this->nodalExternalForces) { - reducedModelNodalForces[thisToOtherViMap.at(fullModelNodalForce.first)] - = fullModelNodalForce.second; + std::unordered_map destination_nodalForces; + for (const auto &source_nodalForces : this->nodalExternalForces) { + destination_nodalForces[sourceToDestinationViMap.at(source_nodalForces.first)] + = source_nodalForces.second; } - std::unordered_map reducedNodalForcedDisplacements; - for (auto fullForcedDisplacement : this->nodalForcedDisplacements) { - reducedNodalForcedDisplacements[thisToOtherViMap.at(fullForcedDisplacement.first)] - = fullForcedDisplacement.second; + std::unordered_map destination_forcedDisplacements; + for (const auto &source_forcedDisplacements : this->nodalForcedDisplacements) { + destination_forcedDisplacements[sourceToDestinationViMap.at( + source_forcedDisplacements.first)] + = source_forcedDisplacements.second; } - simulationJobMapped.constrainedVertices = reducedModelFixedVertices; - simulationJobMapped.nodalExternalForces = reducedModelNodalForces; - simulationJobMapped.label = this->getLabel(); - simulationJobMapped.nodalForcedDisplacements = reducedNodalForcedDisplacements; + destination_simulationJob.constrainedVertices = destination_fixedVertices; + destination_simulationJob.nodalExternalForces = destination_nodalForces; + destination_simulationJob.label = this->getLabel(); + destination_simulationJob.nodalForcedDisplacements = destination_forcedDisplacements; } SimulationJob getCopy() const { SimulationJob jobCopy; @@ -377,7 +377,7 @@ json[jsonLabels.meshFilename]= std::filesystem::relative(std::filesystem::path(m return; } std::vector> nodeColors(pMesh->VN()); - for (auto fixedVertex : constrainedVertices) { + for (const auto &fixedVertex : constrainedVertices) { const bool hasRotationalDoFConstrained = fixedVertex.second.contains(3) || fixedVertex.second.contains(4) || fixedVertex.second.contains(5); @@ -393,7 +393,8 @@ json[jsonLabels.meshFilename]= std::filesystem::relative(std::filesystem::path(m } } if (!nodalForcedDisplacements.empty()) { - for (std::pair viDisplPair : nodalForcedDisplacements) { + for (const std::pair &viDisplPair : + nodalForcedDisplacements) { const VertexIndex vi = viDisplPair.first; nodeColors[vi][0] += 1; nodeColors[vi][0] /= 2; @@ -712,21 +713,21 @@ struct SimulationResults // /*polyscope::getNextUniqueColor()*/ glm::vec3(0, 0, 1)); // //} pJob->registerForDrawing(getLabel()); - static bool wasExecuted = false; - if (!wasExecuted) { - std::function callback = [&]() { - static bool showFrames = shouldShowFrames; + // static bool wasExecuted =false; + // if (!wasExecuted) { + // std::function callback = [&]() { + // static bool showFrames = shouldShowFrames; - if (ImGui::Checkbox("Show Frames", &showFrames) && showFrames) { - polyscopeHandle_frameX->setEnabled(showFrames); - polyscopeHandle_frameY->setEnabled(showFrames); - polyscopeHandle_frameZ->setEnabled(showFrames); - } - }; + // if (ImGui::Checkbox("Show Frames", &showFrames) && showFrames) { + // polyscopeHandle_frameX->setEnabled(showFrames); + // polyscopeHandle_frameY->setEnabled(showFrames); + // polyscopeHandle_frameZ->setEnabled(showFrames); + // } + // }; - PolyscopeInterface::addUserCallback(callback); - wasExecuted = true; - } + // PolyscopeInterface::addUserCallback(callback); + // wasExecuted = true; + // } return polyscopeHandle_deformedEdmeMesh; } #endif @@ -738,7 +739,7 @@ private: //load job //Use the first .eigenBin file for loading the displacements for (auto const &entry : std::filesystem::recursive_directory_iterator(loadFromPath)) { - if (filesystem::is_regular_file(entry) && entry.path().extension() == ".eigenBin") { + if (std::filesystem::is_regular_file(entry) && entry.path().extension() == ".eigenBin") { Eigen::MatrixXd displacements_eigen; Eigen::read_binary(entry.path().string(), displacements_eigen); displacements = Utilities::fromEigenMatrix(displacements_eigen); diff --git a/simulationhistoryplotter.hpp b/simulationhistoryplotter.hpp index 1c803f3..f49612b 100755 --- a/simulationhistoryplotter.hpp +++ b/simulationhistoryplotter.hpp @@ -15,27 +15,25 @@ struct SimulationResultsReporter { void writeStatistics(const SimulationResults &results, const std::string &reportFolderPath) { + std::ofstream file; + file.open(std::filesystem::path(reportFolderPath).append("results.txt").string()); - ofstream file; - file.open( - std::filesystem::path(reportFolderPath).append("results.txt").string()); + const size_t numberOfSteps = results.history.numberOfSteps; + file << "Number of steps " << numberOfSteps << "\n"; - const size_t numberOfSteps = results.history.numberOfSteps; - file << "Number of steps " << numberOfSteps << "\n"; + // file << "Force threshold used " << 1000 << "\n"; - // file << "Force threshold used " << 1000 << "\n"; - - // assert(numberOfSteps == results.history.potentialEnergy.size() && - // numberOfSteps == results.history.residualForces.size()); - // Write kinetic energies - const SimulationHistory &history = results.history; - if (!history.kineticEnergy.empty()) { - file << "Kinetic energies" - << "\n"; - for (size_t step = 0; step < numberOfSteps; step++) { - file << history.kineticEnergy[step] << "\n"; - } - file << "\n"; + // assert(numberOfSteps == results.history.potentialEnergy.size() && + // numberOfSteps == results.history.residualForces.size()); + // Write kinetic energies + const SimulationHistory &history = results.history; + if (!history.kineticEnergy.empty()) { + file << "Kinetic energies" + << "\n"; + for (size_t step = 0; step < numberOfSteps; step++) { + file << history.kineticEnergy[step] << "\n"; + } + file << "\n"; } if (!history.logResidualForces.empty()) { diff --git a/simulationmesh.cpp b/simulationmesh.cpp index 446a3eb..425a196 100755 --- a/simulationmesh.cpp +++ b/simulationmesh.cpp @@ -284,85 +284,85 @@ std::vector SimulationMesh::getBeamMaterial() { return beamMaterial; } -bool SimulationMesh::load(const string &plyFilename) { - this->Clear(); - // assert(plyFileHasAllRequiredFields(plyFilename)); - // Load the ply file - // VCGEdgeMesh::PerEdgeAttributeHandle handleBeamDimensions = - // vcg::tri::Allocator::AddPerEdgeAttribute< - // CrossSectionType>(*this, plyPropertyBeamDimensionsID); - // VCGEdgeMesh::PerEdgeAttributeHandle handleBeamMaterial = - // vcg::tri::Allocator::AddPerEdgeAttribute( - // *this, plyPropertyBeamMaterialID); - // nanoply::NanoPlyWrapper::CustomAttributeDescriptor - // customAttrib; - // customAttrib.GetMeshAttrib(plyFilename); - // customAttrib.AddEdgeAttribDescriptor( - // plyPropertyBeamDimensionsID, nanoply::NNP_LIST_INT8_FLOAT64, nullptr); - // /*FIXME: Since I allow CrossSectionType to take two types I should export the - // * type as well such that that when loaded the correct type of cross section - // * is used. - // */ - // customAttrib.AddEdgeAttribDescriptor( - // plyPropertyBeamMaterialID, nanoply::NNP_LIST_INT8_FLOAT64, nullptr); +bool SimulationMesh::load(const std::string &plyFilename) +{ + this->Clear(); + // assert(plyFileHasAllRequiredFields(plyFilename)); + // Load the ply file + // VCGEdgeMesh::PerEdgeAttributeHandle handleBeamDimensions = + // vcg::tri::Allocator::AddPerEdgeAttribute< + // CrossSectionType>(*this, plyPropertyBeamDimensionsID); + // VCGEdgeMesh::PerEdgeAttributeHandle handleBeamMaterial = + // vcg::tri::Allocator::AddPerEdgeAttribute( + // *this, plyPropertyBeamMaterialID); + // nanoply::NanoPlyWrapper::CustomAttributeDescriptor + // customAttrib; + // customAttrib.GetMeshAttrib(plyFilename); + // customAttrib.AddEdgeAttribDescriptor( + // plyPropertyBeamDimensionsID, nanoply::NNP_LIST_INT8_FLOAT64, nullptr); + // /*FIXME: Since I allow CrossSectionType to take two types I should export the + // * type as well such that that when loaded the correct type of cross section + // * is used. + // */ + // customAttrib.AddEdgeAttribDescriptor( + // plyPropertyBeamMaterialID, nanoply::NNP_LIST_INT8_FLOAT64, nullptr); - // VCGEdgeMesh::PerEdgeAttributeHandle> - // handleBeamProperties = - // vcg::tri::Allocator::AddPerEdgeAttribute< - // std::array>(*this, plyPropertyBeamProperties); - // customAttrib.AddEdgeAttribDescriptor, double, 6>( - // plyPropertyBeamProperties, nanoply::NNP_LIST_INT8_FLOAT64, nullptr); + // VCGEdgeMesh::PerEdgeAttributeHandle> + // handleBeamProperties = + // vcg::tri::Allocator::AddPerEdgeAttribute< + // std::array>(*this, plyPropertyBeamProperties); + // customAttrib.AddEdgeAttribDescriptor, double, 6>( + // plyPropertyBeamProperties, nanoply::NNP_LIST_INT8_FLOAT64, nullptr); - // VCGEdgeMesh::PerEdgeAttributeHandle - // handleBeamRigidityContants; - // customAttrib.AddEdgeAttribDescriptor( - // plyPropertyBeamRigidityConstantsID, nanoply::NNP_LIST_INT8_FLOAT32, - // nullptr); - // unsigned int mask = 0; - // mask |= nanoply::NanoPlyWrapper::IO_VERTCOORD; - // mask |= nanoply::NanoPlyWrapper::IO_VERTNORMAL; - // mask |= nanoply::NanoPlyWrapper::IO_EDGEINDEX; - // mask |= nanoply::NanoPlyWrapper::IO_EDGEATTRIB; - // mask |= nanoply::NanoPlyWrapper::IO_MESHATTRIB; - if (!VCGEdgeMesh::load(plyFilename)) { - return false; - } - - // elements = vcg::tri::Allocator::GetPerEdgeAttribute( - // *this, std::string("Elements")); - // nodes = vcg::tri::Allocator::GetPerVertexAttribute( - // *this, std::string("Nodes")); - vcg::tri::UpdateTopology::VertexEdge(*this); - initializeNodes(); - initializeElements(); - setBeamMaterial(0.3, 1 * 1e9); - updateEigenEdgeAndVertices(); - - // if (!handleBeamProperties._handle->data.empty()) { - // for (size_t ei = 0; ei < EN(); ei++) { - // elements[ei] = - // Element::Properties(handleBeamProperties[ei]); - // elements[ei].updateRigidity(); - // } - // } - // for (size_t ei = 0; ei < EN(); ei++) { - // elements[ei].setDimensions(handleBeamDimensions[ei]); - // elements[ei].setMaterial(handleBeamMaterial[ei]); - // elements[ei].updateRigidity(); - // } - - bool normalsAreAbsent = vert[0].cN().Norm() < 0.000001; - if (normalsAreAbsent) { - CoordType normalVector(0, 0, 1); - std::cout << "Warning: Normals are missing from " << plyFilename - << ". Added normal vector:" << toString(normalVector) - << std::endl; - for (auto &v : vert) { - v.N() = normalVector; + // VCGEdgeMesh::PerEdgeAttributeHandle + // handleBeamRigidityContants; + // customAttrib.AddEdgeAttribDescriptor( + // plyPropertyBeamRigidityConstantsID, nanoply::NNP_LIST_INT8_FLOAT32, + // nullptr); + // unsigned int mask = 0; + // mask |= nanoply::NanoPlyWrapper::IO_VERTCOORD; + // mask |= nanoply::NanoPlyWrapper::IO_VERTNORMAL; + // mask |= nanoply::NanoPlyWrapper::IO_EDGEINDEX; + // mask |= nanoply::NanoPlyWrapper::IO_EDGEATTRIB; + // mask |= nanoply::NanoPlyWrapper::IO_MESHATTRIB; + if (!VCGEdgeMesh::load(plyFilename)) { + return false; } - } - return true; + // elements = vcg::tri::Allocator::GetPerEdgeAttribute( + // *this, std::string("Elements")); + // nodes = vcg::tri::Allocator::GetPerVertexAttribute( + // *this, std::string("Nodes")); + vcg::tri::UpdateTopology::VertexEdge(*this); + initializeNodes(); + initializeElements(); + setBeamMaterial(0.3, 1 * 1e9); + updateEigenEdgeAndVertices(); + + // if (!handleBeamProperties._handle->data.empty()) { + // for (size_t ei = 0; ei < EN(); ei++) { + // elements[ei] = + // Element::Properties(handleBeamProperties[ei]); + // elements[ei].updateRigidity(); + // } + // } + // for (size_t ei = 0; ei < EN(); ei++) { + // elements[ei].setDimensions(handleBeamDimensions[ei]); + // elements[ei].setMaterial(handleBeamMaterial[ei]); + // elements[ei].updateRigidity(); + // } + + bool normalsAreAbsent = vert[0].cN().Norm() < 0.000001; + if (normalsAreAbsent) { + CoordType normalVector(0, 0, 1); + std::cout << "Warning: Normals are missing from " << plyFilename + << ". Added normal vector:" << toString(normalVector) << std::endl; + for (auto &v : vert) { + v.N() = normalVector; + } + } + + return true; } bool SimulationMesh::save(const std::string &plyFilename) diff --git a/simulationmesh.hpp b/simulationmesh.hpp index 791e9ea..ecdcb52 100755 --- a/simulationmesh.hpp +++ b/simulationmesh.hpp @@ -34,7 +34,7 @@ public: std::vector getIncidentElements(const VCGEdgeMesh::VertexType &v); - virtual bool load(const string &plyFilename); + virtual bool load(const std::string &plyFilename); std::vector getBeamDimensions(); std::vector getBeamMaterial(); diff --git a/trianglepatterngeometry.cpp b/trianglepatterngeometry.cpp index 74d67be..9d30d9a 100755 --- a/trianglepatterngeometry.cpp +++ b/trianglepatterngeometry.cpp @@ -145,7 +145,7 @@ PatternGeometry::PatternGeometry(PatternGeometry &other) { vcg::tri::UpdateTopology::EdgeEdge(*this); } -bool PatternGeometry::load(const string &plyFilename) +bool PatternGeometry::load(const std::string &plyFilename) { if (!VCGEdgeMesh::load(plyFilename)) { return false; diff --git a/utilities.hpp b/utilities.hpp index ad42f77..c0fee9a 100755 --- a/utilities.hpp +++ b/utilities.hpp @@ -4,9 +4,11 @@ #include #include #include +#include #include #include #include +#include #include #define GET_VARIABLE_NAME(Variable) (#Variable) @@ -140,16 +142,37 @@ struct Vector6d : public std::array { }; namespace Utilities { -//inline std::vector split(std::string text, char delim) -//{ -// std::string line; -// std::vector vec; -// std::stringstream ss(text); -// while (std::getline(ss, line, delim)) { -// vec.push_back(line); -// } -// return vec; -//} +inline bool compareNat(const std::string &a, const std::string &b) +{ + if (a.empty()) + return true; + if (b.empty()) + return false; + if (std::isdigit(a[0]) && !std::isdigit(b[0])) + return true; + if (!std::isdigit(a[0]) && std::isdigit(b[0])) + return false; + if (!std::isdigit(a[0]) && !std::isdigit(b[0])) { + if (std::toupper(a[0]) == std::toupper(b[0])) + return compareNat(a.substr(1), b.substr(1)); + return (std::toupper(a[0]) < std::toupper(b[0])); + } + + // Both strings begin with digit --> parse both numbers + std::istringstream issa(a); + std::istringstream issb(b); + int ia, ib; + issa >> ia; + issb >> ib; + if (ia != ib) + return ia < ib; + + // Numbers are the same --> remove numbers and recurse + std::string anew, bnew; + std::getline(issa, anew); + std::getline(issb, bnew); + return (compareNat(anew, bnew)); +} inline std::string_view leftTrimSpaces(const std::string_view& str) { @@ -175,6 +198,21 @@ std::string_view trimmedString=str; return trimmedString; } +template +inline void normalize(InputIt itBegin, InputIt itEnd) +{ + const auto squaredSumOfElements = std::accumulate(itBegin, + itEnd, + 0.0, + [](const auto &sum, const auto &el) { + return sum + el * el; + }); + assert(squaredSumOfElements != 0); + std::transform(itBegin, itEnd, itBegin, [&](auto &element) { + return element / std::sqrt(squaredSumOfElements); + }); +} + inline std::vector split(const std::string& text, std::string delim) { std::vector vec; @@ -273,6 +311,19 @@ inline std::vector fromEigenMatrix(const Eigen::MatrixXd &m) // return true; //} +inline std::filesystem::path getFilepathWithExtension(const std::filesystem::path &folderPath, + const std::string &extension) +{ + for (const std::filesystem::directory_entry &dirEntry : + std::filesystem::directory_iterator(folderPath)) { + if (dirEntry.is_regular_file() && std::filesystem::path(dirEntry).extension() == extension) { + return std::filesystem::path(dirEntry); + } + } + + return ""; +} + } // namespace Utilities #ifdef POLYSCOPE_DEFINED @@ -293,7 +344,7 @@ inline void mainCallback() // instead of full width. Must have // matching PopItemWidth() below. - for (std::function userCallback : globalPolyscopeData.userCallbacks) { + for (std::function &userCallback : globalPolyscopeData.userCallbacks) { userCallback(); } @@ -418,17 +469,4 @@ inline size_t computeHashOrdered(const std::vector &v) return std::hash{}(elementsString); } -//inline std::filesystem::path getFilepathWithExtension(const std::filesystem::path &folderPath, -// const std::string &extension) -//{ -// for (const std::filesystem::directory_entry &dirEntry : -// std::filesystem::directory_iterator(folderPath)) { -// if (dirEntry.is_regular_file() && std::filesystem::path(dirEntry).extension() == extension) { -// return std::filesystem::path(dirEntry); -// } -// } - -// return ""; -//} - #endif // UTILITIES_H