Loading [MathJax]/extensions/tex2jax.js
PolyFEM
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DiffCache.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <polyfem/Common.hpp>
5#include <ipc/ipc.hpp>
6#include <ipc/collisions/normal/normal_collisions.hpp>
7#include <ipc/collisions/tangential/tangential_collisions.hpp>
8
9namespace polyfem::solver
10{
11 enum class CacheLevel
12 {
13 None,
16 };
17
19 {
20 public:
21 void init(const int dimension, const int ndof, const int n_time_steps = 0)
22 {
23 cur_size_ = 0;
24 n_time_steps_ = n_time_steps;
25
26 u_.setZero(ndof, n_time_steps + 1);
27 disp_grad_.assign(n_time_steps + 1, Eigen::MatrixXd::Zero(dimension,dimension));
28 if (n_time_steps_ > 0)
29 {
30 bdf_order_.setZero(n_time_steps + 1);
31 v_.setZero(ndof, n_time_steps + 1);
32 acc_.setZero(ndof, n_time_steps + 1);
33 // gradu_h_prev_.resize(n_time_steps + 1);
34 }
35 gradu_h_.resize(n_time_steps + 1);
36 collision_set_.resize(n_time_steps + 1);
37 friction_collision_set_.resize(n_time_steps + 1);
38 normal_adhesion_collision_set_.resize(n_time_steps + 1);
39 tangential_adhesion_collision_set_.resize(n_time_steps + 1);
40 }
41
43 const Eigen::MatrixXd &u,
45 const ipc::NormalCollisions &contact_set,
46 const ipc::TangentialCollisions &friction_constraint_set,
47 const ipc::NormalCollisions &normal_adhesion_set,
48 const ipc::TangentialCollisions &tangential_adhesion_set,
49 const Eigen::MatrixXd &disp_grad)
50 {
51 u_ = u;
52
53 gradu_h_[0] = gradu_h;
54 collision_set_[0] = contact_set;
55 friction_collision_set_[0] = friction_constraint_set;
56 normal_adhesion_collision_set_[0] = normal_adhesion_set;
57 tangential_adhesion_collision_set_[0] = tangential_adhesion_set;
59
60 cur_size_ = 1;
61 }
62
64 const int cur_step,
65 const int cur_bdf_order,
66 const Eigen::MatrixXd &u,
67 const Eigen::MatrixXd &v,
68 const Eigen::MatrixXd &acc,
70 // const StiffnessMatrix &gradu_h_prev,
71 const ipc::NormalCollisions &collision_set,
72 const ipc::TangentialCollisions &friction_collision_set)
73 {
74 bdf_order_(cur_step) = cur_bdf_order;
75
76 u_.col(cur_step) = u;
77 v_.col(cur_step) = v;
78 acc_.col(cur_step) = acc;
79
80 gradu_h_[cur_step] = gradu_h;
81 // gradu_h_prev_[cur_step] = gradu_h_prev;
82
83 collision_set_[cur_step] = collision_set;
85
86 cur_size_++;
87 }
88
90 const int cur_step,
91 const Eigen::MatrixXd &u,
93 const ipc::NormalCollisions &contact_set,
94 const ipc::NormalCollisions &normal_adhesion_set,
95 const Eigen::MatrixXd &disp_grad)
96 {
97 u_.col(cur_step) = u;
98 gradu_h_[cur_step] = gradu_h;
99 collision_set_[cur_step] = contact_set;
100 normal_adhesion_collision_set_[cur_step] = normal_adhesion_set;
101 disp_grad_[cur_step] = disp_grad;
102
103 cur_size_++;
104 }
105
106 void cache_adjoints(const Eigen::MatrixXd &adjoint_mat) { adjoint_mat_ = adjoint_mat; }
107 const Eigen::MatrixXd &adjoint_mat() const { return adjoint_mat_; }
108
109 inline int size() const { return cur_size_; }
110 inline int bdf_order(int step) const
111 {
112 assert(step < size());
113 if (step < 0)
114 step += bdf_order_.size();
115 return bdf_order_(step);
116 }
117
118 Eigen::MatrixXd disp_grad(int step = 0) const { assert(step < size()); if (step < 0) step += disp_grad_.size(); return disp_grad_[step]; }
119
120 Eigen::VectorXd u(int step) const
121 {
122 assert(step < size());
123 if (step < 0)
124 step += u_.cols();
125 return u_.col(step);
126 }
127 Eigen::VectorXd v(int step) const
128 {
129 assert(step < size());
130 if (step < 0)
131 step += v_.cols();
132 return v_.col(step);
133 }
134 Eigen::VectorXd acc(int step) const
135 {
136 assert(step < size());
137 if (step < 0)
138 step += acc_.cols();
139 return acc_.col(step);
140 }
141
142 const StiffnessMatrix &gradu_h(int step) const
143 {
144 assert(step < size());
145 if (step < 0)
146 step += gradu_h_.size();
147 return gradu_h_[step];
148 }
149 // const StiffnessMatrix &gradu_h_prev(const int step) const { assert(step < size()); return gradu_h_prev_[step]; }
150
151 const ipc::NormalCollisions &collision_set(int step) const
152 {
153 assert(step < size());
154 if (step < 0)
155 step += collision_set_.size();
156 return collision_set_[step];
157 }
158 const ipc::TangentialCollisions &friction_collision_set(int step) const
159 {
160 assert(step < size());
161 if (step < 0)
162 step += friction_collision_set_.size();
163 return friction_collision_set_[step];
164 }
165 const ipc::NormalCollisions &normal_adhesion_collision_set(int step) const
166 {
167 assert(step < size());
168 if (step < 0)
169 step += normal_adhesion_collision_set_.size();
171 }
172 const ipc::TangentialCollisions &tangential_adhesion_collision_set(int step) const
173 {
174 assert(step < size());
175 if (step < 0)
178 }
179
180 private:
182 int cur_size_ = 0;
183
184 std::vector<Eigen::MatrixXd> disp_grad_; // macro linear displacement in homogenization
185 Eigen::MatrixXd u_; // PDE solution
186 Eigen::MatrixXd v_; // velocity in transient elastic simulations
187 Eigen::MatrixXd acc_; // acceleration in transient elastic simulations
188
189 Eigen::VectorXi bdf_order_; // BDF orders used at each time step in forward simulation
190
191 std::vector<StiffnessMatrix> gradu_h_; // gradient of force at time T wrt. u at time T
192 // std::vector<StiffnessMatrix> gradu_h_prev_; // gradient of force at time T wrt. u at time (T-1) in transient simulations
193
194 std::vector<ipc::NormalCollisions> collision_set_;
195 std::vector<ipc::TangentialCollisions> friction_collision_set_;
196
197 std::vector<ipc::NormalCollisions> normal_adhesion_collision_set_;
198 std::vector<ipc::TangentialCollisions> tangential_adhesion_collision_set_;
199
200 Eigen::MatrixXd adjoint_mat_;
201 };
202} // namespace polyfem::solver
const ipc::NormalCollisions & collision_set(int step) const
const StiffnessMatrix & gradu_h(int step) const
int bdf_order(int step) const
void cache_quantities_transient(const int cur_step, const int cur_bdf_order, const Eigen::MatrixXd &u, const Eigen::MatrixXd &v, const Eigen::MatrixXd &acc, const StiffnessMatrix &gradu_h, const ipc::NormalCollisions &collision_set, const ipc::TangentialCollisions &friction_collision_set)
Definition DiffCache.hpp:63
std::vector< ipc::TangentialCollisions > tangential_adhesion_collision_set_
Eigen::MatrixXd adjoint_mat_
const ipc::TangentialCollisions & friction_collision_set(int step) const
Eigen::MatrixXd disp_grad(int step=0) const
const ipc::TangentialCollisions & tangential_adhesion_collision_set(int step) const
void init(const int dimension, const int ndof, const int n_time_steps=0)
Definition DiffCache.hpp:21
std::vector< Eigen::MatrixXd > disp_grad_
std::vector< ipc::NormalCollisions > normal_adhesion_collision_set_
void cache_quantities_quasistatic(const int cur_step, const Eigen::MatrixXd &u, const StiffnessMatrix &gradu_h, const ipc::NormalCollisions &contact_set, const ipc::NormalCollisions &normal_adhesion_set, const Eigen::MatrixXd &disp_grad)
Definition DiffCache.hpp:89
const Eigen::MatrixXd & adjoint_mat() const
const ipc::NormalCollisions & normal_adhesion_collision_set(int step) const
std::vector< StiffnessMatrix > gradu_h_
Eigen::VectorXd acc(int step) const
void cache_adjoints(const Eigen::MatrixXd &adjoint_mat)
std::vector< ipc::TangentialCollisions > friction_collision_set_
void cache_quantities_static(const Eigen::MatrixXd &u, const StiffnessMatrix &gradu_h, const ipc::NormalCollisions &contact_set, const ipc::TangentialCollisions &friction_constraint_set, const ipc::NormalCollisions &normal_adhesion_set, const ipc::TangentialCollisions &tangential_adhesion_set, const Eigen::MatrixXd &disp_grad)
Definition DiffCache.hpp:42
std::vector< ipc::NormalCollisions > collision_set_
Eigen::VectorXi bdf_order_
Eigen::VectorXd v(int step) const
Eigen::VectorXd u(int step) const
Eigen::SparseMatrix< double, Eigen::ColMajor > StiffnessMatrix
Definition Types.hpp:22