AlgEqSystem.h

Go to the documentation of this file.

// $Id$
//==============================================================================
//==============================================================================
 
#ifndef _ALG_EQ_SYSTEM_H
#define _ALG_EQ_SYSTEM_H
 
#include "GlobalIntegral.h"
#include "LinAlgenums.h"
#include "MatVec.h"
 
class SAM;
class ProcessAdm;
class LinSolParams;
class SystemMatrix;
class SystemVector;
 
 
class AlgEqSystem : public GlobalIntegral
{
public:
  explicit AlgEqSystem(const SAM& s, const ProcessAdm* a = nullptr);
  AlgEqSystem(const AlgEqSystem& a) : sam(a.sam), adm(a.adm) { this->copy(a); }
  virtual ~AlgEqSystem() { this->clear(); }
 
  AlgEqSystem& operator=(const AlgEqSystem& s) { return this->copy(s); }
 
  AlgEqSystem& copy(const AlgEqSystem& that);
  AlgEqSystem& add(const AlgEqSystem& that);
 
  bool init(LinAlg::MatrixType mtype, const LinSolParams* spar = nullptr,
            size_t nmat = 1, size_t nvec = 1, size_t nscl = 0,
            bool reactions = false, int num_threads_SLU = 1,
            bool forcePreAss = false);
 
  void clear();
 
  bool setAssociatedVector(size_t imat, size_t ivec);
 
  virtual void initialize(char initLHS);
  virtual bool finalize(bool newLHS);
 
  virtual bool assemble(const LocalIntegral* elmObj, int elmId);
 
  size_t getNoRHS() const { return b.size(); }
  const SAM& getSAM() const { return sam; }
 
  SystemMatrix* getMatrix(size_t i = 0) { return i < A.size() ? A[i]._A : 0; }
  SystemVector* getVector(size_t i = 0) { return i < b.size() ? b[i] : 0; }
  double getScalar(size_t i = 0) { return i < c.size() ? c[i] : 0.0; }
  void addScalar(double s, size_t i = 0) { if (i < c.size()) c[i] += s; }
  const std::vector<double>& getScalars() const { return c; }
 
  const std::vector<double>* getReactions() const { return R.empty() ? 0 : &R; }
 
  bool staticCondensation(Matrix& Ared, Vector& bred,
                          const std::vector<int>& extNodes,
                          size_t imat = 0,
                          const char* recmatFile = nullptr) const;
 
  static bool readRecoveryMatrix(Matrix& Rmat, const char* recmatFile);
  bool recoverInternals(const Matrix& Rmat,
                        const std::vector<int>& extNodes,
                        const Vector& xe, Vector& xFull) const;
 
private:
  struct SysMatrixPair
  {
    SystemMatrix* _A; 
    SystemVector* _b; 
 
    SysMatrixPair() : _A(nullptr), _b(nullptr) {}
  };
 
  std::vector<SysMatrixPair> A; 
  std::vector<SystemVector*> b; 
  std::vector<double>        c; 
  std::vector<double>*       d; 
  std::vector<double>        R; 
 
  const SAM&        sam; 
  const ProcessAdm* adm; 
};
 
#endif