Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F63673834
parabolic_driver.hpp
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Tue, May 21, 18:43
Size
3 KB
Mime Type
text/x-c++
Expires
Thu, May 23, 18:43 (1 d, 23 h)
Engine
blob
Format
Raw Data
Handle
17805080
Attached To
rSPECMICP SpecMiCP / ReactMiCP
parabolic_driver.hpp
View Options
#ifndef SPECMICP_DFPMSOLVER_PARABOLICDRIVER_HPP
#define SPECMICP_DFPMSOLVER_PARABOLICDRIVER_HPP
#include "common.hpp"
#include <Eigen/Sparse>
#include "driver.hpp"
#include "parabolic_structs.hpp"
namespace specmicp {
namespace dfpmsolver {
//! \brief The parabolic driver
//!
//! Parabolic driver for finite element programs
template <class Program>
class ParabolicDriver: public Driver<Program, ParabolicDriverOptions, ParabolicDriverPerformance>
{
public:
using base=Driver<Program, ParabolicDriverOptions, ParabolicDriverPerformance>;
using base::program;
using base::get_neq;
using base::get_options;
using base::get_perfs;
ParabolicDriver(Program& the_program):
Driver<Program, ParabolicDriverOptions, ParabolicDriverPerformance>(the_program)
{}
//! \brief Solve a timestep of length dt
ParabolicDriverReturnCode solve_timestep(scalar_t dt, Vector& displacement);
//! \brief Restart the current timestep
ParabolicDriverReturnCode restart_timestep(Vector& displacement);
//! \brief Check if the solution has been found
ParabolicDriverReturnCode check_convergence(scalar_t norm_update);
void compute_residuals(const Vector& displacement,
Vector& residual)
{
compute_residuals(displacement, m_velocity, residual);
}
//! \brief Compute the residuals
void compute_residuals(const Vector& displacement,
const Vector& velocity,
Vector& residual)
{
program().compute_residuals(displacement, velocity, residual);
get_perfs().nb_call_residuals += 1;
}
//! \brief Compute the jacobian
void compute_jacobian(Vector& displacement,
Vector& velocity,
Eigen::SparseMatrix<scalar_t>& jacobian
)
{
program().compute_jacobian(displacement, velocity, jacobian, get_options().alpha*m_current_dt);
get_perfs().nb_call_jacobian += 1;
}
//! \brief Return the norm of the current residuals
double norm_residuals() {return m_residuals.norm();}
Eigen::VectorXd& get_velocity() {return m_velocity;}
private:
void update_variable(
const Vector& update,
scalar_t lambda,
Vector& displacement
);
//! \brief Set the predictor
void set_predictor(Vector& displacement);
//! \brief solve the linear system
int solve_linear_system(Eigen::VectorXd& update);
//! \brief perform the linesearch
int linesearch(Vector &update,
Vector &displacements);
double compute_residuals_linesearch(
Vector& update,
scalar_t lambda,
Vector& displacement
);
//! \brief Initialize the computation
void initialize_timestep(scalar_t dt, Eigen::VectorXd& displacement);
double m_norm_0;
double m_current_dt;
Eigen::VectorXd m_gradient;
Eigen::VectorXd m_velocity;
Eigen::VectorXd m_residuals;
Eigen::VectorXd m_predictor;
Eigen::SparseMatrix<double> m_jacobian;
};
} // end namespace dfpmsolver
} // end namespace specmicp
// implementation
// ==============
#include "parabolic_driver.inl"
#endif // SPECMICP_DFPMSOLVER_PARABOLICDRIVER_HPP
Event Timeline
Log In to Comment