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