Loading [MathJax]/extensions/tex2jax.js
PolyFEM
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SolveData.hpp
Go to the documentation of this file.
1#pragma once
2
10
11#include <ipc/collision_mesh.hpp>
12#include <ipc/broad_phase/broad_phase.hpp>
13
14#include <Eigen/Core>
15
16#include <memory>
17#include <string>
18#include <unordered_map>
19
21
23{
24 class ImplicitTimeIntegrator;
25} // namespace polyfem::time_integrator
26
27namespace polyfem::utils
28{
29 class PeriodicBoundary;
30}
31
32namespace polyfem::assembler
33{
34 class ViscousDamping;
35 class MacroStrainValue;
36 class PressureAssembler;
37} // namespace polyfem::assembler
38
39namespace polyfem::solver
40{
41 class NLProblem;
42 class Form;
43 class ContactForm;
44 class PeriodicContactForm;
45 class MacroStrainALForm;
46 class FrictionForm;
47 class BodyForm;
48 class AugmentedLagrangianForm;
49 class MacroStrainLagrangianForm;
50 class MacroStrainALForm;
51 class InertiaForm;
52 class ElasticForm;
53 class PressureForm;
54 class NormalAdhesionForm;
55 class TangentialAdhesionForm;
56
59 {
60 public:
63 std::vector<std::shared_ptr<Form>> init_forms(
64 // General
65 const Units &units,
66 const int dim,
67 const double t,
68 const Eigen::VectorXi &in_node_to_node,
69
70 // Elastic form
71 const int n_bases,
72 std::vector<basis::ElementBases> &bases,
73 const std::vector<basis::ElementBases> &geom_bases,
74 const assembler::Assembler &assembler,
75 assembler::AssemblyValsCache &ass_vals_cache,
76 const assembler::AssemblyValsCache &mass_ass_vals_cache,
77 const double jacobian_threshold,
78 const solver::ElementInversionCheck check_inversion,
79
80 // Body form
81 const int n_pressure_bases,
82 const std::vector<int> &boundary_nodes,
83 const std::vector<mesh::LocalBoundary> &local_boundary,
84 const std::vector<mesh::LocalBoundary> &local_neumann_boundary,
85 const int n_boundary_samples,
86 const Eigen::MatrixXd &rhs,
87 const Eigen::MatrixXd &sol,
88 const assembler::Density &density,
89
90 // Pressure form
91 const std::vector<mesh::LocalBoundary> &local_pressure_boundary,
92 const std::unordered_map<int, std::vector<mesh::LocalBoundary>> &local_pressure_cavity,
93 const std::shared_ptr<assembler::PressureAssembler> pressure_assembler,
94
95 // Inertia form
96 const bool ignore_inertia,
97 const StiffnessMatrix &mass,
98 const std::shared_ptr<assembler::ViscousDamping> damping_assembler,
99
100 // Lagged regularization form
101 const double lagged_regularization_weight,
102 const int lagged_regularization_iterations,
103
104 // Constraint forms
105 const size_t obstacle_ndof,
106 const std::vector<std::string> &hard_constraint_files,
107 const std::vector<json> &soft_constraint_files,
108
109 // Contact form
110 const bool contact_enabled,
111 const ipc::CollisionMesh &collision_mesh,
112 const double dhat,
113 const double avg_mass,
114 const bool use_area_weighting,
115 const bool use_improved_max_operator,
116 const bool use_physical_barrier,
117 const json &barrier_stiffness,
118 const ipc::BroadPhaseMethod broad_phase,
119 const double ccd_tolerance,
120 const long ccd_max_iterations,
121 const bool enable_shape_derivatives,
122
123 // Normal Adhesion Form
124 const bool adhesion_enabled,
125 const double dhat_p,
126 const double dhat_a,
127 const double Y,
128
129 // Tangential Adhesion Form
130 const double tangential_adhesion_coefficient,
131 const double epsa,
132 const int tangential_adhesion_iterations,
133
134 // Homogenization
135 const assembler::MacroStrainValue &macro_strain_constraint,
136
137 // Periodic contact
138 const bool periodic_contact,
139 const Eigen::VectorXi &tiled_to_single,
140 const std::shared_ptr<utils::PeriodicBoundary> &periodic_bc,
141
142 // Friction form
143 const double friction_coefficient,
144 const double epsv,
145 const int friction_iterations,
146
147 // Rayleigh damping form
148 const json &rayleigh_damping);
149
152 void update_barrier_stiffness(const Eigen::VectorXd &x);
153
155 void update_dt();
156
157 std::vector<std::pair<std::string, std::shared_ptr<solver::Form>>> named_forms() const;
158
159 public:
160 std::shared_ptr<assembler::RhsAssembler> rhs_assembler;
161 std::shared_ptr<assembler::PressureAssembler> pressure_assembler;
162 std::shared_ptr<solver::NLProblem> nl_problem;
163
164 std::vector<std::shared_ptr<solver::AugmentedLagrangianForm>> al_form;
165 std::shared_ptr<solver::MacroStrainLagrangianForm> strain_al_lagr_form;
166 std::shared_ptr<solver::BodyForm> body_form;
167 std::shared_ptr<solver::ContactForm> contact_form;
168 std::shared_ptr<solver::ElasticForm> damping_form;
169 std::shared_ptr<solver::ElasticForm> elastic_form;
170 std::shared_ptr<solver::FrictionForm> friction_form;
171 std::shared_ptr<solver::InertiaForm> inertia_form;
172 std::shared_ptr<solver::PressureForm> pressure_form;
173 std::shared_ptr<solver::NormalAdhesionForm> normal_adhesion_form;
174 std::shared_ptr<solver::TangentialAdhesionForm> tangential_adhesion_form;
175
176 std::shared_ptr<solver::PeriodicContactForm> periodic_contact_form;
177
178 std::shared_ptr<time_integrator::ImplicitTimeIntegrator> time_integrator;
179 };
180} // namespace polyfem::solver
int x
Caches basis evaluation and geometric mapping at every element.
class to store time stepping data
Definition SolveData.hpp:59
std::shared_ptr< solver::FrictionForm > friction_form
std::shared_ptr< solver::InertiaForm > inertia_form
std::shared_ptr< solver::PeriodicContactForm > periodic_contact_form
std::vector< std::shared_ptr< Form > > init_forms(const Units &units, const int dim, const double t, const Eigen::VectorXi &in_node_to_node, const int n_bases, std::vector< basis::ElementBases > &bases, const std::vector< basis::ElementBases > &geom_bases, const assembler::Assembler &assembler, assembler::AssemblyValsCache &ass_vals_cache, const assembler::AssemblyValsCache &mass_ass_vals_cache, const double jacobian_threshold, const solver::ElementInversionCheck check_inversion, const int n_pressure_bases, const std::vector< int > &boundary_nodes, const std::vector< mesh::LocalBoundary > &local_boundary, const std::vector< mesh::LocalBoundary > &local_neumann_boundary, const int n_boundary_samples, const Eigen::MatrixXd &rhs, const Eigen::MatrixXd &sol, const assembler::Density &density, const std::vector< mesh::LocalBoundary > &local_pressure_boundary, const std::unordered_map< int, std::vector< mesh::LocalBoundary > > &local_pressure_cavity, const std::shared_ptr< assembler::PressureAssembler > pressure_assembler, const bool ignore_inertia, const StiffnessMatrix &mass, const std::shared_ptr< assembler::ViscousDamping > damping_assembler, const double lagged_regularization_weight, const int lagged_regularization_iterations, const size_t obstacle_ndof, const std::vector< std::string > &hard_constraint_files, const std::vector< json > &soft_constraint_files, const bool contact_enabled, const ipc::CollisionMesh &collision_mesh, const double dhat, const double avg_mass, const bool use_area_weighting, const bool use_improved_max_operator, const bool use_physical_barrier, const json &barrier_stiffness, const ipc::BroadPhaseMethod broad_phase, const double ccd_tolerance, const long ccd_max_iterations, const bool enable_shape_derivatives, const bool adhesion_enabled, const double dhat_p, const double dhat_a, const double Y, const double tangential_adhesion_coefficient, const double epsa, const int tangential_adhesion_iterations, const assembler::MacroStrainValue &macro_strain_constraint, const bool periodic_contact, const Eigen::VectorXi &tiled_to_single, const std::shared_ptr< utils::PeriodicBoundary > &periodic_bc, const double friction_coefficient, const double epsv, const int friction_iterations, const json &rayleigh_damping)
Initialize the forms and return a vector of pointers to them.
Definition SolveData.cpp:34
void update_dt()
updates the dt inside the different forms
std::shared_ptr< solver::PressureForm > pressure_form
std::shared_ptr< solver::BodyForm > body_form
std::shared_ptr< solver::NLProblem > nl_problem
std::shared_ptr< solver::MacroStrainLagrangianForm > strain_al_lagr_form
std::shared_ptr< solver::NormalAdhesionForm > normal_adhesion_form
std::shared_ptr< solver::ContactForm > contact_form
std::vector< std::pair< std::string, std::shared_ptr< solver::Form > > > named_forms() const
std::shared_ptr< solver::ElasticForm > damping_form
std::shared_ptr< solver::ElasticForm > elastic_form
std::vector< std::shared_ptr< solver::AugmentedLagrangianForm > > al_form
std::shared_ptr< time_integrator::ImplicitTimeIntegrator > time_integrator
void update_barrier_stiffness(const Eigen::VectorXd &x)
update the barrier stiffness for the forms
std::shared_ptr< solver::TangentialAdhesionForm > tangential_adhesion_form
std::shared_ptr< assembler::PressureAssembler > pressure_assembler
std::shared_ptr< assembler::RhsAssembler > rhs_assembler
Used for test only.
nlohmann::json json
Definition Common.hpp:9
Eigen::SparseMatrix< double, Eigen::ColMajor > StiffnessMatrix
Definition Types.hpp:22