Public Member Functions | Private Types | Private Member Functions | Private Attributes | Friends | List of all members
oomph::PitchForkHandler Class Reference

A class that is used to assemble the augmented system that defines a pitchfork (symmetry-breaking) bifurcation. The "standard" problem must be a function of a global parameter $ \lambda $ and a solution is $R(u,\lambda) = 0$, where $u$ are the unknowns in the problem. A pitchfork bifurcation may be specified by the augmented system of size $2N+2$. More...

#include <assembly_handler.h>

+ Inheritance diagram for oomph::PitchForkHandler:

Public Member Functions

 PitchForkHandler (Problem *const &problem_pt, AssemblyHandler *const &assembly_handler_pt, double *const &parameter_pt, const DoubleVector &symmetry_vector)
 Constructor, initialise the systems. More...
 
 ~PitchForkHandler ()
 Destructor, return the problem to its original state, before the augmented system was added. More...
 
unsigned ndof (GeneralisedElement *const &elem_pt)
 Get the number of elemental degrees of freedom. More...
 
unsigned long eqn_number (GeneralisedElement *const &elem_pt, const unsigned &ieqn_local)
 Get the global equation number of the local unknown. More...
 
void get_residuals (GeneralisedElement *const &elem_pt, Vector< double > &residuals)
 Get the residuals. More...
 
