oomph::ExplicitTimeStepHandler Class Reference

#include <assembly_handler.h>

+ Inheritance diagram for oomph::ExplicitTimeStepHandler:

Public Member Functions

 ExplicitTimeStepHandler ()
 Empty Constructor. More...
 
unsigned ndof (GeneralisedElement *const &elem_pt)
 Return the number of degrees of freedom in the element elem_pt. More...
 
unsigned long eqn_number (GeneralisedElement *const &elem_pt, const unsigned &ieqn_local)
 
void get_residuals (GeneralisedElement *const &elem_pt, Vector< double > &residuals)
 Call the element's residuals. More...
 
void get_jacobian (GeneralisedElement *const &elem_pt, Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Replace get jacobian with the call to get the mass matrix. More...
 
void get_all_vectors_and_matrices (GeneralisedElement *const &elem_pt, Vector< Vector< double >> &vec, Vector< DenseMatrix< double >> &matrix)
 
 ~ExplicitTimeStepHandler ()
 Empty virtual destructor. More...
 
- Public Member Functions inherited from oomph::AssemblyHandler
 AssemblyHandler ()
 Empty constructor. More...
 
virtual void dof_vector (GeneralisedElement *const &elem_pt, const unsigned &t, Vector< double > &dof)
 Return vector of dofs at time level t in the element elem_pt. More...
 
virtual void dof_pt_vector (GeneralisedElement *const &elem_pt, Vector< double * > &dof_pt)
 Return vector of pointers to dofs in the element elem_pt. More...
 
virtual doublelocal_problem_dof (Problem *const &problem_pt, const unsigned &t, const unsigned &i)
 
virtual void get_dresiduals_dparameter (GeneralisedElement *const &elem_pt, double *const &parameter_pt, Vector< double > &dres_dparam)
 
virtual void get_djacobian_dparameter (GeneralisedElement *const &elem_pt, double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
 
virtual void get_hessian_vector_products (GeneralisedElement *const &elem_pt, Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
 
virtual int bifurcation_type () const
 
virtual doublebifurcation_parameter_pt () const
 
virtual void get_eigenfunction (Vector< DoubleVector > &eigenfunction)
 
virtual void get_inner_products (GeneralisedElement *const &elem_pt, Vector< std::pair< unsigned, unsigned >> const &history_index, Vector< double > &inner_product)
 
virtual void get_inner_product_vectors (GeneralisedElement *const &elem_pt, Vector< unsigned > const &history_index, Vector< Vector< double >> &inner_product_vector)
 
virtual ~AssemblyHandler ()
 Empty virtual destructor. More...
 

Detailed Description

A class that is used to define the functions used to assemble and invert the mass matrix when taking an explicit timestep. The idea is simply to replace the jacobian matrix with the mass matrix and then our standard linear solvers will solve the required system

Constructor & Destructor Documentation

◆ ExplicitTimeStepHandler()

oomph::ExplicitTimeStepHandler::ExplicitTimeStepHandler ( )
inline

Empty Constructor.

184 {}

◆ ~ExplicitTimeStepHandler()

oomph::ExplicitTimeStepHandler::~ExplicitTimeStepHandler ( )
inline

Empty virtual destructor.

212 {}

Member Function Documentation

◆ eqn_number()

unsigned long oomph::ExplicitTimeStepHandler::eqn_number ( GeneralisedElement *const &  elem_pt,
const unsigned ieqn_local 
)
virtual

Return the global equation number of the local unknown ieqn_local in elem_pt.

Get the global equation number of the local unknown. Direct call to the function in the element.

Reimplemented from oomph::AssemblyHandler.

250  {
251  return elem_pt->eqn_number(ieqn_local);
252  }

References oomph::GeneralisedElement::eqn_number().

◆ get_all_vectors_and_matrices()

void oomph::ExplicitTimeStepHandler::get_all_vectors_and_matrices ( GeneralisedElement *const &  elem_pt,
Vector< Vector< double >> &  vec,
Vector< DenseMatrix< double >> &  matrix 
)
virtual

Calculate all desired vectors and matrices provided by the element elem_pt.

Calculate all desired vectors and matrices that are required by the problem by calling those of the underlying element.

Reimplemented from oomph::AssemblyHandler.

282  {
283 #ifdef PARANOID
284  // Check dimension
285  if (matrix.size() != 1)
286  {
287  throw OomphLibError("ExplicitTimeSteps should return one matrix",
290  }
291 #endif
292  // Get the residuals and <mass matrix
293  elem_pt->get_mass_matrix(vec[0], matrix[0]);
294  }
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::ColMajor >, 0, Eigen::OuterStride<> > matrix(T *data, int rows, int cols, int stride)
Definition: common.h:85
#define OOMPH_EXCEPTION_LOCATION
Definition: oomph_definitions.h:61
#define OOMPH_CURRENT_FUNCTION
Definition: oomph_definitions.h:86

References oomph::GeneralisedElement::get_mass_matrix(), matrix(), OOMPH_CURRENT_FUNCTION, and OOMPH_EXCEPTION_LOCATION.

◆ get_jacobian()

void oomph::ExplicitTimeStepHandler::get_jacobian ( GeneralisedElement *const &  elem_pt,
Vector< double > &  residuals,
DenseMatrix< double > &  jacobian 
)
virtual

Replace get jacobian with the call to get the mass matrix.

Calculate the elemental Jacobian matrix "d equation / d variable" for elem_pt. Again deliberately broken in the eigenproblem

Reimplemented from oomph::AssemblyHandler.

269  {
270  elem_pt->get_mass_matrix(residuals, jacobian);
271  }

References oomph::GeneralisedElement::get_mass_matrix().

◆ get_residuals()

void oomph::ExplicitTimeStepHandler::get_residuals ( GeneralisedElement *const &  elem_pt,
Vector< double > &  residuals 
)
virtual

Call the element's residuals.

Return the contribution to the residuals of the element elem_pt This is deliberately broken in our eigenproblem

Reimplemented from oomph::AssemblyHandler.

259  {
260  elem_pt->get_residuals(residuals);
261  }

References oomph::GeneralisedElement::get_residuals().

◆ ndof()

unsigned oomph::ExplicitTimeStepHandler::ndof ( GeneralisedElement *const &  elem_pt)
virtual

Return the number of degrees of freedom in the element elem_pt.

Get the number of elemental degrees of freedom. Direct call to the function in the element.

Reimplemented from oomph::AssemblyHandler.

240  {
241  return elem_pt->ndof();
242  }

References oomph::GeneralisedElement::ndof().


The documentation for this class was generated from the following files: