Solver.h

Go to the documentation of this file.

/* **********************************************************************************************/
/*                                        stand-alone                                          */
/*                               Derived from RESEARCH VERSION 9.0                             */
/* **********************************************************************************************/
/* **********************************************************************************/
/*  Copyright WSL Institute for Snow and Avalanche Research    SLF-DAVOS           */
/* **********************************************************************************/
/* This file is part of Snowpack.
    Snowpack is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Snowpack is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with Snowpack.  If not, see <http://www.gnu.org/licenses/>.
*/

#ifndef  SOLVER_H
#define  SOLVER_H

#include <cstddef> //needed for size_t

typedef struct
{
        int     nChunks;
        int     pChunksSize;
        char   **pChunks;
        int     TotChunkSize;
} SD_CHUNK_DATA;

typedef struct
{
        int                Row0;
        int                Row1;
        size_t             nCol;
        int                nColBlock;
        int                iColBlock;
        int                iFloat;
}  SD_ROW_BLOCK_DATA;

typedef struct
{
        size_t              Dim;
        int                *pPerm;
        int                 nRowBlock;
        SD_ROW_BLOCK_DATA  *pRowBlock;
        int                 nColBlock;
        int                *pFirstColBlock;
        int                *pSizeColBlock;
        int                 SizeBlockJump;
        int                *pBlockJump;
        int                 SizeUpper;
        double              *pUpper;
}  SD_BLOCK_MATRIX_DATA;

typedef struct SD_COL_DATA
{
        size_t               Col;
        struct SD_COL_DATA  *Next;
} SD_COL_DATA;

typedef struct SD_ROW_DATA
{
        SD_COL_DATA  *Col;
}  SD_ROW_DATA;

typedef struct
{
        size_t                 nRow;
        int                   *pPerm;
        int                   *pPermInv;
        int                    nSupernode;
        int                   *pSupernode;

        SD_ROW_DATA           *pRow;
        SD_CHUNK_DATA          PoolCol;
        SD_COL_DATA           *FreeCol;
        int                    nFreeCol;
        size_t                 nCol;
}  SD_CON_MATRIX_DATA;

typedef struct SD_COL_BLOCK_DATA
{
        size_t Col0, Col1;
        struct SD_COL_BLOCK_DATA  *Next;
} SD_COL_BLOCK_DATA;

typedef union
{
        SD_ROW_BLOCK_DATA     UnusedData;
        struct
        {
                int                Row0;
                int                Row1;
        }  Any;
        struct
        {
                int                Row0;
                int                Row1;
                SD_COL_BLOCK_DATA *ColBlock;
        }  Data;
}  SD_TMP_ROW_BLOCK_DATA;

typedef struct
{
        size_t                 nRow;
        int                   *pPerm;

        int                    nRowBlock;
        SD_TMP_ROW_BLOCK_DATA *pRowBlock;
        int                    nColBlock;
        SD_CHUNK_DATA          PoolColBlock;
        SD_COL_BLOCK_DATA     *FreeColBlock;
}  SD_TMP_CON_MATRIX_DATA;

typedef enum StateType {ConMatrix, BlockConMatrix, BlockMatrix}  StateType;

typedef  struct
{
        int   nEq;
        int   nDeletedEq;

        StateType State;
        union
        {  SD_CON_MATRIX_DATA      Con;
        SD_TMP_CON_MATRIX_DATA  TmpCon;
        SD_BLOCK_MATRIX_DATA    Block;
        }  Mat;
}  SD_MATRIX_DATA;

typedef enum SD_MATRIX_WHAT
{
   /*
    * For symbolically factorizing and optimizing the structure of the [L] and [U] matrices
    * after the connectivity of the matrix [A] has been defined by calling
    * ds_DefineConnectivity (...)
    */
   SymbolicFactorize = 1,
   /*
    * For numerically factorizing the matrix [A] = [L][U] in to the lower resp. upper
    * triangular matrices [L] resp. [U] after [A] has been symbolically factorized by calling
    * ds_Solve(SymbolicFactorize,...) and assembled by calling ds_AssemblLocalMatrix(...) for
    * each element.
    */
   NumericFactorize  = 2,
   /*
    * For solving for a new right hand side load vector {B} after the matrix [A] has been
    * numerically factorized by calling ds_Solve(NumericFactorize,...)
    */
   BackForwardSubst  = 4,
   /*
    * For both numerically factorizing the matrix [A] and solving for a 1st right hand side
    * vector {B} i.e. is equivalent to the calls ds_Solve(NumericFactorize,...) &
    * ds_Solve(BackForwardSubst,...)
    */
   ComputeSolution   = NumericFactorize | BackForwardSubst,
   /*
    * For resetting all real coefficients of the matrix [A] to zero before reassembling a new
    * matrix [A] with identical connectivity by calling ds_AssembleLocalMatrix(...) with
    * changed local element matrices [ElMat] but the same element equations list Eq[...]
    */
   ResetMatrixData   = 8,
   /*
    * For releasing all storage space needed for solving the current problem defined by the
    * data pointed to by *pMat
    */
   ReleaseMatrixData = 16

} SD_MATRIX_WHAT;

//functions that are really used

int ds_DefineConnectivity( SD_MATRIX_DATA *const pMat0, const int& nEq, int Eq[], const int& nEl, const int& Dim );

int ds_Initialize( const size_t& MatDim, SD_MATRIX_DATA **ppMat );

int ds_AssembleMatrix( SD_MATRIX_DATA *pMat0, const int& nEq, int Eq[], const int& Dim, const double *ElMat );

int ds_Solve( const SD_MATRIX_WHAT& Code, SD_MATRIX_DATA *pMat, double *pX );

int ReleaseConMatrix( SD_CON_MATRIX_DATA * pMat );
int ReleaseBlockMatrix( SD_BLOCK_MATRIX_DATA * pMat );
#endif