CppNoddy  0.90
Public Member Functions | List of all members
CppNoddy::ODE_EVP< _Type > Class Template Reference

A templated object for real/complex vector system of first-order ordinary differential equations. More...

#include <ODE_EVP.h>

Inheritance diagram for CppNoddy::ODE_EVP< _Type >:
CppNoddy::Uncopyable

Public Member Functions

 ODE_EVP (Equation_2matrix< _Type > *equation_ptr, const DenseVector< double > &nodes, Residual< _Type > *ptr_to_left_residual, Residual< _Type > *ptr_to_right_residual)
 The class is defined by a vector function for the system. More...
 
 ~ODE_EVP ()
 Destructor. More...
 
void eigensolve ()
 Formulate and solve the global eigenvalue problem for a linear system. More...
 
LinearEigenSystem_basep_eigensystem ()
 Allow access to the underlying dense linear eigensystem through a pointer to the private member data. More...
 
void add_tagged_to_mesh ()
 
OneD_Node_Mesh< D_complexget_mesh (const unsigned &i) const
 

Detailed Description

template<typename _Type>
class CppNoddy::ODE_EVP< _Type >

A templated object for real/complex vector system of first-order ordinary differential equations.

Definition at line 36 of file ODE_EVP.h.

Constructor & Destructor Documentation

template<typename _Type >
CppNoddy::ODE_EVP< _Type >::ODE_EVP ( Equation_2matrix< _Type > *  equation_ptr,
const DenseVector< double > &  nodes,
Residual< _Type > *  ptr_to_left_residual,
Residual< _Type > *  ptr_to_right_residual 
)

The class is defined by a vector function for the system.

Parameters
equation_ptrA pointer to an equation with 2 associated matrices; matrix1 will define the eigenvalue problem.
nodesA vector of nodal points.
ptr_to_left_residualA pointer to a residual object that defines the LHS boundary conditions.
ptr_to_right_residualA pointer to a residual object that defines the RHS boundary conditions.

Definition at line 27 of file ODE_EVP.cpp.

References CppNoddy::DenseVector< _Type >::size().

30  :
31  p_EQUATION( ptr_to_equation ),
32  p_LEFT_RESIDUAL( ptr_to_left_residual ),
33  p_RIGHT_RESIDUAL( ptr_to_right_residual ),
34  NODES( nodes )
35  {
36  unsigned n( nodes.size() );
37  unsigned order( p_EQUATION -> get_order() );
38  // first make the dense matrix equivalent of the problem
39  // using the Dense( banded ) constructor
40  p_A_DENSE = new DenseMatrix<_Type>( n * order, n * order, 0.0 );
41  p_B_DENSE = new DenseMatrix<_Type>( n * order, n * order, 0.0 );
42  // point the system pointer to the dense problem
43  p_SYSTEM = new DenseLinearEigenSystem<_Type>( p_A_DENSE, p_B_DENSE );
44  }
std::size_t size() const
A pass-thru definition to get the size of the vector.
Definition: DenseVector.h:365
template<typename _Type >
CppNoddy::ODE_EVP< _Type >::~ODE_EVP ( )

Destructor.

Definition at line 47 of file ODE_EVP.cpp.

48  {
49  // clean up the matrices & evp system
50  delete p_SYSTEM;
51  delete p_A_DENSE;
52  delete p_B_DENSE;
53  }

Member Function Documentation

template<typename _Type>
void CppNoddy::ODE_EVP< _Type >::add_tagged_to_mesh ( )
inline

Definition at line 62 of file ODE_EVP.h.

References CppNoddy::OneD_Node_Mesh< _Type, _Xtype >::set_nodes_vars(), and CppNoddy::DenseVector< _Type >::size().

Referenced by main(), and CppNoddy::Example::Neutral_residual::Neutral_residual().

