CppNoddy  0.85
 All Classes Namespaces Files Functions Variables Typedefs Enumerator Friends Macros Groups Pages
Public Member Functions | List of all members
CppNoddy::PDE_double_IBVP< _Type > Class Template Reference

A templated object for real/complex vector system of unsteady equations. More...

#include <PDE_double_IBVP.h>

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

Public Member Functions

 PDE_double_IBVP (Equation_3matrix< _Type > *equation_ptr, const DenseVector< double > &xnodes, const DenseVector< double > &ynodes, Residual_with_coords< _Type > *ptr_to_bottom_residual, Residual_with_coords< _Type > *ptr_to_top_residual)
 The class is defined by a vector function for the system.
 ~PDE_double_IBVP ()
 Destructor.
void update_previous_solution ()
 Copy the current solution to the previous solution.
void step2 (const double &dt)
 A Crank-Nicolson 'time' stepper.
double & t ()
 Return a reference to the current value of the 'timelike' coordinate.
double & tolerance ()
 Return a reference to the convergence tolerance.
TwoD_Node_Mesh< _Type > & solution ()

Additional Inherited Members

Detailed Description

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

A templated object for real/complex vector system of unsteady equations.

Definition at line 45 of file PDE_double_IBVP.h.

Constructor & Destructor Documentation

template<typename _Type >
CppNoddy::PDE_double_IBVP< _Type >::PDE_double_IBVP ( Equation_3matrix< _Type > *  equation_ptr,
const DenseVector< double > &  xnodes,
const DenseVector< double > &  ynodes,
Residual_with_coords< _Type > *  ptr_to_bottom_residual,
Residual_with_coords< _Type > *  ptr_to_top_residual 
)

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

Parameters
equation_ptrA pointer to an inherited Equation object.
xnodesA vector that defines the nodal x-positions.
ynodesA vector that defines the nodal y-positions.
ptr_to_bottom_residualA pointer to a residual object that defines the y=y1 boundary conditions.
ptr_to_top_residualA pointer to a residual object that defines the y=y2 boundary conditions.

Definition at line 19 of file PDE_double_IBVP.cpp.

:
TOL( 1.e-8 ),
T( 0.0 ),
MAX_ITERATIONS( 12 ),
p_EQUATION( ptr_to_equation ),
p_BOTTOM_RESIDUAL( ptr_to_bottom_residual ),
p_TOP_RESIDUAL( ptr_to_top_residual ),
UPDATED( false )
{
// create the 2D mesh for the current time level and previous time level
SOLN = TwoD_Node_Mesh<_Type>( xnodes, ynodes, p_EQUATION -> get_order() );
// copy construct the previous solution's storage
PREV_SOLN = SOLN;
// initialise the eqn object
// "x"
p_EQUATION -> coord(2) = xnodes[0];
#ifdef TIME
// timers
T_ASSEMBLE = Timer( "Assembling of the matrix (incl. equation updates):" );
T_SOLVE = Timer( "Solving of the matrix:" );
#endif
}
template<typename _Type >
CppNoddy::PDE_double_IBVP< _Type >::~PDE_double_IBVP ( )

Destructor.

Definition at line 48 of file PDE_double_IBVP.cpp.

{
#ifdef TIME
std::cout << "\n";
T_ASSEMBLE.stop();
T_ASSEMBLE.print();
T_SOLVE.stop();
T_SOLVE.print();
#endif
}

Member Function Documentation

template<typename _Type >
TwoD_Node_Mesh< _Type > & CppNoddy::PDE_double_IBVP< _Type >::solution ( )
inline
Returns
A handle to the solution mesh

Definition at line 127 of file PDE_double_IBVP.h.

Referenced by main().

{
return SOLN;
}
template<typename _Type >
void CppNoddy::PDE_double_IBVP< _Type >::step2 ( const double &  dt)

A Crank-Nicolson 'time' stepper.

Definition at line 60 of file PDE_double_IBVP.cpp.

References CppNoddy::DenseVector< _Type >::inf_norm(), and CppNoddy::BandedLinearSystem< _Type >::solve().

Referenced by main().

