Added linear and non linear chronos euler simulation model classes in order to be able to out-of-the-box use them with a factory class

This commit is contained in:
iasonmanolas 2022-08-08 12:02:14 +03:00
parent 9fd536ccf3
commit 68d5655214
6 changed files with 264 additions and 229 deletions

View File

@ -1,5 +1,4 @@
#include "chronoseulersimulationmodel.hpp" #include "chronoseulersimulationmodel.hpp"
#include "chrono/physics/ChLoadContainer.h"
#include <chrono/fea/ChBeamSectionEuler.h> #include <chrono/fea/ChBeamSectionEuler.h>
#include <chrono/fea/ChBuilderBeam.h> #include <chrono/fea/ChBuilderBeam.h>
#include <chrono/fea/ChLoadsBeam.h> #include <chrono/fea/ChLoadsBeam.h>
@ -8,6 +7,7 @@
#include <chrono/physics/ChBody.h> #include <chrono/physics/ChBody.h>
#include <chrono/physics/ChSystemSMC.h> #include <chrono/physics/ChSystemSMC.h>
#include <chrono/solver/ChIterativeSolverLS.h> #include <chrono/solver/ChIterativeSolverLS.h>
#include "chrono/physics/ChLoadContainer.h"
#include <chrono/assets/ChVisualization.h> #include <chrono/assets/ChVisualization.h>
#include <chrono/fea/ChElementBeamEuler.h> #include <chrono/fea/ChElementBeamEuler.h>
@ -16,242 +16,268 @@
using namespace chrono; using namespace chrono;
using namespace chrono::fea; using namespace chrono::fea;
std::shared_ptr<ChMesh> ChronosEulerSimulationModel::convertToChronosMesh_Euler( std::shared_ptr<ChMesh> ChronosEulerSimulationModel::convertToChronosMesh_Euler(
const std::shared_ptr<SimulationMesh> &pMesh, const std::shared_ptr<SimulationEdgeMesh>& pMesh,
std::vector<std::shared_ptr<ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes) std::vector<std::shared_ptr<ChNodeFEAxyzrot>>&
{ edgeMeshVertsToChronosNodes) {
auto mesh_chronos = chrono_types::make_shared<ChMesh>(); auto mesh_chronos = chrono_types::make_shared<ChMesh>();
edgeMeshVertsToChronosNodes.clear(); edgeMeshVertsToChronosNodes.clear();
edgeMeshVertsToChronosNodes.resize(pMesh->VN(), nullptr); edgeMeshVertsToChronosNodes.resize(pMesh->VN(), nullptr);
//add nodes // add nodes
for (int vi = 0; vi < pMesh->VN(); vi++) { for (int vi = 0; vi < pMesh->VN(); vi++) {
const auto &vertex = pMesh->vert[vi]; const auto& vertex = pMesh->vert[vi];
ChVector<> vertexPos(vertex.cP()[0], vertex.cP()[1], vertex.cP()[2]); ChVector<> vertexPos(vertex.cP()[0], vertex.cP()[1], vertex.cP()[2]);
edgeMeshVertsToChronosNodes[vi] = chrono_types::make_shared<ChNodeFEAxyzrot>( edgeMeshVertsToChronosNodes[vi] =
ChFrame<>(vertexPos)); chrono_types::make_shared<ChNodeFEAxyzrot>(ChFrame<>(vertexPos));
mesh_chronos->AddNode(edgeMeshVertsToChronosNodes[vi]); mesh_chronos->AddNode(edgeMeshVertsToChronosNodes[vi]);
} }
//add elements // add elements
ChBuilderBeamEuler builder; ChBuilderBeamEuler builder;
for (int ei = 0; ei < pMesh->EN(); ei++) { for (int ei = 0; ei < pMesh->EN(); ei++) {
const SimulationMesh::EdgeType &edge = pMesh->edge[ei]; const SimulationEdgeMesh::EdgeType& edge = pMesh->edge[ei];
//define end-points // define end-points
const auto vi0 = pMesh->getIndex(edge.cV(0)); const auto vi0 = pMesh->getIndex(edge.cV(0));
const auto vi1 = pMesh->getIndex(edge.cV(1)); const auto vi1 = pMesh->getIndex(edge.cV(1));
//define cross section // define cross section
const Element &element = pMesh->elements[ei]; const Element& element = pMesh->elements[ei];
const double beam_wz = element.dimensions.getDim1(); const double beam_wz = element.dimensions.getDim1();
const double beam_wy = element.dimensions.getDim2(); const double beam_wy = element.dimensions.getDim2();
const double E = element.material.youngsModulus; const double E = element.material.youngsModulus;
// const double poisson = element.material.poissonsRatio; // const double poisson = element.material.poissonsRatio;
const double density = 1e0; const double density = 1e0;
// auto msection = chrono_types::make_shared<ChBeamSectionEulerAdvanced>(); // auto msection =
auto msection = chrono_types::make_shared<ChBeamSectionEulerEasyRectangular>( // chrono_types::make_shared<ChBeamSectionEulerAdvanced>();
auto msection =
chrono_types::make_shared<ChBeamSectionEulerEasyRectangular>(
beam_wy, beam_wz, E, element.material.G, density); beam_wy, beam_wz, E, element.material.G, density);
// msection->SetDensity(density); msection->SetArea(element.dimensions.A);
// msection->SetYoungModulus(E); msection->SetIyy(element.dimensions.inertia.I2);
// msection->SetGwithPoissonRatio(poisson); msection->SetIzz(element.dimensions.inertia.I3);
// // msection->SetBeamRaleyghDamping(0.0); msection->SetJ(element.dimensions.inertia.J);
// msection->SetAsRectangularSection(beam_wy, beam_wz); // msection->SetDensity(density);
builder // msection->SetYoungModulus(E);
.BuildBeam(mesh_chronos, // the mesh where to put the created nodes and elements // msection->SetGwithPoissonRatio(poisson);
msection, // the ChBeamSectionEuler to use for the ChElementBeamEuler elements // // msection->SetBeamRaleyghDamping(0.0);
1, // the number of ChElementBeamEuler to create // msection->SetAsRectangularSection(beam_wy, beam_wz);
edgeMeshVertsToChronosNodes[vi0], // the 'A' point in space (beginning of beam) builder.BuildBeam(
edgeMeshVertsToChronosNodes[vi1], // the 'B' point in space (end of beam) mesh_chronos, // the mesh where to put the created nodes and elements
ChVector<>(0, 0, 1) msection, // the ChBeamSectionEuler to use for the ChElementBeamEuler
// ChVector<>(0, cos(rot_rad), sin(rot_rad)) // elements
); // the 'Y' up direction of the section for the beam 1, // the number of ChElementBeamEuler to create
} edgeMeshVertsToChronosNodes[vi0], // the 'A' point in space (beginning
// of beam)
edgeMeshVertsToChronosNodes[vi1], // the 'B' point in space (end of
// beam)
ChVector<>(0, 0, 1)
// ChVector<>(0, cos(rot_rad), sin(rot_rad))
); // the 'Y' up direction of the section for the beam
}
return mesh_chronos; return mesh_chronos;
} }
ChronosEulerSimulationModel::ChronosEulerSimulationModel() {} ChronosEulerSimulationModel::ChronosEulerSimulationModel() {}
//SimulationResults ChronosEulerSimulationModel::executeSimulation( // SimulationResults ChronosEulerSimulationModel::executeSimulation(
// const std::shared_ptr<SimulationJob> &pJob) // const std::shared_ptr<SimulationJob> &pJob)
//{} //{}
//chrono::ChSystemSMC convertToChronosSystem(const std::shared_ptr<SimulationJob> &pJob) // chrono::ChSystemSMC convertToChronosSystem(const
// std::shared_ptr<SimulationJob> &pJob)
//{ //{
// chrono::ChSystemSMC my_system; // chrono::ChSystemSMC my_system;
//} //}
void ChronosEulerSimulationModel::parseForces( void ChronosEulerSimulationModel::parseForces(
const std::shared_ptr<chrono::fea::ChMesh> &mesh_chronos, const std::shared_ptr<chrono::fea::ChMesh>& mesh_chronos,
const std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>> &edgeMeshVertsToChronoNodes, const std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>>&
const std::unordered_map<VertexIndex, Vector6d> &nodalExternalForces) edgeMeshVertsToChronoNodes,
{ const std::unordered_map<VertexIndex, Vector6d>& nodalExternalForces) {
mesh_chronos->SetAutomaticGravity(false); mesh_chronos->SetAutomaticGravity(false);
for (const std::pair<VertexIndex, Vector6d> &externalForce : nodalExternalForces) { for (const std::pair<VertexIndex, Vector6d>& externalForce :
const int &forceVi = externalForce.first; nodalExternalForces) {
const Vector6d &force = externalForce.second; const int& forceVi = externalForce.first;
edgeMeshVertsToChronoNodes[forceVi]->SetForce(ChVector<>(force[0], force[1], force[2])); const Vector6d& force = externalForce.second;
edgeMeshVertsToChronoNodes[forceVi]->SetTorque(ChVector<>(force[3], force[4], force[5])); edgeMeshVertsToChronoNodes[forceVi]->SetForce(
} ChVector<>(force[0], force[1], force[2]));
edgeMeshVertsToChronoNodes[forceVi]->SetTorque(
ChVector<>(force[3], force[4], force[5]));
}
} }
void ChronosEulerSimulationModel::parseConstrainedVertices( void ChronosEulerSimulationModel::parseConstrainedVertices(
const std::shared_ptr<const SimulationJob> &pJob, const std::shared_ptr<const SimulationJob>& pJob,
const std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>> &edgeMeshVertsToChronoNodes, const std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>>&
chrono::ChSystemSMC &my_system) edgeMeshVertsToChronoNodes,
{ chrono::ChSystemSMC& my_system) {
assert(!edgeMeshVertsToChronoNodes.empty()); assert(!edgeMeshVertsToChronoNodes.empty());
for (const std::pair<VertexIndex, std::unordered_set<int>> &constrainedVertex : for (const std::pair<VertexIndex, std::unordered_set<int>>&
pJob->constrainedVertices) { constrainedVertex : pJob->constrainedVertices) {
const int &constrainedVi = constrainedVertex.first; const int& constrainedVi = constrainedVertex.first;
const std::unordered_set<int> &constrainedDoF = constrainedVertex.second; const std::unordered_set<int>& constrainedDoF = constrainedVertex.second;
// Create a truss // Create a truss
auto truss = chrono_types::make_shared<ChBody>(); auto truss = chrono_types::make_shared<ChBody>();
truss->SetBodyFixed(true); truss->SetBodyFixed(true);
my_system.Add(truss); my_system.Add(truss);
const auto &constrainedChronoNode = edgeMeshVertsToChronoNodes[constrainedVi]; const auto& constrainedChronoNode =
// Create a constraint at the end of the beam edgeMeshVertsToChronoNodes[constrainedVi];
auto constr_a = chrono_types::make_shared<ChLinkMateGeneric>(); // Create a constraint at the end of the beam
constr_a->SetConstrainedCoords(constrainedDoF.contains(0), auto constr_a = chrono_types::make_shared<ChLinkMateGeneric>();
constrainedDoF.contains(1), constr_a->SetConstrainedCoords(
constrainedDoF.contains(2), constrainedDoF.contains(0), constrainedDoF.contains(1),
constrainedDoF.contains(3), constrainedDoF.contains(2), constrainedDoF.contains(3),
constrainedDoF.contains(4), constrainedDoF.contains(4), constrainedDoF.contains(5));
constrainedDoF.contains(5)); constr_a->Initialize(constrainedChronoNode, truss, false,
constr_a->Initialize(constrainedChronoNode, constrainedChronoNode->Frame(),
truss, constrainedChronoNode->Frame());
false, // const auto frameNode = constrainedChronoNode->Frame();
constrainedChronoNode->Frame(), my_system.Add(constr_a);
constrainedChronoNode->Frame());
// const auto frameNode = constrainedChronoNode->Frame();
my_system.Add(constr_a);
// edgeMeshVertsToChronoNodes[constrainedVi]->SetFixed(true); // edgeMeshVertsToChronoNodes[constrainedVi]->SetFixed(true);
// if (vertexIsFullyConstrained) { // if (vertexIsFullyConstrained) {
// } else { // } else {
// std::cerr << "Currently only rigid vertices are handled." << std::endl; // std::cerr << "Currently only rigid vertices are handled." <<
// // SimulationResults simulationResults; // std::endl;
// // simulationResults.converged = false; // // SimulationResults simulationResults;
// // assert(false); // // simulationResults.converged = false;
// // return simulationResults; // // assert(false);
// } // // return simulationResults;
} // }
}
} }
SimulationResults ChronosEulerSimulationModel::executeSimulation( SimulationResults ChronosEulerSimulationModel::executeSimulation(
const std::shared_ptr<SimulationJob> &pJob) const std::shared_ptr<SimulationJob>& pJob) {
{ assert(pJob->pMesh->VN() != 0);
// assert(pJob->pMesh->VN() != 0); // const bool structureInitialized = mesh_chronos != nullptr;
// const bool structureInitialized = mesh_chronos != nullptr; // const bool wasInitializedWithDifferentStructure =
// const bool wasInitializedWithDifferentStructure = structureInitialized // structureInitialized &&
// && (pJob->pMesh->EN() // (pJob->pMesh->EN() != mesh_chronos->GetNelements() ||
// != mesh_chronos->GetNelements() // pJob->pMesh->VN() != mesh_chronos->GetNnodes());
// || pJob->pMesh->VN() // if (!structureInitialized || wasInitializedWithDifferentStructure) {
// != mesh_chronos->GetNnodes()); setStructure(pJob->pMesh);
// if (!structureInitialized || wasInitializedWithDifferentStructure) { // }
setStructure(pJob->pMesh); chrono::ChSystemSMC my_system;
// } // chrono::irrlicht::ChIrrApp application(&my_system,
chrono::ChSystemSMC my_system; // L"Irrlicht FEM visualization",
// chrono::irrlicht::ChIrrApp application(&my_system, // irr::core::dimension2d<irr::u32>(800,
// L"Irrlicht FEM visualization", // 600),
// irr::core::dimension2d<irr::u32>(800, 600), // chrono::irrlicht::VerticalDir::Z,
// chrono::irrlicht::VerticalDir::Z, // false,
// false, // true);
// true); // const std::string chronoDataFolderPath =
// const std::string chronoDataFolderPath = "/home/iason/Coding/build/external " // "/home/iason/Coding/build/external "
// "dependencies/CHRONO-src/data/"; // "dependencies/CHRONO-src/data/";
// application.AddTypicalLogo(chronoDataFolderPath + "logo_chronoengine_alpha.png"); // application.AddTypicalLogo(chronoDataFolderPath +
// application.AddTypicalSky(chronoDataFolderPath + "skybox/"); // "logo_chronoengine_alpha.png");
// application.AddTypicalLights(); // application.AddTypicalSky(chronoDataFolderPath + "skybox/");
// application.AddTypicalCamera(irr::core::vector3df(0, (irr::f32) 0.6, -1)); // application.AddTypicalLights();
// my_system.SetTimestepperType(chrono::ChTimestepper::Type::EULER_IMPLICIT_LINEARIZED); // application.AddTypicalCamera(irr::core::vector3df(0, (irr::f32) 0.6,
//parse forces // -1));
parseForces(mesh_chronos, edgeMeshVertsToChronoNodes, pJob->nodalExternalForces); // my_system.SetTimestepperType(chrono::ChTimestepper::Type::EULER_IMPLICIT_LINEARIZED);
//parse constrained vertices // parse forces
parseConstrainedVertices(pJob, edgeMeshVertsToChronoNodes, my_system); parseForces(mesh_chronos, edgeMeshVertsToChronoNodes,
// std::dynamic_pointer_cast<std::shared_ptr<ChNodeFEAxyzrot>>(mesh_chronos->GetNode(1)) pJob->nodalExternalForces);
// ->SetFixed(true); // parse constrained vertices
// and load containers must be added to your system parseConstrainedVertices(pJob, edgeMeshVertsToChronoNodes, my_system);
// auto mvisualizemesh = chrono_types::make_shared<ChVisualizationFEAmesh>(*(mesh_chronos.get())); // std::dynamic_pointer_cast<std::shared_ptr<ChNodeFEAxyzrot>>(mesh_chronos->GetNode(1))
// mvisualizemesh->SetFEMdataType(ChVisualizationFEAmesh::E_PLOT_NODE_DISP_NORM); // ->SetFixed(true);
// mvisualizemesh->SetColorscaleMinMax(0.0, 5.50); // and load containers must be added to your system
// mvisualizemesh->SetShrinkElements(false, 0.85); // auto mvisualizemesh =
// mvisualizemesh->SetSmoothFaces(false); // chrono_types::make_shared<ChVisualizationFEAmesh>(*(mesh_chronos.get()));
// mesh_chronos->AddAsset(mvisualizemesh); // mvisualizemesh->SetFEMdataType(ChVisualizationFEAmesh::E_PLOT_NODE_DISP_NORM);
// mvisualizemesh->SetColorscaleMinMax(0.0, 5.50);
// mvisualizemesh->SetShrinkElements(false, 0.85);
// mvisualizemesh->SetSmoothFaces(false);
// mesh_chronos->AddAsset(mvisualizemesh);
// application.AssetBindAll(); // application.AssetBindAll();
// application.AssetUpdateAll(); // application.AssetUpdateAll();
my_system.Add(mesh_chronos); my_system.Add(mesh_chronos);
auto solver = chrono_types::make_shared<ChSolverMINRES>(); auto solver = chrono_types::make_shared<ChSolverMINRES>();
my_system.SetSolver(solver); my_system.SetSolver(solver);
solver->SetMaxIterations(1e5); // solver->SetMaxIterations(1e5);
// solver->SetTolerance(1e-12); // solver->SetTolerance(1e-12);
solver->EnableWarmStart(true); // IMPORTANT for convergence when using EULER_IMPLICIT_LINEARIZED solver->EnableWarmStart(
solver->EnableDiagonalPreconditioner(true); true); // IMPORTANT for convergence when using EULER_IMPLICIT_LINEARIZED
my_system.SetSolverForceTolerance(1e-9); solver->EnableDiagonalPreconditioner(true);
solver->SetVerbose(false); my_system.SetSolverForceTolerance(1e-9);
solver->SetVerbose(false);
SimulationResults simulationResults; SimulationResults simulationResults;
if (settings.analysisType == Settings::AnalysisType::Linear) {
simulationResults.converged = my_system.DoStaticLinear();
} else {
simulationResults.converged = my_system.DoStaticNonlinear(); simulationResults.converged = my_system.DoStaticNonlinear();
// simulationResults.converged = my_system.DoStaticLinear(); }
if (!simulationResults.converged) { if (!simulationResults.converged) {
std::cerr << "Simulation failed" << std::endl; std::cerr << "Simulation failed" << std::endl;
assert(false); assert(false);
return simulationResults;
}
// my_system.SetTimestepperType(ChTimestepper::Type::EULER_IMPLICIT_LINEARIZED);
// application.SetTimestep(0.001);
// while (application.GetDevice()->run()) {
// application.BeginScene();
// application.DrawAll();
// application.EndScene();
// }
simulationResults.pJob = pJob;
simulationResults.displacements.resize(pJob->pMesh->VN());
simulationResults.rotationalDisplacementQuaternion.resize(pJob->pMesh->VN());
for (size_t vi = 0; vi < pJob->pMesh->VN(); vi++) {
const auto node_chronos = edgeMeshVertsToChronoNodes[vi];
const auto posDisplacement = node_chronos->Frame().GetPos()
- node_chronos->GetX0().GetPos();
// std::cout << "Node " << vi << " coordinate x= " << node_chronos->Frame().GetPos().x()
// << " y=" << node_chronos->Frame().GetPos().y()
// << " z=" << node_chronos->Frame().GetPos().z() << "\n";
//Translations
simulationResults.displacements[vi][0] = posDisplacement[0];
simulationResults.displacements[vi][1] = posDisplacement[1];
simulationResults.displacements[vi][2] = posDisplacement[2];
//Rotations
chrono::ChQuaternion<double> rotQuat = node_chronos->GetRot();
simulationResults.rotationalDisplacementQuaternion[vi]
= Eigen::Quaternion<double>(rotQuat.e0(), rotQuat.e1(), rotQuat.e2(), rotQuat.e3());
const ChVector<double> eulerAngles = rotQuat.Q_to_Euler123();
// std::cout << "Euler angles:" << eulerAngles << std::endl;
simulationResults.displacements[vi][3] = eulerAngles[0];
simulationResults.displacements[vi][4] = eulerAngles[1];
simulationResults.displacements[vi][5] = eulerAngles[2];
}
simulationResults.simulationModelUsed = label;
return simulationResults; return simulationResults;
}
// VCGEdgeMesh deformedMesh; // my_system.SetTimestepperType(ChTimestepper::Type::EULER_IMPLICIT_LINEARIZED);
// deformedMesh.copy(*(pJob->pMesh)); // application.SetTimestep(0.001);
// for (size_t vi = 0; vi < pJob->pMesh->VN(); vi++) {
// const std::shared_ptr<ChNodeFEAxyzrot> node_chronos = edgeMeshVertsToChronosNodes[vi];
// deformedMesh.vert[vi].P() = CoordType(node_chronos->GetPos()[0],
// node_chronos->GetPos()[1],
// node_chronos->GetPos()[2]);
// }
// deformedMesh.updateEigenEdgeAndVertices(); // while (application.GetDevice()->run()) {
// deformedMesh.setLabel("deformed"); // application.BeginScene();
// deformedMesh.registerForDrawing(); // application.DrawAll();
// polyscope::show(); // application.EndScene();
// return simulationResults; // }
simulationResults.pJob = pJob;
simulationResults.displacements.resize(pJob->pMesh->VN());
simulationResults.rotationalDisplacementQuaternion.resize(pJob->pMesh->VN());
for (size_t vi = 0; vi < pJob->pMesh->VN(); vi++) {
const auto node_chronos = edgeMeshVertsToChronoNodes[vi];
const auto posDisplacement =
node_chronos->Frame().GetPos() - node_chronos->GetX0().GetPos();
// std::cout << "Node " << vi << " coordinate x= " <<
// node_chronos->Frame().GetPos().x()
// << " y=" << node_chronos->Frame().GetPos().y()
// << " z=" << node_chronos->Frame().GetPos().z() <<
// "\n";
// Translations
simulationResults.displacements[vi][0] = posDisplacement[0];
simulationResults.displacements[vi][1] = posDisplacement[1];
simulationResults.displacements[vi][2] = posDisplacement[2];
// Rotations
chrono::ChQuaternion<double> rotQuat = node_chronos->GetRot();
simulationResults.rotationalDisplacementQuaternion[vi] =
Eigen::Quaternion<double>(rotQuat.e0(), rotQuat.e1(), rotQuat.e2(),
rotQuat.e3());
const Eigen::Vector3d eulerAngles =
simulationResults.rotationalDisplacementQuaternion[vi]
.toRotationMatrix()
.eulerAngles(0, 1, 2);
// std::cout << "Euler angles:" << eulerAngles << std::endl;
simulationResults.displacements[vi][3] = eulerAngles[0];
simulationResults.displacements[vi][4] = eulerAngles[1];
simulationResults.displacements[vi][5] = eulerAngles[2];
}
simulationResults.simulationModelUsed = label;
return simulationResults;
// VCGEdgeMesh deformedMesh;
// deformedMesh.copy(*(pJob->pMesh));
// for (size_t vi = 0; vi < pJob->pMesh->VN(); vi++) {
// const std::shared_ptr<ChNodeFEAxyzrot> node_chronos =
// edgeMeshVertsToChronosNodes[vi]; deformedMesh.vert[vi].P() =
// CoordType(node_chronos->GetPos()[0],
// node_chronos->GetPos()[1],
// node_chronos->GetPos()[2]);
// }
// deformedMesh.updateEigenEdgeAndVertices();
// deformedMesh.setLabel("deformed");
// deformedMesh.registerForDrawing();
// polyscope::show();
// return simulationResults;
} }
void ChronosEulerSimulationModel::setStructure(const std::shared_ptr<SimulationMesh> &pMesh) void ChronosEulerSimulationModel::setStructure(
{ const std::shared_ptr<SimulationEdgeMesh>& pMesh) {
mesh_chronos = convertToChronosMesh_Euler(pMesh, edgeMeshVertsToChronoNodes); mesh_chronos = convertToChronosMesh_Euler(pMesh, edgeMeshVertsToChronoNodes);
} }

View File

@ -28,12 +28,19 @@ class ChronosEulerSimulationModel : public SimulationModel
public: public:
ChronosEulerSimulationModel(); ChronosEulerSimulationModel();
SimulationResults executeSimulation(const std::shared_ptr<SimulationJob> &pJob) override; SimulationResults executeSimulation(const std::shared_ptr<SimulationJob> &pJob) override;
void setStructure(const std::shared_ptr<SimulationMesh> &pMesh) override; void setStructure(const std::shared_ptr<SimulationEdgeMesh> &pMesh) override;
static std::shared_ptr<chrono::fea::ChMesh> convertToChronosMesh_Euler( static std::shared_ptr<chrono::fea::ChMesh> convertToChronosMesh_Euler(
const std::shared_ptr<SimulationMesh> &pMesh, const std::shared_ptr<SimulationEdgeMesh> &pMesh,
std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes); std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes);
inline const static std::string label{"Chronos_nonLinear_Euler"}; inline const static std::string label{"Chronos_Euler"};
struct Settings
{
enum AnalysisType { Linear = 0, NonLinear };
AnalysisType analysisType{NonLinear};
};
Settings settings;
}; };
#endif // CHRONOSEULERSIMULATIONMODEL_HPP #endif // CHRONOSEULERSIMULATIONMODEL_HPP

