1 // =============================================================================
2 // PROJECT CHRONO - http://projectchrono.org
3 //
4 // Copyright (c) 2014 projectchrono.org
5 // All rights reserved.
6 //
7 // Use of this source code is governed by a BSD-style license that can be found
8 // in the LICENSE file at the top level of the distribution and at
9 // http://projectchrono.org/license-chrono.txt.
10 //
11 // =============================================================================
12 // Authors: Alessandro Tasora, Radu Serban
13 // =============================================================================
14 //
15 // FEA for 3D beams of 'cable' type (ANCF gradient-deficient beams)
16 //
17 // =============================================================================
18
19 #include "chrono/physics/ChSystemSMC.h"
20 #include "chrono/solver/ChDirectSolverLS.h"
21 #include "chrono/solver/ChIterativeSolverLS.h"
22 #include "chrono/timestepper/ChTimestepper.h"
23 #include "chrono_irrlicht/ChIrrApp.h"
24
25 #include "FEAcables.h"
26
27 using namespace chrono;
28 using namespace chrono::fea;
29 using namespace chrono::irrlicht;
30
31 using namespace irr;
32
33 // Select solver type (SPARSE_QR, SPARSE_LU, or MINRES).
34 ChSolver::Type solver_type = ChSolver::Type::SPARSE_QR;
35
main(int argc,char * argv[])36 int main(int argc, char* argv[]) {
37 GetLog() << "Copyright (c) 2017 projectchrono.org\nChrono version: " << CHRONO_VERSION << "\n\n";
38
39 // Create a Chrono::Engine physical system
40 ChSystemSMC my_system;
41
42 // Create the Irrlicht visualization (open the Irrlicht device,
43 // bind a simple user interface, etc. etc.)
44 ChIrrApp application(&my_system, L"Cables FEM", core::dimension2d<u32>(800, 600));
45
46 // Easy shortcuts to add camera, lights, logo and sky in Irrlicht scene:
47 application.AddTypicalLogo();
48 application.AddTypicalSky();
49 application.AddTypicalLights();
50 application.AddTypicalCamera(core::vector3df(0.f, 0.6f, -1.f));
51
52 // Create a mesh, that is a container for groups of elements and
53 // their referenced nodes.
54 auto my_mesh = chrono_types::make_shared<ChMesh>();
55
56 // Create one of the available models (defined in FEAcables.h)
57 ////auto model = Model1(my_system, my_mesh);
58 ////auto model = Model2(my_system, my_mesh);
59 auto model = Model3(my_system, my_mesh);
60
61 // Remember to add the mesh to the system!
62 my_system.Add(my_mesh);
63
64 // ==Asset== attach a visualization of the FEM mesh.
65 // This will automatically update a triangle mesh (a ChTriangleMeshShape
66 // asset that is internally managed) by setting proper
67 // coordinates and vertex colors as in the FEM elements.
68 // Such triangle mesh can be rendered by Irrlicht or POVray or whatever
69 // postprocessor that can handle a colored ChTriangleMeshShape).
70 // Do not forget AddAsset() at the end!
71
72 auto mvisualizebeamA = chrono_types::make_shared<ChVisualizationFEAmesh>(*(my_mesh.get()));
73 mvisualizebeamA->SetFEMdataType(ChVisualizationFEAmesh::E_PLOT_ELEM_BEAM_MZ);
74 mvisualizebeamA->SetColorscaleMinMax(-0.4, 0.4);
75 mvisualizebeamA->SetSmoothFaces(true);
76 mvisualizebeamA->SetWireframe(false);
77 my_mesh->AddAsset(mvisualizebeamA);
78
79 auto mvisualizebeamC = chrono_types::make_shared<ChVisualizationFEAmesh>(*(my_mesh.get()));
80 mvisualizebeamC->SetFEMglyphType(ChVisualizationFEAmesh::E_GLYPH_NODE_DOT_POS); // E_GLYPH_NODE_CSYS
81 mvisualizebeamC->SetFEMdataType(ChVisualizationFEAmesh::E_PLOT_NONE);
82 mvisualizebeamC->SetSymbolsThickness(0.006);
83 mvisualizebeamC->SetSymbolsScale(0.01);
84 mvisualizebeamC->SetZbufferHide(false);
85 my_mesh->AddAsset(mvisualizebeamC);
86
87 // ==IMPORTANT!== Use this function for adding a ChIrrNodeAsset to all items
88 // in the system. These ChIrrNodeAsset assets are 'proxies' to the Irrlicht meshes.
89 // If you need a finer control on which item really needs a visualization proxy in
90 // Irrlicht, just use application.AssetBind(myitem); on a per-item basis.
91 application.AssetBindAll();
92
93 // ==IMPORTANT!== Use this function for 'converting' into Irrlicht meshes the assets
94 // that you added to the bodies into 3D shapes, they can be visualized by Irrlicht!
95 application.AssetUpdateAll();
96
97 // Set solver and solver settings
98 switch (solver_type) {
99 case ChSolver::Type::SPARSE_QR: {
100 std::cout << "Using SparseQR solver" << std::endl;
101 auto solver = chrono_types::make_shared<ChSolverSparseQR>();
102 my_system.SetSolver(solver);
103 solver->UseSparsityPatternLearner(true);
104 solver->LockSparsityPattern(true);
105 solver->SetVerbose(false);
106 break;
107 }
108 case ChSolver::Type::SPARSE_LU: {
109 std::cout << "Using SparseLU solver" << std::endl;
110 auto solver = chrono_types::make_shared<ChSolverSparseLU>();
111 my_system.SetSolver(solver);
112 solver->UseSparsityPatternLearner(true);
113 solver->LockSparsityPattern(true);
114 solver->SetVerbose(false);
115 break;
116 }
117 case ChSolver::Type::MINRES: {
118 std::cout << "Using MINRES solver" << std::endl;
119 auto solver = chrono_types::make_shared<ChSolverMINRES>();
120 my_system.SetSolver(solver);
121 solver->SetMaxIterations(200);
122 solver->SetTolerance(1e-10);
123 solver->EnableDiagonalPreconditioner(true);
124 solver->EnableWarmStart(true); // IMPORTANT for convergence when using EULER_IMPLICIT_LINEARIZED
125 solver->SetVerbose(false);
126 break;
127 }
128 default: {
129 std::cout << "Solver type not supported." << std::endl;
130 break;
131 }
132 }
133
134 my_system.SetSolverForceTolerance(1e-13);
135
136 // Set integrator
137 my_system.SetTimestepperType(ChTimestepper::Type::EULER_IMPLICIT_LINEARIZED);
138
139 // SIMULATION LOOP
140 application.SetTimestep(0.01);
141
142 while (application.GetDevice()->run()) {
143 application.BeginScene();
144 application.DrawAll();
145 application.DoStep();
146 application.EndScene();
147 ////model.PrintBodyPositions();
148 }
149
150 return 0;
151 }
152