PolyFEM
|
Form of the lagged friction disapative potential and forces. More...
#include <FrictionForm.hpp>
Public Member Functions | |
FrictionForm (const ipc::CollisionMesh &collision_mesh, const std::shared_ptr< time_integrator::ImplicitTimeIntegrator > time_integrator, const double epsv, const double mu, const ipc::BroadPhaseMethod broad_phase_method, const ContactForm &contact_form, const int n_lagging_iters) | |
Construct a new Friction Form object. | |
std::string | name () const override |
void | force_shape_derivative (const Eigen::MatrixXd &prev_solution, const Eigen::MatrixXd &solution, const Eigen::MatrixXd &adjoint, const ipc::FrictionCollisions &friction_constraints_set, Eigen::VectorXd &term) |
void | init_lagging (const Eigen::VectorXd &x) override |
Initialize lagged fields. | |
void | update_lagging (const Eigen::VectorXd &x, const int iter_num) override |
Update lagged fields. | |
void | update_lagging (const Eigen::VectorXd &x) |
Update lagged fields. | |
int | max_lagging_iterations () const override |
Get the maximum number of lagging iteration allowable. | |
bool | uses_lagging () const override |
Does this form require lagging? | |
Eigen::MatrixXd | compute_displaced_surface (const Eigen::VectorXd &x) const |
Compute the displaced positions of the surface nodes. | |
Eigen::MatrixXd | compute_surface_velocities (const Eigen::VectorXd &x) const |
Compute the surface velocities. | |
double | dv_dx () const |
Compute the derivative of the velocities wrt x. | |
double | mu () const |
double | epsv () const |
const ipc::FrictionCollisions & | friction_collision_set () const |
const ipc::FrictionPotential & | friction_potential () const |
Public Member Functions inherited from polyfem::solver::Form | |
virtual | ~Form () |
virtual void | init (const Eigen::VectorXd &x) |
Initialize the form. | |
virtual void | finish () |
virtual double | value (const Eigen::VectorXd &x) const |
Compute the value of the form multiplied with the weigth. | |
Eigen::VectorXd | value_per_element (const Eigen::VectorXd &x) const |
Compute the value of the form multiplied with the weigth. | |
virtual void | first_derivative (const Eigen::VectorXd &x, Eigen::VectorXd &gradv) const |
Compute the first derivative of the value wrt x multiplied with the weigth. | |
void | second_derivative (const Eigen::VectorXd &x, StiffnessMatrix &hessian) const |
Compute the second derivative of the value wrt x multiplied with the weigth. | |
virtual bool | is_step_valid (const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) const |
Determine if a step from solution x0 to solution x1 is allowed. | |
virtual double | max_step_size (const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) const |
Determine the maximum step size allowable between the current and next solution. | |
virtual void | line_search_begin (const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) |
Initialize variables used during the line search. | |
virtual void | line_search_end () |
Clear variables used during the line search. | |
virtual void | post_step (const polysolve::nonlinear::PostStepData &data) |
Update fields after a step in the optimization. | |
virtual void | solution_changed (const Eigen::VectorXd &new_x) |
Update cached fields upon a change in the solution. | |
virtual void | update_quantities (const double t, const Eigen::VectorXd &x) |
Update time-dependent fields. | |
void | set_project_to_psd (bool val) |
Set project to psd. | |
bool | is_project_to_psd () const |
Get if the form's second derivative is projected to psd. | |
void | enable () |
Enable the form. | |
void | disable () |
Disable the form. | |
void | set_enabled (const bool enabled) |
Set if the form is enabled. | |
bool | enabled () const |
Determine if the form is enabled. | |
virtual double | weight () const |
Get the form's multiplicative constant weight. | |
void | set_weight (const double weight) |
Set the form's multiplicative constant weight. | |
virtual bool | is_step_collision_free (const Eigen::VectorXd &x0, const Eigen::VectorXd &x1) const |
Checks if the step is collision free. | |
void | set_output_dir (const std::string &output_dir) |
Protected Member Functions | |
double | value_unweighted (const Eigen::VectorXd &x) const override |
Compute the value of the form. | |
void | first_derivative_unweighted (const Eigen::VectorXd &x, Eigen::VectorXd &gradv) const override |
Compute the first derivative of the value wrt x. | |
void | second_derivative_unweighted (const Eigen::VectorXd &x, StiffnessMatrix &hessian) const override |
Compute the second derivative of the value wrt x. | |
Protected Member Functions inherited from polyfem::solver::Form | |
std::string | resolve_output_path (const std::string &path) const |
virtual Eigen::VectorXd | value_per_element_unweighted (const Eigen::VectorXd &x) const |
Compute the value of the form multiplied per element. | |
Private Attributes | |
const ipc::CollisionMesh & | collision_mesh_ |
Reference to the collision mesh. | |
const std::shared_ptr< time_integrator::ImplicitTimeIntegrator > | time_integrator_ |
Pointer to the time integrator. | |
const double | epsv_ |
Smoothing factor between static and dynamic friction. | |
const double | mu_ |
Global coefficient of friction. | |
const ipc::BroadPhaseMethod | broad_phase_method_ |
Broad-phase method used for distance computation and collision detection. | |
const int | n_lagging_iters_ |
Number of lagging iterations. | |
ipc::FrictionCollisions | friction_collision_set_ |
Lagged friction constraint set. | |
const ContactForm & | contact_form_ |
necessary to have the barrier stiffnes, maybe clean me | |
const ipc::FrictionPotential | friction_potential_ |
Additional Inherited Members | |
Protected Attributes inherited from polyfem::solver::Form | |
bool | project_to_psd_ = false |
If true, the form's second derivative is projected to be positive semidefinite. | |
double | weight_ = 1 |
weight of the form (e.g., AL penalty weight or Δt²) | |
bool | enabled_ = true |
If true, the form is enabled. | |
std::string | output_dir_ |
Form of the lagged friction disapative potential and forces.
Definition at line 18 of file FrictionForm.hpp.
polyfem::solver::FrictionForm::FrictionForm | ( | const ipc::CollisionMesh & | collision_mesh, |
const std::shared_ptr< time_integrator::ImplicitTimeIntegrator > | time_integrator, | ||
const double | epsv, | ||
const double | mu, | ||
const ipc::BroadPhaseMethod | broad_phase_method, | ||
const ContactForm & | contact_form, | ||
const int | n_lagging_iters | ||
) |
Construct a new Friction Form object.
collision_mesh | Reference to the collision mesh |
time_integrator | Pointer to the time integrator |
epsv | Smoothing factor between static and dynamic friction |
mu | Global coefficient of friction |
dhat | Barrier activation distance |
broad_phase_method | Broad-phase method used for distance computation and collision detection |
contact_form | Pointer to contact form; necessary to have the barrier stiffnes, maybe clean me |
n_lagging_iters | Number of lagging iterations |
Definition at line 9 of file FrictionForm.cpp.
References epsv_.
Eigen::MatrixXd polyfem::solver::FrictionForm::compute_displaced_surface | ( | const Eigen::VectorXd & | x | ) | const |
Compute the displaced positions of the surface nodes.
Definition at line 88 of file FrictionForm.cpp.
References polyfem::solver::ContactForm::compute_displaced_surface(), contact_form_, and x.
Referenced by update_lagging().
Eigen::MatrixXd polyfem::solver::FrictionForm::compute_surface_velocities | ( | const Eigen::VectorXd & | x | ) | const |
Compute the surface velocities.
Definition at line 93 of file FrictionForm.cpp.
References collision_mesh_, time_integrator_, polyfem::utils::unflatten(), and x.
Referenced by first_derivative_unweighted(), second_derivative_unweighted(), and value_unweighted().
double polyfem::solver::FrictionForm::dv_dx | ( | ) | const |
Compute the derivative of the velocities wrt x.
Definition at line 100 of file FrictionForm.cpp.
References time_integrator_.
Referenced by second_derivative_unweighted(), and value_unweighted().
|
inline |
Definition at line 87 of file FrictionForm.hpp.
References epsv_.
|
overrideprotectedvirtual |
Compute the first derivative of the value wrt x.
[in] | x | Current solution |
[out] | gradv | Output gradient of the value wrt x |
Implements polyfem::solver::Form.
Definition at line 110 of file FrictionForm.cpp.
References collision_mesh_, compute_surface_velocities(), friction_collision_set_, friction_potential_, and x.
void polyfem::solver::FrictionForm::force_shape_derivative | ( | const Eigen::MatrixXd & | prev_solution, |
const Eigen::MatrixXd & | solution, | ||
const Eigen::MatrixXd & | adjoint, | ||
const ipc::FrictionCollisions & | friction_constraints_set, | ||
Eigen::VectorXd & | term | ||
) |
Definition at line 29 of file FrictionForm.cpp.
References polyfem::solver::ContactForm::barrier_potential(), polyfem::solver::ContactForm::barrier_stiffness(), collision_mesh_, contact_form_, friction_potential_, time_integrator_, and polyfem::utils::unflatten().
|
inline |
Definition at line 88 of file FrictionForm.hpp.
References friction_collision_set_.
|
inline |
Definition at line 89 of file FrictionForm.hpp.
References friction_potential_.
|
inlineoverridevirtual |
Initialize lagged fields.
x | Current solution |
Reimplemented from polyfem::solver::Form.
Definition at line 62 of file FrictionForm.hpp.
References update_lagging(), and x.
|
inlineoverridevirtual |
Get the maximum number of lagging iteration allowable.
Reimplemented from polyfem::solver::Form.
Definition at line 73 of file FrictionForm.hpp.
References n_lagging_iters_.
|
inline |
Definition at line 86 of file FrictionForm.hpp.
References mu_.
|
inlineoverridevirtual |
Implements polyfem::solver::Form.
Definition at line 39 of file FrictionForm.hpp.
|
overrideprotectedvirtual |
Compute the second derivative of the value wrt x.
[in] | x | Current solution |
[out] | hessian | Output Hessian of the value wrt x |
Implements polyfem::solver::Form.
Definition at line 117 of file FrictionForm.cpp.
References collision_mesh_, compute_surface_velocities(), dv_dx(), friction_collision_set_, friction_potential_, POLYFEM_SCOPED_TIMER, polyfem::solver::Form::project_to_psd_, and x.
|
inline |
Update lagged fields.
x | Current solution |
Definition at line 70 of file FrictionForm.hpp.
References update_lagging(), and x.
Referenced by update_lagging().
|
overridevirtual |
Update lagged fields.
x | Current solution |
Reimplemented from polyfem::solver::Form.
Definition at line 127 of file FrictionForm.cpp.
References polyfem::solver::ContactForm::barrier_potential(), polyfem::solver::ContactForm::barrier_stiffness(), broad_phase_method_, collision_mesh_, compute_displaced_surface(), contact_form_, polyfem::solver::ContactForm::dhat(), polyfem::solver::ContactForm::enable_shape_derivatives(), friction_collision_set_, mu_, polyfem::solver::ContactForm::use_convergent_formulation(), and x.
Referenced by init_lagging().
|
inlineoverridevirtual |
Does this form require lagging?
Reimplemented from polyfem::solver::Form.
Definition at line 77 of file FrictionForm.hpp.
|
overrideprotectedvirtual |
Compute the value of the form.
x | Current solution |
Implements polyfem::solver::Form.
Definition at line 105 of file FrictionForm.cpp.
References collision_mesh_, compute_surface_velocities(), dv_dx(), friction_collision_set_, friction_potential_, and x.
|
private |
Broad-phase method used for distance computation and collision detection.
Definition at line 100 of file FrictionForm.hpp.
Referenced by update_lagging().
|
private |
Reference to the collision mesh.
Definition at line 93 of file FrictionForm.hpp.
Referenced by compute_surface_velocities(), first_derivative_unweighted(), force_shape_derivative(), second_derivative_unweighted(), update_lagging(), and value_unweighted().
|
private |
necessary to have the barrier stiffnes, maybe clean me
Definition at line 105 of file FrictionForm.hpp.
Referenced by compute_displaced_surface(), force_shape_derivative(), and update_lagging().
|
private |
Smoothing factor between static and dynamic friction.
Definition at line 98 of file FrictionForm.hpp.
Referenced by epsv(), and FrictionForm().
|
private |
Lagged friction constraint set.
Definition at line 103 of file FrictionForm.hpp.
Referenced by first_derivative_unweighted(), friction_collision_set(), second_derivative_unweighted(), update_lagging(), and value_unweighted().
|
private |
Definition at line 107 of file FrictionForm.hpp.
Referenced by first_derivative_unweighted(), force_shape_derivative(), friction_potential(), second_derivative_unweighted(), and value_unweighted().
|
private |
Global coefficient of friction.
Definition at line 99 of file FrictionForm.hpp.
Referenced by mu(), and update_lagging().
|
private |
Number of lagging iterations.
Definition at line 101 of file FrictionForm.hpp.
Referenced by max_lagging_iterations().
|
private |
Pointer to the time integrator.
Definition at line 96 of file FrictionForm.hpp.
Referenced by compute_surface_velocities(), dv_dx(), and force_shape_derivative().