View File

@ -16,7 +16,7 @@
using namespace chrono; using namespace chrono;
using namespace chrono::fea; using namespace chrono::fea;
std::shared_ptr<ChMesh> ChronosIGASimulationModel::convertToChronosMesh_IGA( std::shared_ptr<ChMesh> ChronosIGASimulationModel::convertToChronosMesh_IGA(
const std::shared_ptr<SimulationMesh> &pMesh, const std::shared_ptr<SimulationEdgeMesh> &pMesh,
std::vector<std::shared_ptr<ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes) std::vector<std::shared_ptr<ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes)
{ {
auto mesh_chronos = chrono_types::make_shared<ChMesh>(); auto mesh_chronos = chrono_types::make_shared<ChMesh>();
@ -34,7 +34,7 @@ std::shared_ptr<ChMesh> ChronosIGASimulationModel::convertToChronosMesh_IGA(
//add elements //add elements
ChBuilderBeamIGA builder; ChBuilderBeamIGA builder;
for (int ei = 0; ei < pMesh->EN(); ei++) { for (int ei = 0; ei < pMesh->EN(); ei++) {
const SimulationMesh::EdgeType &edge = pMesh->edge[ei]; const SimulationEdgeMesh::EdgeType &edge = pMesh->edge[ei];
//define end-points //define end-points
const auto vi0 = pMesh->getIndex(edge.cV(0)); const auto vi0 = pMesh->getIndex(edge.cV(0));
const auto vi1 = pMesh->getIndex(edge.cV(1)); const auto vi1 = pMesh->getIndex(edge.cV(1));
@ -54,15 +54,16 @@ std::shared_ptr<ChMesh> ChronosIGASimulationModel::convertToChronosMesh_IGA(
// msection->SetGwithPoissonRatio(poisson); // msection->SetGwithPoissonRatio(poisson);
// // msection->SetBeamRaleyghDamping(0.0); // // msection->SetBeamRaleyghDamping(0.0);
// msection->SetAsRectangularSection(beam_wy, beam_wz); // msection->SetAsRectangularSection(beam_wy, beam_wz);
builder.BuildBeam( builder
mesh_chronos, // the mesh where to put the created nodes and elements .BuildBeam(mesh_chronos, // the mesh where to put the created nodes and elements
msection, // the ChBeamSectionEuler to use for the ChElementBeamEuler elements msection, // the ChBeamSectionEuler to use for the ChElementBeamEuler elements
4, // the number of ChElementBeamEuler to create 1, // the number of ChElementBeamEuler to create
edgeMeshVertsToChronosNodes[vi0]->GetPos(), // the 'A' point in space (beginning of beam) edgeMeshVertsToChronosNodes[vi0], // the 'A' point in space (beginning of beam)
edgeMeshVertsToChronosNodes[vi1]->GetPos(), // the 'B' point in space (end of beam) edgeMeshVertsToChronosNodes[vi1], // the 'B' point in space (end of beam)
ChVector<>(0, 0, 1) ChVector<>(0, 0, 1),
// ChVector<>(0, cos(rot_rad), sin(rot_rad)) 3
); // the 'Y' up direction of the section for the beam // ChVector<>(0, cos(rot_rad), sin(rot_rad))
); // the 'Y' up direction of the section for the beam
const auto lastBeamNodesVector = builder.GetLastBeamNodes(); const auto lastBeamNodesVector = builder.GetLastBeamNodes();
// assert(lastBeamNodesVector.size() == 4); // assert(lastBeamNodesVector.size() == 4);
edgeMeshVertsToChronosNodes[vi0] = std::dynamic_pointer_cast<ChNodeFEAxyzrot>( edgeMeshVertsToChronosNodes[vi0] = std::dynamic_pointer_cast<ChNodeFEAxyzrot>(
@ -195,7 +196,8 @@ SimulationResults ChronosIGASimulationModel::executeSimulation(
//parse constrained vertices //parse constrained vertices
// std::cout << node_0->GetCoord().pos[0] << " " << node_0->GetCoord().pos[1] << " " // std::cout << node_0->GetCoord().pos[0] << " " << node_0->GetCoord().pos[1] << " "
// << node_0->GetCoord().pos[2] << std::endl; // << node_0->GetCoord().pos[2] << std::endl;
parseConstrainedVertices(pJob, edgeMeshVertsToChronoNodes, my_system); // parseConstrainedVertices(pJob, edgeMeshVertsToChronoNodes, my_system);
// std::cout << node_0->GetCoord().pos[0] << " " << node_0->GetCoord().pos[1] << " " // std::cout << node_0->GetCoord().pos[0] << " " << node_0->GetCoord().pos[1] << " "
// << node_0->GetCoord().pos[2] << std::endl; // << node_0->GetCoord().pos[2] << std::endl;
// std::dynamic_pointer_cast<std::shared_ptr<ChNodeFEAxyzrot>>(mesh_chronos->GetNode(1)) // std::dynamic_pointer_cast<std::shared_ptr<ChNodeFEAxyzrot>>(mesh_chronos->GetNode(1))
@ -290,7 +292,7 @@ SimulationResults ChronosIGASimulationModel::executeSimulation(
// return simulationResults; // return simulationResults;
} }
void ChronosIGASimulationModel::setStructure(const std::shared_ptr<SimulationMesh> &pMesh) void ChronosIGASimulationModel::setStructure(const std::shared_ptr<SimulationEdgeMesh> &pMesh)
{ {
mesh_chronos = convertToChronosMesh_IGA(pMesh, edgeMeshVertsToChronoNodes); mesh_chronos = convertToChronosMesh_IGA(pMesh, edgeMeshVertsToChronoNodes);
} }

View File

@ -28,9 +28,9 @@ class ChronosIGASimulationModel : public SimulationModel
public: public:
ChronosIGASimulationModel(); ChronosIGASimulationModel();
SimulationResults executeSimulation(const std::shared_ptr<SimulationJob> &pJob) override; SimulationResults executeSimulation(const std::shared_ptr<SimulationJob> &pJob) override;
void setStructure(const std::shared_ptr<SimulationMesh> &pMesh) override; void setStructure(const std::shared_ptr<SimulationEdgeMesh> &pMesh) override;
static std::shared_ptr<chrono::fea::ChMesh> convertToChronosMesh_IGA( static std::shared_ptr<chrono::fea::ChMesh> convertToChronosMesh_IGA(
const std::shared_ptr<SimulationMesh> &pMesh, const std::shared_ptr<SimulationEdgeMesh> &pMesh,
std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes); std::vector<std::shared_ptr<chrono::fea::ChNodeFEAxyzrot>> &edgeMeshVertsToChronosNodes);
inline const static std::string label{"Chronos_linear_IGA"}; inline const static std::string label{"Chronos_linear_IGA"};

View File

@ -4210,7 +4210,7 @@ Eigen::Vector3d DER_leimer::getPerpendicularVector(const Eigen::Vector3d &t) con
} }
DER_leimer::RodConfig *DER_leimer::readRodGrid( DER_leimer::RodConfig *DER_leimer::readRodGrid(
const std::shared_ptr<SimulationMesh> &pMesh, const std::shared_ptr<SimulationEdgeMesh> &pMesh,
const std::vector<int> &numVertsPerRod, const std::vector<int> &numVertsPerRod,
const std::vector<int> & const std::vector<int> &
vInd, //index after the physVerts in the verts vector? Contains "global" indices of all centerline vertices eg at intersections all endpoints vertices get the same "global" index. Do interior centerline vertices get a global index? vInd, //index after the physVerts in the verts vector? Contains "global" indices of all centerline vertices eg at intersections all endpoints vertices get the same "global" index. Do interior centerline vertices get a global index?
@ -4837,7 +4837,7 @@ int DER_leimer::RodConfig::getNumDoF() const
return ndofs; return ndofs;
} }
void DER_leimer::setStructure(const std::shared_ptr<SimulationMesh> &pMesh) void DER_leimer::setStructure(const std::shared_ptr<SimulationEdgeMesh> &pMesh)
{ {
std::cout << "This function is currently not implemented" << std::endl; std::cout << "This function is currently not implemented" << std::endl;
assert(false); assert(false);

View File

@ -140,7 +140,7 @@ public:
int getNumDoF() const; int getNumDoF() const;
}; };
void setStructure(const std::shared_ptr<SimulationMesh> &pMesh) override; void setStructure(const std::shared_ptr<SimulationEdgeMesh> &pMesh) override;
private: private:
bool simulateOneStepGrid(RodConfig *config); bool simulateOneStepGrid(RodConfig *config);
@ -159,7 +159,7 @@ private:
std::vector<Eigen::Triplet<int>> *Hu_index = NULL, std::vector<Eigen::Triplet<int>> *Hu_index = NULL,
std::vector<double> *Hu_value = NULL); std::vector<double> *Hu_value = NULL);
RodConfig *readRodGrid(const std::shared_ptr<SimulationMesh> &pMesh, RodConfig *readRodGrid(const std::shared_ptr<SimulationEdgeMesh> &pMesh,
const std::vector<int> &numVertsPerRod, const std::vector<int> &numVertsPerRod,
const std::vector<int> &vInd, const std::vector<int> &vInd,
const std::vector<Eigen::Vector3d> &verts, const std::vector<Eigen::Vector3d> &verts,