Merged potentialEnergy branch

This commit is contained in:
iasonmanolas 2022-01-24 11:46:33 +02:00
commit 872d3b1144
4 changed files with 60 additions and 55 deletions

3
.gitmodules vendored Normal file
View File

@ -0,0 +1,3 @@
[submodule "deps/MySources"]
path = deps/MySources
url = https://gitea-s2i2s.isti.cnr.it/manolas/MySources

1
deps/MySources vendored Submodule

@ -0,0 +1 @@
Subproject commit f79e0acb4f3361939cc495956b7afb52a318d02b

View File

@ -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;
}
@ -140,27 +140,27 @@ double ReducedModelOptimizer::computeError(
const double &scenarioWeight,
const Settings::ObjectiveWeights &objectiveWeights)
{
const double translationalError
= computeDisplacementError(simulationResults_fullPattern.displacements,
simulationResults_reducedPattern.displacements,
reducedToFullInterfaceViMap,
normalizationFactor_translationalDisplacement);
// const double translationalError
// = computeRawTranslationalError(simulationResults_fullPattern.displacements,
// = computeDisplacementError(simulationResults_fullPattern.displacements,
// simulationResults_reducedPattern.displacements,
// reducedToFullInterfaceViMap);
// reducedToFullInterfaceViMap,
// normalizationFactor_translationalDisplacement);
const double translationalError
= computeRawTranslationalError(simulationResults_fullPattern.displacements,
simulationResults_reducedPattern.displacements,
reducedToFullInterfaceViMap);
// std::cout << "normalization factor:" << normalizationFactor_rotationalDisplacement << std::endl;
// std::cout << "trans error:" << translationalError << std::endl;
const double rotationalError
= computeRotationalError(simulationResults_fullPattern.rotationalDisplacementQuaternion,
simulationResults_reducedPattern.rotationalDisplacementQuaternion,
reducedToFullInterfaceViMap,
normalizationFactor_rotationalDisplacement);
// const double rotationalError
// = computeRawRotationalError(simulationResults_fullPattern.rotationalDisplacementQuaternion,
// = computeRotationalError(simulationResults_fullPattern.rotationalDisplacementQuaternion,
// simulationResults_reducedPattern.rotationalDisplacementQuaternion,
// reducedToFullInterfaceViMap);
// reducedToFullInterfaceViMap,
// normalizationFactor_rotationalDisplacement);
const double rotationalError
= computeRawRotationalError(simulationResults_fullPattern.rotationalDisplacementQuaternion,
simulationResults_reducedPattern.rotationalDisplacementQuaternion,
reducedToFullInterfaceViMap);
// std::cout << "rot error:" << rotationalError<< std::endl;
const double simulationError = objectiveWeights.translational * translationalError
+ objectiveWeights.rotational * rotationalError;
@ -186,18 +186,16 @@ struct EnsmallenOptimizationObjective
//#ifdef POLYSCOPE_DEFINED
// std::cout << "Out of range" << std::endl;
//#endif
// return std::numeric_limits<double>::max();
// return 1e7;
// return 1000;
x[xi] = global.xMax[xi];
return std::numeric_limits<double>::max();
// return 1e10;
// x[xi] = global.xMax[xi];
} else if (x[xi] < global.xMin[xi]) {
//#ifdef POLYSCOPE_DEFINED
// std::cout << "Out of range" << std::endl;
//#endif
// return std::numeric_limits<double>::max();
// return 1e7;
// return 1000;
x[xi] = global.xMin[xi];
return std::numeric_limits<double>::max();
// return 1e10;
// x[xi] = global.xMin[xi];
}
}
return ReducedModelOptimizer::objective(x);
@ -314,7 +312,7 @@ double ReducedModelOptimizer::objective(const std::vector<double> &x)
// double totalError = std::accumulate(simulationErrorsPerScenario.begin(),
// simulationErrorsPerScenario.end(),
// 0);
double totalError = std::reduce(std::execution::par_unseq,
const double totalError = std::reduce(std::execution::par_unseq,
simulationErrorsPerScenario.begin(),
simulationErrorsPerScenario.end());
// std::cout << totalError << std::endl;
@ -1099,12 +1097,12 @@ void ReducedModelOptimizer::runOptimization(const Settings &settings,
// ens::CNE optimizer;
// ens::DE optimizer;
// ens::SPSA optimizer;
// arma::mat xMin_arma(global.xMin);
// arma::mat xMax_arma(global.xMax);
arma::mat xMin_arma(global.xMin);
arma::mat xMax_arma(global.xMax);
// ens::LBestPSO optimizer(64, xMin_arma, xMax_arma, 1000);
ens::LBestPSO optimizer(64,
1,
1,
ens::LBestPSO optimizer(200,
xMin_arma,
xMax_arma,
settings.numberOfFunctionCalls,
500,
settings.solverAccuracy,
@ -1112,13 +1110,13 @@ void ReducedModelOptimizer::runOptimization(const Settings &settings,
2.35);
// ens::LBestPSO optimizer;
const double minima = optimizer.Optimize(optimizationFunction, x);
for (int xIndex = 0; xIndex < parameterGroup.size(); xIndex++) {
if (x[xIndex] > global.xMax[xIndex]) {
x[xIndex] = global.xMax[xIndex];
} else if (x[xIndex] < global.xMin[xIndex]) {
x[xIndex] = global.xMin[xIndex];
}
}
// for (int xIndex = 0; xIndex < parameterGroup.size(); xIndex++) {
// if (x[xIndex] > global.xMax[xIndex]) {
// x[xIndex] = global.xMax[xIndex];
// } else if (x[xIndex] < global.xMin[xIndex]) {
// x[xIndex] = global.xMin[xIndex];
// }
// }
#ifdef POLYSCOPE_DEFINED
std::cout << "optimal x:"
@ -1327,14 +1325,14 @@ void ReducedModelOptimizer::constructDomeSimulationScenario(
= Eigen::Vector3d(-interfaceVector[0],
-interfaceVector[1],
0)
* 0.005
* 0.001
* std::abs(
forceMagnitude); //NOTE:Should the forced displacement change relatively to the magnitude?
// * std::abs(forceMagnitude / maxForceMagnitude_dome);
job.nodalForcedDisplacements[viPair.second] = Eigen::Vector3d(interfaceVector[0],
interfaceVector[1],
0)
* 0.005 * std::abs(forceMagnitude);
* 0.001 * std::abs(forceMagnitude);
// * std::abs(forceMagnitude / maxForceMagnitude_dome);
// CoordType v = (pMesh->vert[viPair.first].cP() - pMesh->vert[viPair.second].cP())
// ^ CoordType(0, 0, -1).Normalize();
@ -1565,19 +1563,22 @@ double ReducedModelOptimizer::computeFullPatternMaxSimulationForce(
global.desiredMaxDisplacementValue = 0.03 * baseTriangleHeight;
break;
case Shear:
global.desiredMaxDisplacementValue = 0.03 * baseTriangleHeight;
global.desiredMaxDisplacementValue = 0.04 * baseTriangleHeight;
break;
case Bending:
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first;
global.desiredMaxDisplacementValue = 0.05 * baseTriangleHeight;
break;
case Dome:
global.desiredMaxRotationAngle = vcg::math::ToRad(35.0);
global.desiredMaxRotationAngle = vcg::math::ToRad(20.0);
objectiveFunction = &fullPatternMaxSimulationForceRotationalObjective;
forceMagnitudeEpsilon *= 1e-2;
objectiveEpsilon = vcg::math::ToRad(3.0);
objectiveEpsilon = vcg::math::ToRad(1.0);
forceMagnitude = 0.6;
break;
case Saddle:
global.interfaceViForComputingScenarioError = global.fullPatternInterfaceViPairs[0].first;
global.desiredMaxDisplacementValue = 0.05 * baseTriangleHeight;
break;
}
@ -1801,9 +1802,9 @@ void ReducedModelOptimizer::computeObjectiveValueNormalizationFactors()
continue;
}
}
angularDistanceSum += global.fullPatternResults[simulationScenarioIndex]
angularDistanceSum += std::abs(global.fullPatternResults[simulationScenarioIndex]
.rotationalDisplacementQuaternion[fullPatternVi]
.angularDistance(Eigen::Quaterniond::Identity());
.angularDistance(Eigen::Quaterniond::Identity()));
}
fullPatternTranslationalDisplacementNormSum[simulationScenarioIndex]

View File

@ -66,7 +66,7 @@ public:
// inline constexpr static ParameterLabels parameterLabels();
inline static std::array<std::string, ReducedModelOptimization::NumberOfOptimizationVariables>
parameterLabels = {"R", "A", "I2", "I3", "J", "Theta", "R"};
parameterLabels = {"E", "A", "I2", "I3", "J", "Theta", "R"};
constexpr static std::array<double, ReducedModelOptimization::NumberOfBaseSimulationScenarios>
simulationScenariosResolution = {11, 11, 20, 20, 20};
constexpr static std::array<double, ReducedModelOptimization::NumberOfBaseSimulationScenarios>