Virtual base class. Abstracts a tridiagonal system to solve in a column of ice and/or bedrock. More...

#include <ColumnSystem.hh>

Public Member Functions
	TridiagonalSystem (unsigned int max_size, const std::string &prefix)
	Allocate a tridiagonal system of maximum size max_system_size. More...

double	norm1 (unsigned int system_size) const
	Compute 1-norm, which is max sum of absolute values of columns. More...

double	ddratio (unsigned int system_size) const
	Compute diagonal-dominance ratio. If this is less than one then the matrix is strictly diagonally-dominant. More...

void	reset ()
	Zero all entries. More...

void	solve (unsigned int system_size, std::vector< double > &result)

void	solve (unsigned int system_size, double *result)
	The actual code for solving a tridiagonal system. More...

void	save_system_with_solution (const std::string &filename, unsigned int system_size, const std::vector< double > &solution)

void	save_system (std::ostream &output, unsigned int system_size) const
	View the tridiagonal system A x = b to an output stream, both A as a full matrix and b as a vector. More...

void	save_matrix (std::ostream &output, unsigned int system_size, const std::string &variable) const
	View the tridiagonal matrix. More...

std::string	prefix () const

double &	L (size_t i)

double &	D (size_t i)

double &	U (size_t i)

double &	RHS (size_t i)

Static Public Member Functions
static void	save_vector (std::ostream &output, const std::vector< double > &v, unsigned int system_size, const std::string &variable)
	Utility for simple ascii view of a vector (one-dimensional column) quantity. More...

Private Attributes
unsigned int	m_max_system_size

std::vector< double >	m_L

std::vector< double >	m_D

std::vector< double >	m_U

std::vector< double >	m_rhs

std::vector< double >	m_work

std::string	m_prefix

Detailed Description

Virtual base class. Abstracts a tridiagonal system to solve in a column of ice and/or bedrock.

Because both the age evolution and conservation of energy equations require us to set up and solve a tridiagonal system of equations, this is structure is worth abstracting.

This base class just holds the tridiagonal system and the ability to solve it, but does not insert entries into the relevant matrix locations. Derived classes will actually set up instances of the system.

The sequence requires setting the column-independent (public) data members, calling the initAllColumns() routine, and then setting up and solving the system in each column.

The tridiagonal algorithm here comes from Numerical Recipes in C Section 2.4, page 50. It solves the system:

  [b_1  c_1   0   ...                           ] [  u_1  ]   [   r_1   ]
  [a_2  b_2  c_2  ...                           ] [  u_2  ]   [   r_2   ]
  [               ...                           ].[  ...  ] = [   ...   ]
  [               ...  a_{N-1}  b_{N-1}  c_{N-1}] [u_{N-1}]   [ r_{N-1} ]
  [               ...     0     a_N      b_N    ] [  u_N  ]   [   r_N   ]

HOWEVER... the version in this code is different from Numerical Recipes in two ways:

Indexing is zero-based
Variables have been renamed.

  NR      PISM
  ==================
  a       L       "Lower Diagonal" (L doesn't use index 0)
  b       D       "Diagonal"
  c       U       "Upper Diagonal"
  u       x
  r       rhs
  bet     b
  j       k
  n       n
  gam     work

Therefore... this version of the code solves the following problem:

  [D_0  U_0   0   ...                           ] [  x_0  ]   [   r_0   ]
  [L_1  D_1  U_1  ...                           ] [  x_1  ]   [   r_1   ]
  [               ...                           ].[  ...  ] = [   ...   ]
  [               ...  L_{N-2}  D_{N-2}  U_{N-2}] [x_{N-2}]   [ r_{N-2} ]
  [               ...     0     L_{N-1}  D_{N-1}] [x_{N-1}]   [ r_{N-1} ]

Definition at line 86 of file ColumnSystem.hh.

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

src/util/ColumnSystem.hh
src/util/ColumnSystem.cc

Public Member Functions

Static Public Member Functions

Private Attributes

Detailed Description