63  {
64  // clear the existing data if any
65  MESHES.clear();
66  // order of the equation
67  unsigned order = p_EQUATION -> get_order();
68  // get the eigenvalues
69  DenseVector<D_complex> vals( p_SYSTEM -> get_tagged_eigenvalues() );
70  // get the eigenvectors
71  DenseMatrix<D_complex> vecs( p_SYSTEM -> get_tagged_eigenvectors() );
72  // loop through the eigenvectors
73  for ( unsigned ivec = 0; ivec < vals.size(); ++ivec )
74  {
75  // make a mesh with the right node distribution
76  // we'll increase the order by 1 to allow the eigenvalue to be
77  // stored at each nodal point -- this is wasteful but very useful
78  // for feeding this into a BVP for local refinement.
79  OneD_Node_Mesh<D_complex> eigfn( NODES, order + 1 );
80  // loop through all nodes
81  for ( unsigned node = 0; node < NODES.size(); ++node )
82  {
83  // complex vector of the dof at this node ( + 1 for the eigenvalue)
84  DenseVector<D_complex> vars_at_node( order + 1, 0.0 );
85  // get the dof from the eigenvector
86  for ( unsigned var = 0; var < order; ++var )
87  {
88  vars_at_node[ var ] = vecs[ ivec ][ node * order + var ];
89  }
90  // the last variable at each node is the corresponding eigenvalue
91  vars_at_node[ order ] = vals[ ivec ];
92  // set the first 'order' dof to be those from the eigenvector
93  eigfn.set_nodes_vars( node, vars_at_node );
94  }
95  //// store the eigenvalue in the mesh at each node ... wasteful, but useful
96  //// a complex vector filled with the same value N times
97  //DenseVector<D_complex> c( NODES.size(), vals[ ivec ] );
98  //// add it to the mesh -- for use in nonlinear local refinement via ODE_BVP
99  //eigfn.push_var( c );
100  // add the eigenfunction to the vector of meshes
101  MESHES.push_back( eigfn );
102  }
103  }
std::size_t size() const
A pass-thru definition to get the size of the vector.
Definition: DenseVector.h:365
template<typename _Type >
void CppNoddy::ODE_EVP< _Type >::eigensolve ( )

Formulate and solve the global eigenvalue problem for a linear system.

Definition at line 62 of file ODE_EVP.cpp.

References h, and CppNoddy::DenseVector< _Type >::size().

Referenced by main(), and CppNoddy::Example::Neutral_residual::Neutral_residual().

63  {
64  // NOTE: moving construct_banded_system to the constructor is a bad idea, because
65  // sometimes we will want to construct an EVP that depends on an
66  // underlying baseflow solution that is unknown at the time of construction.
67  //
68  // construct the banded matrix eigenvalue problem that corresponds
69  // to the equation & boundary conditions specified in the constructor.
70  assemble_dense_problem();
71  // solve the system
72  p_SYSTEM -> eigensolve();
73  }
void eigensolve()
Formulate and solve the global eigenvalue problem for a linear system.
Definition: ODE_EVP.cpp:62
template<typename _Type>
OneD_Node_Mesh<D_complex> CppNoddy::ODE_EVP< _Type >::get_mesh ( const unsigned &  i) const
inline

Definition at line 105 of file ODE_EVP.h.

Referenced by main(), and CppNoddy::Example::Neutral_residual::Neutral_residual().

106  {
107 #ifdef PARANOID
108  if ( i > MESHES.size() )
109  {
110  std::string problem;
111  problem = "You have tried to extract an eigenfunction from the ODE_EVP class\n";
112  problem += "whose index is outside the range of stored meshes.\n";
113  throw ExceptionRange( problem, MESHES.size(), i );
114  }
115 #endif
116  return MESHES[ i ];
117  }
template<typename _Type >
LinearEigenSystem_base * CppNoddy::ODE_EVP< _Type >::p_eigensystem ( )

Allow access to the underlying dense linear eigensystem through a pointer to the private member data.

Definition at line 56 of file ODE_EVP.cpp.

Referenced by main(), and CppNoddy::Example::Neutral_residual::Neutral_residual().

57  {
58  return p_SYSTEM;
59  }

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

© 2012

R.E. Hewitt