Skip to content
Snippets Groups Projects
Select Git revision
  • d4c7210267c15e0eda2e8ee431c9f9e35d2b8a3d
  • master default protected
  • dof-renumbering
  • gdemesy-master-patch-30528
  • eval-space-time
  • oscillating_multiharm
  • MH_movement
  • axisqu
  • write_vtu_and_ensight_formats
  • movingband
  • CP_1972_add_vtu_file_writing
  • mortar
  • fast_freq_sweep_Resolution
  • applyresolvent_again
  • marteaua-master-patch-54323
  • patch-1
  • binde-master-patch-08072
  • binde-master-patch-52461
  • BCGSL
  • resolvent
  • TreeElementsOf
  • getdp_3_5_0
  • getdp_3_4_0
  • getdp_3_3_0
  • getdp_3_2_0
  • getdp_3_1_0
  • getdp_3_0_4
  • getdp_3_0_3
  • getdp_3_0_2
  • getdp_3_0_1
  • getdp_3_0_0
  • onelab_mobile_2.1.0
  • getdp_2_11_3 protected
  • getdp_2_11_2 protected
  • getdp_2_11_1 protected
  • getdp_2_11_0 protected
  • getdp_2_10_0 protected
  • getdp_2_9_2 protected
  • getdp_2_9_1 protected
  • getdp_2_9_0 protected
  • getdp_2_8_0 protected
41 results

ProParser.h

Blame
  • linearSystemGMM.h 3.96 KiB
    // Gmsh - Copyright (C) 1997-2019 C. Geuzaine, J.-F. Remacle
    //
    // See the LICENSE.txt file for license information. Please report all
    // issues on https://gitlab.onelab.info/gmsh/gmsh/issues.
    
    #ifndef _LINEAR_SYSTEM_GMM_H_
    #define _LINEAR_SYSTEM_GMM_H_
    
    // Interface to GMM++
    
    #include "GmshConfig.h"
    #include "GmshMessage.h"
    #include "linearSystem.h"
    
    #if defined(HAVE_GMM)
    #include <gmm.h>
    
    template <class scalar> class linearSystemGmm : public linearSystem<scalar> {
    protected:
      std::vector<scalar>
        *_x; // the nonLinearSystemGmm has to access to this vector
      std::vector<scalar> *_b; // idem
      gmm::row_matrix<gmm::wsvector<scalar> > *_a;
    
    private:
      double _prec;
      int _noisy, _gmres;
    
    public:
      linearSystemGmm() : _x(0), _b(0), _a(0), _prec(1.e-8), _noisy(0), _gmres(0) {}
      virtual bool isAllocated() const { return _a != 0; }
      virtual void allocate(int nbRows)
      {
        clear();
        _a = new gmm::row_matrix<gmm::wsvector<scalar> >(nbRows, nbRows);
        _b = new std::vector<scalar>(nbRows);
        _x = new std::vector<scalar>(nbRows);
      }
      virtual ~linearSystemGmm() { clear(); }
      virtual void clear()
      {
        if(_a) {
          delete _a;
          delete _b;
          delete _x;
        }
        _a = 0;
      }
    
      virtual void addToMatrix(int row, int col, const scalar &val)
      {
        if(val != 0.0) (*_a)(row, col) += val;
      }
      virtual void getFromMatrix(int row, int col, scalar &val) const
      {
        val = (*_a)(row, col);
      }
    
      virtual void addToRightHandSide(int row, const scalar &val, int ith = 0)
      {
        if(val != 0.0) (*_b)[row] += val;
      }
    
      virtual void getFromRightHandSide(int row, scalar &val) const
      {
        val = (*_b)[row];
      }
    
      virtual void addToSolution(int row, const scalar &val)
      {
        if(val != 0.0) (*_x)[row] += val;
      }
    
      virtual void getFromSolution(int row, scalar &val) const { val = (*_x)[row]; }
      virtual void zeroMatrix() { gmm::clear(*_a); }
      virtual void zeroRightHandSide()
      {
        for(unsigned int i = 0; i < _b->size(); i++) (*_b)[i] = 0.;
      }
      virtual void zeroSolution()
      {
        for(unsigned int i = 0; i < _x->size(); i++) (*_x)[i] = 0.;
      }
    
      virtual double normInfRightHandSide() const
      {
        double nor = 0.;
        double temp;
        for(unsigned int i = 0; i < _b->size(); i++) {
          temp = abs((*_b)[i]); // this is valid also for complex
          // if(temp<0) temp = -temp;
          if(nor < temp) nor = temp;
        }
        return nor;
      }
    
      void setPrec(double p) { _prec = p; }
      void setNoisy(int n) { _noisy = n; }
      void setGmres(int n) { _gmres = n; }
      virtual int systemSolve()
      {
    #if defined(HAVE_MUMPS)
        gmm::MUMPS_solve(*_a, *_x, *_b);
    #else
        // gmm::ilutp_precond<gmm::row_matrix<gmm::wsvector<scalar> > > P(*_a, 25,
        // 0.);
        gmm::ildltt_precond<gmm::row_matrix<gmm::wsvector<scalar> > > P(*_a, 30,
                                                                        1.e-10);
        gmm::iteration iter(_prec);
        iter.set_noisy(_noisy);
        if(_gmres)
          gmm::gmres(*_a, *_x, *_b, P, 100, iter);
        else
          gmm::cg(*_a, *_x, *_b, P, iter);
    #endif
        return 1;
      }
    };
    
    #else
    
    template <class scalar> class linearSystemGmm : public linearSystem<scalar> {
    public:
      linearSystemGmm()
      {
        Msg::Error("Gmm++ is not available in this version of Gmsh");
      }
      virtual bool isAllocated() const { return false; }
      virtual void allocate(int nbRows) {}
      virtual void addToMatrix(int row, int col, const scalar &val) {}
      virtual void getFromMatrix(int row, int col, scalar &val) const {}
      virtual void addToRightHandSide(int row, const scalar &val, int ith = 0) {}
      virtual void getFromRightHandSide(int row, scalar &val) const {}
      virtual void addToSolution(int row, const scalar &val) {}
      virtual void getFromSolution(int row, scalar &val) const {}
      virtual void zeroMatrix() {}
      virtual void zeroRightHandSide() {}
      virtual void zeroSolution() {}
      virtual int systemSolve() { return 0; }
      virtual double normInfRightHandSide() const { return 0.; }
      void setPrec(double p) {}
      virtual void clear() {}
      void setNoisy(int n) {}
      void setGmres(int n) {}
    };
    
    #endif
    
    #endif