{
DT = dt;
if ( !UPDATED )
{
std::string problem;
problem = " You have called the PDE_double_IBVP::step2 method without calling\n";
problem += " the PDE_double_IBVP::update_previous_solution method first. This\n";
problem += " method is required to copy/store the previous time level's solution\n";
problem += " and then allow you to specify/alter the x-initial condition by\n";
problem += " for the current time level.";
throw ExceptionRuntime( problem );
}
// the order of the problem
unsigned order( p_EQUATION -> get_order() );
// get the number of nodes in the mesh
// -- this may have been refined by the user since the last call.
unsigned nx( SOLN.get_nnodes().first );
unsigned ny( SOLN.get_nnodes().second );
// measure of maximum residual
double max_residual( 1.0 );
// the current solution moves to the previous solution
// ANY LARGE STORAGE USED IN THE MAIN LOOP IS
// DEFINED HERE TO AVOID REPEATED CONSTRUCTION.
// Note we blank the A matrix after every iteration.
//
// Banded LHS matrix - max obove diagonal band width is
// from first variable at node i to last variable at node i+1
BandedMatrix<_Type> a( ny * order, 2 * order - 1, 0.0 );
// RHS
DenseVector<_Type> b( ny * order, 0.0 );
// linear solver definition
#ifdef LAPACK
BandedLinearSystem<_Type> system( &a, &b, "lapack" );
#else
BandedLinearSystem<_Type> system( &a, &b, "native" );
#endif
for ( std::size_t j = 0; j < nx - 1; ++j )
{
//// next x-soln guess is the previous x-soln
//for ( std::size_t var = 0; var < order; ++var )
//{
//for ( std::size_t i = 0; i < ny; ++i )
//{
//SOLN( j + 1, i, var ) = SOLN( j, i, var );
//}
//}
//
// iteration counter
int counter = 0;
// loop until converged or too many iterations
do
{
// iteration counter
++counter;
// time the assemble phase
#ifdef TIME
T_ASSEMBLE.start();
#endif
assemble_matrix_problem( a, b, j );
max_residual = b.inf_norm();
#ifdef DEBUG
std::cout.precision( 12 );
std::cout << " PDE_double_IBVP.step2 : x_j = " << SOLN.coord(j,0).first << " Residual_max = " << max_residual << " tol = " << TOL << " counter = " << counter << "\n";
#endif
#ifdef TIME
T_ASSEMBLE.stop();
T_SOLVE.start();
#endif
// linear solver
system.solve();
// keep the solution in a OneD_GenMesh object
for ( std::size_t i = 0; i < ny; ++i )
{
for ( std::size_t var = 0; var < order; ++var )
{
SOLN( j + 1, i, var ) += b[ i * order + var ];
}
}
#ifdef TIME
T_SOLVE.stop();
#endif
}
while ( ( max_residual > TOL ) && ( counter < MAX_ITERATIONS ) );
if ( counter >= MAX_ITERATIONS )
{
std::string problem( "\n The PDE_double_IBVP.step2 method took too many iterations. \n" );
throw ExceptionItn( problem, counter, max_residual );
}
}
// set the time to the updated level
T += DT;
#ifdef DEBUG
std::cout << "[DEBUG] solved at t = " << T << "\n";
#endif
UPDATED = false;
}
template<typename _Type >
double & CppNoddy::PDE_double_IBVP< _Type >::t ( )
inline

Return a reference to the current value of the 'timelike' coordinate.

Returns
A handle to the current timelinke coordinate stored in the object

Definition at line 133 of file PDE_double_IBVP.h.

Referenced by main().

{
return T;
}
template<typename _Type>
double& CppNoddy::PDE_double_IBVP< _Type >::tolerance ( )
inline

Return a reference to the convergence tolerance.

Definition at line 79 of file PDE_double_IBVP.h.

{
return TOL;
}
template<typename _Type>
void CppNoddy::PDE_double_IBVP< _Type >::update_previous_solution ( )
inline

Copy the current solution to the previous solution.

Definition at line 65 of file PDE_double_IBVP.h.

Referenced by main().

{
PREV_SOLN = SOLN;
UPDATED = true;
}

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

© 2012

R.E. Hewitt