CppNoddy  0.85
 All Classes Namespaces Files Functions Variables Typedefs Enumerator Friends Macros Groups Pages
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.
 ~ODE_EVP ()
 Destructor.
void eigensolve ()
 Formulate and solve the global eigenvalue problem for a linear system.
LinearEigenSystem_basep_eigensystem ()
 Allow access to the underlying dense linear eigensystem through a pointer to the private member data.
void add_tagged_to_mesh ()
OneD_Node_Mesh< D_complexget_mesh (const unsigned &i) const

Additional Inherited Members

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 29 of file ODE_EVP.cpp.

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

:
p_EQUATION( ptr_to_equation ),
p_LEFT_RESIDUAL( ptr_to_left_residual ),
p_RIGHT_RESIDUAL( ptr_to_right_residual ),
NODES( nodes )
{
unsigned n( nodes.size() );
unsigned order( p_EQUATION -> get_order() );
// first make the dense matrix equivalent of the problem
// using the Dense( banded ) constructor
p_A_DENSE = new DenseMatrix<_Type>( n * order, n * order, 0.0 );
p_B_DENSE = new DenseMatrix<_Type>( n * order, n * order, 0.0 );
// point the system pointer to the dense problem
p_SYSTEM = new DenseLinearEigenSystem<_Type>( p_A_DENSE, p_B_DENSE );
}
template<typename _Type >
CppNoddy::ODE_EVP< _Type >::~ODE_EVP ( )

Destructor.

Definition at line 49 of file ODE_EVP.cpp.

{
// clean up the matrices & evp system
delete p_SYSTEM;
delete p_A_DENSE;
delete p_B_DENSE;
}

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().

{
// clear the existing data if any
MESHES.clear();
// order of the equation
unsigned order = p_EQUATION -> get_order();
// get the eigenvalues
DenseVector<D_complex> vals( p_SYSTEM -> get_tagged_eigenvalues() );
// get the eigenvectors
DenseMatrix<D_complex> vecs( p_SYSTEM -> get_tagged_eigenvectors() );
// loop through the eigenvectors
for ( unsigned ivec = 0; ivec < vals.size(); ++ivec )
{
// make a mesh with the right node distribution
// we'll increase the order by 1 to allow the eigenvalue to be
// stored at each nodal point -- this is wasteful but very useful
// for feeding this into a BVP for local refinement.
OneD_Node_Mesh<D_complex> eigfn( NODES, order + 1 );
// loop through all nodes
for ( unsigned node = 0; node < NODES.size(); ++node )
{
// complex vector of the dof at this node ( + 1 for the eigenvalue)
DenseVector<D_complex> vars_at_node( order + 1, 0.0 );
// get the dof from the eigenvector
for ( unsigned var = 0; var < order; ++var )
{
vars_at_node[ var ] = vecs[ ivec ][ node * order + var ];
}
// the last variable at each node is the corresponding eigenvalue
vars_at_node[ order ] = vals[ ivec ];
// set the first 'order' dof to be those from the eigenvector
eigfn.set_nodes_vars( node, vars_at_node );
}
//// store the eigenvalue in the mesh at each node ... wasteful, but useful
//// a complex vector filled with the same value N times
//DenseVector<D_complex> c( NODES.size(), vals[ ivec ] );
//// add it to the mesh -- for use in nonlinear local refinement via ODE_BVP
//eigfn.push_var( c );
// add the eigenfunction to the vector of meshes
MESHES.push_back( eigfn );
}
}
template<typename _Type >
void CppNoddy::ODE_EVP< _Type >::eigensolve ( )

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

Definition at line 64 of file ODE_EVP.cpp.

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

{
// NOTE: moving construct_banded_system to the constructor is a bad idea, because
// sometimes we will want to construct an EVP that depends on an
// underlying baseflow solution that is unknown at the time of construction.
//
// construct the banded matrix eigenvalue problem that corresponds
// to the equation & boundary conditions specified in the constructor.
assemble_dense_problem();
// solve the system
p_SYSTEM -> eigensolve();
}
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().

{
#ifdef PARANOID
if ( i > MESHES.size() )
{
std::string problem;
problem = "You have tried to extract an eigenfunction from the ODE_EVP class\n";
problem += "whose index is outside the range of stored meshes.\n";
throw ExceptionRange( problem, MESHES.size(), i );
}
#endif
return MESHES[ i ];
}
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 58 of file ODE_EVP.cpp.

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

{
return p_SYSTEM;
}

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

© 2012

R.E. Hewitt