sss.h

高效的c++科学算法库

/*                                                           -*- c++ -*-
 ***********************************************************************
 * Scientific Library (GNU Public Licence)
 *
 * module Sparse Symmetrix Skyline Matrix
 *
 * Author: Pierre Aubert [e-mail]
 *
 * $Id: sss.h,v 1.3 1998/11/16 17:22:13 paubert Exp $
 *
 * Suggestions: [e-mail]
 * Bugs:   [e-mail] 
 *
 * For more information, please see the sl++ Home Page:
 * http://wwwinfo.cern.ch/~ldeniau/sl.html 
 *
 ***********************************************************************
 */

#ifndef SL_MATRIX_STRUCTURES_SSS_SSS_H
#define SL_MATRIX_STRUCTURES_SSS_SSS_H

#ifndef SL_MATRIX_STRUCTURES_SSS_H
#error <sl/matrix/structures/sss/sss.h> must be included via <sl/matrix/structures/sss.h>
#endif

#ifdef HAVE_NAMESPACE
namespace sl 
{
#endif

    class SymmetricSkyline : public MatrixStructure
    {

          /** @name Explanation
              
              This sparse storage is based on the following idea:
              
              $\dots$
          */
        
    public:
        
        typedef SymmetricSkylineIterator iterator_t;

          /** @name Constructor(s), Destructors(s)
           */
          //@{
          /// default constructor
        SymmetricSkyline ();
          /// build a $n\times n$ matrix
        explicit SymmetricSkyline (Index const& n);
          /// build a $n\times p$ matrix, with $n=p$.
        SymmetricSkyline (Index const& n, Index const& p);
          /** build a $n$ by $n$ matrix and store at least nb\_store datas.
          */
        SymmetricSkyline (Index const& n, Index const& p, UInt const& nb_store);
          /// destructor
        ~SymmetricSkyline();
          // WARNING: Do NOT declare the copy constructor !
          //@}

          /** @name operations on profil
           */
          //@{
          /// get size
        size_t size () const;
          /** resize: almost impossible
          */
        void resize (Index const n);
          /** resize: almost impossible
          */
        void resize (Index const n, Index const p);
          /// if you already have computed the profil
        void compute_profil( VectorUInt const& haut );
          /// from finite element method
        void compute_profil( MatrixUInt const& noda,
                             MatrixUInt const& connex );
          /// display informations on profil
        ostream& display( ostream& os, bool const verbose=false ) const;
          //@}
    
          /** @name Access to element
           */
          //@{
          /// get number of rows
        Index rows() const;
          /// get number of cols
        Index cols() const;
          /// get number of rows
        Index n() const {return my_n;};
          /// get number of cols
        Index p() const {return my_n;};
          /// get first element in the i-th row
        Index firstInRow (Index const i) const;
          /// get first element in the j-th column
        Index firstInCol (Index const j) const;
          /// get last element in the i-th row
        Index lastInRow  (Index const i) const;
          /// get last element in the i-th column
        Index lastInCol  (Index const j) const;
          /// true if the storage is row major (C-like)
        static Bool isRowMajor () {return false;};
          /// true if the storage is column major (fortran-like)
        static Bool isColumnMajor () {return true;};
          /// get value 
        template <typename T_value>
        T_value& get (T_value* const data, Index const i, Index const j);
          /// get value 
        template <typename T_value>
        T_value get (T_value const* const data, Index const i, Index const j) const;
          /// get value fast
        template <typename T_value>
        T_value& get_fast (T_value* const data, Index const i, Index const j);
          /// get value fast
        template <typename T_value>
        T_value get_fast (T_value const* const data, Index const i, Index const j) const;
          /// get value of haut (should not be public)
        UInt haut(Index const& i) const;
          /// get value of diag (should not be public)
        UInt diag(Index const& i) const;
          /// get haut (should not be public)
        VectorUInt const& haut() const {return my_haut;};
          ///get diag (should not be public)
        VectorUInt const& diag() const {return my_diag;};
          /// true if matrix is factored
        Bool isFactored() const;
          /// set isFactored to true
        void setIsFactored();
          /// set isFactored to false
        void setIsNotFactored();
          //@}
    
    protected:
    
          // map(MinIndex,MinIndex)
        ptrdiff_t offset () const { return MinIndex; } // ???

    private:
          // number of rows or columns
        UInt                   my_n;
          /// number of elements asked by user
        size_t                 my_data_length;
          /// is profil set ?
        Bool                   my_profil_is_set;
          /// is matrix factored ?
        Bool                   my_matrix_is_factored;
          /** map(i,j) find the position in the data array
          */
        Index                  map (Index const& i, Index const& j) const;
          /** same as map(i,j) but without any control
          */
        inline Index           map_fast (Index const& i, Index const& j) const;
          /// diagonal index
        VectorUInt             my_diag;
          /// haut index 
        VectorUInt             my_haut;
          // copy constructor is private to prevent unwanted allocation
          // MUST BE supress
          // SymmetricSkyline( SymmetricSkyline const& copy ) {};
    };
  
  
      // ----------------------------------------------------------------------
      //                                                                  get()
      // ----------------------------------------------------------------------
    template <typename T_value>
    T_value&
    SymmetricSkyline::get (T_value* const data, Index const i, Index const j) {
        Index pos = map(i,j);
        if( pos == 0 )
            {
                cerr << "Pierre: debug: hors profil in get pour i=" << i << " j=" << j << endl;
                return zero_ref<T_value>();
            }
        else
            return data[pos];
    }
  
      // ----------------------------------------------------------------------
      //                                                                  get()
      // ----------------------------------------------------------------------
    template <typename T_value>
    T_value
    SymmetricSkyline::get (T_value const* const data, Index const i, Index const j) const {
        Index pos=map(i,j);
        if( pos == 0 )
            return T_value();
        else
            return data[pos];
    }
  
      // ----------------------------------------------------------------------
      //                                                             get_fast()
      // ----------------------------------------------------------------------
    template <typename T_value>
    T_value&
    SymmetricSkyline::get_fast (T_value* const data, Index const i, Index const j) {
        return data[map_fast(i,j)];
    }
  
      // ----------------------------------------------------------------------
      //                                                             get_fast()
      // ----------------------------------------------------------------------
    template <typename T_value>
    T_value
    SymmetricSkyline::get_fast (T_value const* const data, Index const i, Index const j) const{
        return data[map_fast(i,j)];
    }

#ifdef HAVE_NAMESPACE
}
#endif

#endif
// end of SL_MATRIX_STRUCTURES_SSS_SSS_H definition(s)
// do not write anything after this line!