void get_jacobian (GeneralisedElement *const &elem_pt, Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Calculate the elemental Jacobian matrix "d equation / d variable". More...
 
void get_dresiduals_dparameter (GeneralisedElement *const &elem_pt, double *const &parameter_pt, Vector< double > &dres_dparam)
 Overload the derivatives of the residuals with respect to a parameter to apply to the augmented system. More...
 
void get_djacobian_dparameter (GeneralisedElement *const &elem_pt, double *const &parameter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
 Overload the derivative of the residuals and jacobian with respect to a parameter so that it breaks. More...
 
void get_hessian_vector_products (GeneralisedElement *const &elem_pt, Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
 Overload the hessian vector product function so that it breaks. More...
 
int bifurcation_type () const
 Indicate that we are tracking a pitchfork bifurcation by returning 2. More...
 
double * bifurcation_parameter_pt () const
 Return a pointer to the bifurcation parameter in bifurcation tracking problems. More...
 
void get_eigenfunction (Vector< DoubleVector > &eigenfunction)
 Return the eigenfunction(s) associated with the bifurcation that has been detected in bifurcation tracking problems. More...
 
void synchronise ()
 Function that is used to perform any synchronisation required during the solution. More...
 
void solve_augmented_block_system ()
 Set to solve the augmented block system. More...
 
void solve_block_system ()
 Set to solve the block system. More...
 
void solve_full_system ()
 Solve non-block system. 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 double & local_problem_dof (Problem *const &problem_pt, const unsigned &t, const unsigned &i)
 Return the t-th level of storage associated with the i-th (local) dof stored in the problem. More...
 
virtual void get_all_vectors_and_matrices (GeneralisedElement *const &elem_pt, Vector< Vector< double >> &vec, Vector< DenseMatrix< double >> &matrix)
 Calculate all desired vectors and matrices provided by the element elem_pt. More...
 
virtual void get_inner_products (GeneralisedElement *const &elem_pt, Vector< std::pair< unsigned, unsigned >> const &history_index, Vector< double > &inner_product)
 Compute the inner products of the given vector of pairs of history values over the element. More...
 
virtual void get_inner_product_vectors (GeneralisedElement *const &elem_pt, Vector< unsigned > const &history_index, Vector< Vector< double >> &inner_product_vector)
 Compute the vectors that when taken as a dot product with other history values give the inner product over the element. More...
 
virtual ~AssemblyHandler ()
 Empty virtual destructor. More...
 

Private Types

enum  { Full_augmented , Block_J , Block_augmented_J }
 A little private enum to determine whether we are solving the block system or not. More...
 

Private Member Functions

unsigned global_eqn_number (const unsigned &i)
 Function that is used to return map the global equations using the simplistic numbering scheme into the actual distributed scheme. More...
 

Private Attributes

unsigned Solve_which_system
 Integer flag to indicate which system should be assembled. There are three possibilities. The full augmented system (0), the non-augmented jacobian system (1), a system in which the jacobian is augmented by 1 row and column to ensure that it is non-singular (2). See the enum above. More...
 
ProblemProblem_pt
 Pointer to the problem. More...
 
AssemblyHandlerAssembly_handler_pt
 Pointer to the underlying (original) assembly handler. More...
 
unsigned Ndof
 Store the number of degrees of freedom in the non-augmented problem. More...
 
LinearAlgebraDistributionDof_distribution_pt
 Store the original dof distribution. More...
 
LinearAlgebraDistributionAugmented_dof_distribution_pt
 The augmented distribution. More...
 
double Sigma
 A slack variable used to specify the amount of antisymmetry in the solution. More...
 
DoubleVectorWithHaloEntries Y
 Storage for the null vector. More...
 
DoubleVectorWithHaloEntries Psi
 A constant vector that is specifies the symmetry being broken. More...
 
DoubleVectorWithHaloEntries C
 A constant vector used to ensure that the null vector is not trivial. More...
 
DoubleVectorWithHaloEntries Count
 A vector that is used to determine how many elements contribute to a particular equation. It is used to ensure that the global system is correctly formulated. This should really be an integer, but its double so that the distribution can be used. More...
 
Vector< unsigned > Global_eqn_number
 A vector that is used to map the global equations to their actual location in a distributed problem. More...
 
double * Parameter_pt
 Storage for the pointer to the parameter. More...
 
unsigned Nelement
 
bool Distributed
 Boolean to indicate whether the problem is distributed. More...
 

Friends

class BlockPitchForkLinearSolver
 
class AugmentedBlockPitchForkLinearSolver
 

Detailed Description

A class that is used to assemble the augmented system that defines a pitchfork (symmetry-breaking) bifurcation. The "standard" problem must be a function of a global parameter $ \lambda $ and a solution is $R(u,\lambda) = 0$, where $u$ are the unknowns in the problem. A pitchfork bifurcation may be specified by the augmented system of size $2N+2$.

\[ R(u,\lambda) + \sigma \psi = 0,\]

\[ Jy = 0,\]

\[ \langle u, \phi \rangle = 0, \]

\[ \phi \cdot y = 1.\]

In the abovem $J$ is the usual Jacobian matrix, $dR/du$ and $\phi$ is a constant vector that is chosen to ensure that the null vector, $y$, is not trivial. Here $\sigma $ is a slack variable that is used to enforce the constraint $ \langle u, \phi \rangle = 0 $ — the inner product of the solution with the chosen symmetry vector is zero. At the bifurcation $ \sigma $ should be very close to zero. Note that this formulation means that any odd symmetry in the problem must be about zero. Moreover, we use the dot product of two vectors to calculate the inner product, rather than an integral over the domain and so the mesh must be symmetric in the appropriate directions.

Definition at line 772 of file assembly_handler.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
private

A little private enum to determine whether we are solving the block system or not.

Enumerator
Full_augmented 
Block_J 
Block_augmented_J 

Definition at line 779 of file assembly_handler.h.

Constructor & Destructor Documentation

◆ PitchForkHandler()

oomph::PitchForkHandler::PitchForkHandler ( Problem *const &  problem_pt,
AssemblyHandler *const &  assembly_handler_pt,
double *const &  parameter_pt,
const DoubleVector symmetry_vector 
)

◆ ~PitchForkHandler()

oomph::PitchForkHandler::~PitchForkHandler ( )

Member Function Documentation

◆ bifurcation_parameter_pt()

double* oomph::PitchForkHandler::bifurcation_parameter_pt ( ) const
inlinevirtual

Return a pointer to the bifurcation parameter in bifurcation tracking problems.

Reimplemented from oomph::AssemblyHandler.

Definition at line 924 of file assembly_handler.h.

References Parameter_pt.

◆ bifurcation_type()

int oomph::PitchForkHandler::bifurcation_type ( ) const
inlinevirtual

Indicate that we are tracking a pitchfork bifurcation by returning 2.

Reimplemented from oomph::AssemblyHandler.

Definition at line 917 of file assembly_handler.h.

◆ eqn_number()

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

◆ get_djacobian_dparameter()

void oomph::PitchForkHandler::get_djacobian_dparameter ( GeneralisedElement *const &  elem_pt,
double *const &  parameter_pt,
Vector< double > &  dres_dparam,
DenseMatrix< double > &  djac_dparam 
)
virtual

Overload the derivative of the residuals and jacobian with respect to a parameter so that it breaks.

Overload the derivative of the residuals and jacobian with respect to a parameter so that it breaks because it should not be required.

Reimplemented from oomph::AssemblyHandler.

Definition at line 3394 of file assembly_handler.cc.

◆ get_dresiduals_dparameter()

void oomph::PitchForkHandler::get_dresiduals_dparameter ( GeneralisedElement *const &  elem_pt,
double *const &  parameter_pt,
Vector< double > &  dres_dparam 
)
virtual

◆ get_eigenfunction()

void oomph::PitchForkHandler::get_eigenfunction ( Vector< DoubleVector > &  eigenfunction)
virtual

Return the eigenfunction(s) associated with the bifurcation that has been detected in bifurcation tracking problems.

Reimplemented from oomph::AssemblyHandler.

Definition at line 3430 of file assembly_handler.cc.

References Dof_distribution_pt, and Y.

◆ get_hessian_vector_products()

void oomph::PitchForkHandler::get_hessian_vector_products ( GeneralisedElement *const &  elem_pt,
Vector< double > const &  Y,
DenseMatrix< double > const &  C,
DenseMatrix< double > &  product 
)
virtual

Overload the hessian vector product function so that it breaks.

Overload the hessian vector product function so that it calls the underlying element's hessian function.

Reimplemented from oomph::AssemblyHandler.

Definition at line 3416 of file assembly_handler.cc.

References C, oomph::GeneralisedElement::get_hessian_vector_products(), and Y.

◆ get_jacobian()

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

◆ get_residuals()

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

◆ global_eqn_number()

unsigned oomph::PitchForkHandler::global_eqn_number ( const unsigned &  i)
inlineprivate

Function that is used to return map the global equations using the simplistic numbering scheme into the actual distributed scheme.

Definition at line 848 of file assembly_handler.h.

References Distributed, Global_eqn_number, and i.

Referenced by eqn_number().

◆ ndof()

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

◆ solve_augmented_block_system()

void oomph::PitchForkHandler::solve_augmented_block_system ( )

◆ solve_block_system()

void oomph::PitchForkHandler::solve_block_system ( )

Set to solve the block system.

Set to solve the block system. The boolean flag is used to specify whether the block decomposition should use exactly the same jacobian as the original system.

Definition at line 3512 of file assembly_handler.cc.

References Block_J, Dof_distribution_pt, oomph::Problem::Dof_distribution_pt, oomph::Problem::Dof_pt, oomph::LinearAlgebraDistribution::nrow_local(), Problem_pt, Solve_which_system, and oomph::Problem::Sparse_assemble_with_arrays_previous_allocation.

Referenced by oomph::BlockPitchForkLinearSolver::resolve(), and oomph::BlockPitchForkLinearSolver::solve().

◆ solve_full_system()

void oomph::PitchForkHandler::solve_full_system ( )

◆ synchronise()

void oomph::PitchForkHandler::synchronise ( )
virtual

Function that is used to perform any synchronisation required during the solution.

Reimplemented from oomph::AssemblyHandler.

Definition at line 3448 of file assembly_handler.cc.

References oomph::Problem::communicator_pt(), Distributed, Parameter_pt, Problem_pt, Sigma, oomph::DoubleVectorWithHaloEntries::synchronise(), and Y.

Friends And Related Function Documentation

◆ AugmentedBlockPitchForkLinearSolver

Definition at line 775 of file assembly_handler.h.

◆ BlockPitchForkLinearSolver

friend class BlockPitchForkLinearSolver
friend

Definition at line 774 of file assembly_handler.h.

Member Data Documentation

◆ Assembly_handler_pt

AssemblyHandler* oomph::PitchForkHandler::Assembly_handler_pt
private

Pointer to the underlying (original) assembly handler.

Definition at line 797 of file assembly_handler.h.

Referenced by get_jacobian(), and PitchForkHandler().

◆ Augmented_dof_distribution_pt

LinearAlgebraDistribution* oomph::PitchForkHandler::Augmented_dof_distribution_pt
private

◆ C

DoubleVectorWithHaloEntries oomph::PitchForkHandler::C
private

A constant vector used to ensure that the null vector is not trivial.

Definition at line 821 of file assembly_handler.h.

Referenced by get_hessian_vector_products(), get_jacobian(), get_residuals(), and PitchForkHandler().

◆ Count

DoubleVectorWithHaloEntries oomph::PitchForkHandler::Count
private

A vector that is used to determine how many elements contribute to a particular equation. It is used to ensure that the global system is correctly formulated. This should really be an integer, but its double so that the distribution can be used.

Definition at line 828 of file assembly_handler.h.

Referenced by get_jacobian(), get_residuals(), and PitchForkHandler().

◆ Distributed

bool oomph::PitchForkHandler::Distributed
private

Boolean to indicate whether the problem is distributed.

Definition at line 842 of file assembly_handler.h.

Referenced by eqn_number(), global_eqn_number(), PitchForkHandler(), solve_full_system(), and synchronise().

◆ Dof_distribution_pt

LinearAlgebraDistribution* oomph::PitchForkHandler::Dof_distribution_pt
private

◆ Global_eqn_number

Vector<unsigned> oomph::PitchForkHandler::Global_eqn_number
private

A vector that is used to map the global equations to their actual location in a distributed problem.

Definition at line 832 of file assembly_handler.h.

Referenced by global_eqn_number(), and PitchForkHandler().

◆ Ndof

unsigned oomph::PitchForkHandler::Ndof
private

Store the number of degrees of freedom in the non-augmented problem.

Definition at line 801 of file assembly_handler.h.

Referenced by eqn_number(), PitchForkHandler(), and solve_augmented_block_system().

◆ Nelement

unsigned oomph::PitchForkHandler::Nelement
private

Definition at line 838 of file assembly_handler.h.

Referenced by get_residuals(), and PitchForkHandler().

◆ Parameter_pt

double* oomph::PitchForkHandler::Parameter_pt
private

◆ Problem_pt

Problem* oomph::PitchForkHandler::Problem_pt
private

◆ Psi

DoubleVectorWithHaloEntries oomph::PitchForkHandler::Psi
private

◆ Sigma

double oomph::PitchForkHandler::Sigma
private

A slack variable used to specify the amount of antisymmetry in the solution.

Definition at line 811 of file assembly_handler.h.

Referenced by get_jacobian(), get_residuals(), PitchForkHandler(), solve_full_system(), and synchronise().

◆ Solve_which_system

unsigned oomph::PitchForkHandler::Solve_which_system
private

Integer flag to indicate which system should be assembled. There are three possibilities. The full augmented system (0), the non-augmented jacobian system (1), a system in which the jacobian is augmented by 1 row and column to ensure that it is non-singular (2). See the enum above.

Definition at line 791 of file assembly_handler.h.

Referenced by eqn_number(), get_dresiduals_dparameter(), get_jacobian(), get_residuals(), ndof(), solve_augmented_block_system(), solve_block_system(), and solve_full_system().

◆ Y

DoubleVectorWithHaloEntries oomph::PitchForkHandler::Y
private

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