
This class provides an interface to the sparse direct solver MA47 from the Harwell Subroutine Library. MA47 is a direct solver specialized for sparse symmetric indefinite systems of linear equations and uses a frontal elimination method. It is included in the Harwell Subroutine Library and is written in Fortran. The present class only transforms the data stored in SparseMatrix objects into the form which is required by the functions resembling MA47, calls these Fortran functions, and interprets some of the returned values indicating error codes, etc. It also manages allocation of the right amount of temporary storage required by these functions.
Note that this class only works if configuration of the deal.II library has detected the presence of this solver. Please read the README file on what the configure script is looking for and how to provide it.
For the meaning of the three functions initialize(), factorize(), and solve(), as well as for the method used in MA47, please see the documentation of these functions. In practice, one will most often call the second solve() function, which solves the linear system for a given right hand side, but one can as well call the three functions separately if, for example, one would like to solve the same matrix for several right hand side vectors; the MA47 solver can do this efficiently, as it computes a decomposition of the matrix, so that subsequent solves only amount to a forward-backward substitution which is significantly less costly than the decomposition process.
The constructor of this class takes several arguments. Their meaning is equivalent to those of the constructor of the SparseDirectMA27 class; see there for more information.
Due to the use of global variables through COMMON blocks, the calls to the sparse direct solver routines is not multithreading-capable, i.e. at each time there may only be one call to these functions active. You have to synchronise your calls to the functions provided by this class using mutexes (see the Threads namespace for such classes) to avoid multiple active calls at the same time if you use multithreading. Since you may use this class in different parts of your program, and may not want to use a global variable for locking, this class has a lock as static member variable, which may be accessed using the get_synchronisation_lock() function. Note however, that this class does not perform the synchronisation for you within its member functions. The reason is that you will usually want to synchronise over the calls to initialize() and factorize(), since there should probably not be a call to one of these function with another matrix between the calls for one matrix. (The author does not really know whether this is true, but it is probably safe to assume that.) Since such cross-function synchronisation can only be performed from outside, it is left to the user of this class to do so.
A detached mode as for MA27 has not yet been implemented for this class.
Definition at line 681 of file sparse_direct.h.
| SparseDirectMA47::SparseDirectMA47 | ( | const double | LIW_factor_1 = 1.4, |
| const double | LIW_factor_2 = 1, |
||
| const double | LA_factor = 3, |
||
| const double | LIW_increase_factor_1 = 1.2, |
||
| const double | LIW_increase_factor_2 = 1.2, |
||
| const double | LA_increase_factor = 1.2, |
||
| const bool | suppress_output = true |
||
| ) |
Constructor. See the documentation of this class for the meaning of the parameters to this function.
This function already calls the initialization function MA47ID to set up some values.
| void SparseDirectMA47::initialize | ( | const SparseMatrix< double > & | matrix ) |
Initialize some data structures. This function computes symbolically some information based on the sparsity pattern, but does not actually use the values of the matrix, so only the sparsity pattern has to be passed as argument.
Since the MA47 solver requires us to omit zero-entries in the matrix (even if they are in the sparsity pattern), we have to actually use the matrix here, as opposed to the MA27 solver that only required the sparsity pattern.
| void SparseDirectMA47::factorize | ( | const SparseMatrix< double > & | matrix ) |
Actually factorize the matrix. Unlike for the MA27 solver, this function may not be called multiple times for different matrices, since we have eliminated entries from the sparsity pattern where matrix entries happen to be zero. Since this is likely to change between matrices although they have the same sparsity pattern.
If the initialization step has not been performed yet, then the initialize() function is called at the beginning of this function.
Solve for a certain right hand side vector. This function may be called multiple times for different right hand side vectors after the matrix has been factorized. This yields a big saving in computing time, since the actual solution is fast, compared to the factorization of the matrix.
The solution will be returned in place of the right hand side vector.
If the factorization has not happened before, strange things will happen. Note that we can't actually call the factorize() function from here if it has not yet been called, since we have no access to the actual matrix.
| void SparseDirectMA47::solve | ( | const SparseMatrix< double > & | matrix, |
| Vector< double > & | rhs_and_solution | ||
| ) |
Call the three functions initialize, factorize and solve in that order, i.e. perform the whole solution process for the given right hand side vector.
The solution will be returned in place of the right hand side vector.
| std::size_t SparseDirectMA47::memory_consumption | ( | ) | const |
Return an estimate of the memory used by this class.
| Threads::ThreadMutex& SparseDirectMA47::get_synchronisation_lock | ( | ) | const |
Get a reference to the synchronisation lock which can be used for this class. See the general description of this class for more information.
| void SparseDirectMA47::fill_A | ( | const SparseMatrix< double > & | matrix ) | [private] |
Fill the A array from the symmetric part of the given matrix.
Call the ma47id function with the given args.
| void SparseDirectMA47::call_ma47ad | ( | const unsigned int * | n_rows, |
| const unsigned int * | n_nonzero_elements, | ||
| unsigned int * | row_numbers, | ||
| unsigned int * | column_numbers, | ||
| unsigned int * | IW, | ||
| const unsigned int * | LIW, | ||
| unsigned int * | KEEP, | ||
| const unsigned int * | ICNTL, | ||
| int * | INFO | ||
| ) | [private] |
Call the ma47ad function with the given args.
| void SparseDirectMA47::call_ma47bd | ( | const unsigned int * | n_rows, |
| const unsigned int * | n_nonzero_elements, | ||
| const unsigned int * | column_numbers, | ||
| double * | A, | ||
| const unsigned int * | LA, | ||
| unsigned int * | IW, | ||
| const unsigned int * | LIW, | ||
| const unsigned int * | KEEP, | ||
| const double * | CNTL, | ||
| const unsigned int * | ICNTL, | ||
| unsigned int * | IW1, | ||
| int * | INFO | ||
| ) | [private] |
Call the ma47bd function with the given args.
| void SparseDirectMA47::call_ma47cd | ( | const unsigned int * | n_rows, |
| const double * | A, | ||
| const unsigned int * | LA, | ||
| const unsigned int * | IW, | ||
| const unsigned int * | LIW, | ||
| double * | rhs_and_solution, | ||
| unsigned int * | IW1, | ||
| const unsigned int * | ICNTL | ||
| ) | [private] |
Call the ma47bd function with the given args.
const bool SparseDirectMA47::suppress_output [private] |
Store in the constructor whether the MA47 routines shall deliver output to stdout or not.
Definition at line 858 of file sparse_direct.h.
const double SparseDirectMA47::LIW_factor_1 [private] |
Store the three values passed to the cinstructor. See the documentation of this class for the meaning of these variables.
Definition at line 867 of file sparse_direct.h.
const double SparseDirectMA47::LIW_factor_2 [private] |
Definition at line 868 of file sparse_direct.h.
const double SparseDirectMA47::LA_factor [private] |
Definition at line 869 of file sparse_direct.h.
const double SparseDirectMA47::LIW_increase_factor_1 [private] |
Increase factors in case a call to a function fails.
Definition at line 875 of file sparse_direct.h.
const double SparseDirectMA47::LIW_increase_factor_2 [private] |
Definition at line 876 of file sparse_direct.h.
const double SparseDirectMA47::LA_increase_factor [private] |
Definition at line 877 of file sparse_direct.h.
bool SparseDirectMA47::initialize_called [private] |
Flags storing whether the first two functions have already been called.
Definition at line 884 of file sparse_direct.h.
bool SparseDirectMA47::factorize_called [private] |
Definition at line 885 of file sparse_direct.h.
SmartPointer<const SparseMatrix<double>,SparseDirectMA47> SparseDirectMA47::matrix [private] |
Store a pointer to the matrix, to make sure that we use the same thing for all calls.
Definition at line 892 of file sparse_direct.h.
unsigned int SparseDirectMA47::n_nonzero_elements [private] |
Number of nonzero elements in the sparsity pattern on and above the diagonal.
Definition at line 899 of file sparse_direct.h.
double SparseDirectMA47::CNTL[2] [private] |
Control values set by MA47ID.
Definition at line 904 of file sparse_direct.h.
unsigned int SparseDirectMA47::ICNTL[7] [private] |
Definition at line 905 of file sparse_direct.h.
int SparseDirectMA47::INFO[24] [private] |
Info field filled by the MA47 functions and (partially) used for subsequent MA47 calls.
Definition at line 912 of file sparse_direct.h.
std::vector<unsigned int> SparseDirectMA47::row_numbers [private] |
Arrays holding row and column indices.
Definition at line 918 of file sparse_direct.h.
std::vector<unsigned int> SparseDirectMA47::column_numbers [private] |
Definition at line 919 of file sparse_direct.h.
std::vector<double> SparseDirectMA47::A [private] |
Array to hold the matrix elements, and later the elements of the factors.
Definition at line 926 of file sparse_direct.h.
unsigned int SparseDirectMA47::LA [private] |
Length of the A array.
Definition at line 931 of file sparse_direct.h.
unsigned int SparseDirectMA47::LIW [private] |
Scratch arrays and variables used by the MA47 functions. We keep to the names introduced in the documentation of these functions, in all uppercase letters as is usual in Fortran.
Definition at line 942 of file sparse_direct.h.
std::vector<unsigned int> SparseDirectMA47::IW [private] |
Definition at line 943 of file sparse_direct.h.
std::vector<unsigned int> SparseDirectMA47::KEEP [private] |
Definition at line 944 of file sparse_direct.h.
std::vector<unsigned int> SparseDirectMA47::IW1 [private] |
Definition at line 945 of file sparse_direct.h.
Threads::ThreadMutex SparseDirectMA47::synchronisation_lock [static, private] |
Mutex for synchronising access to this class.
Definition at line 951 of file sparse_direct.h.
documentation generated on Fri Feb 3 2012 06:04:11 by
doxygen
1.7.2