PETSc.hh

Go to the documentation of this file.

 
#ifndef op_PETSC_hh
#define op_PETSC_hh
 
extern "C" {
  // forwards for PETSc
  typedef struct _p_SLES* SLES;
  typedef struct _p_Mat* Mat;
  typedef struct _p_Vec* Vec;
}
#include <string>
 
#include "compositions.hh"
#include "basics/typedefs.hh"
#include "matrix.hh"
 
namespace concepts {
 
  // forward declarations
  template<class F>
  class Function;
 
  template<class F>
  class Vector;
 
  template<class F>
  class Space;
 
  template<class F, class G>
  class BilinearForm;
 
  // **************************************************************** PETSc **
 
  class PETSc : public Operator<Real> {
  public:
    PETSc(Operator<Real>& A, const Real maxeps, const std::string ksptype, 
          const std::string pctype = std::string("none"));
    ~PETSc();
    virtual void operator()(const Function<Real>& fncY, Function<Real>& fncX);
 
    uint iterations() const { return iterations_; }
  protected:
    virtual std::ostream& info(std::ostream& os) const;
  private:
    typedef int (*multPointer)(Mat, Vec, Vec);
    typedef void (*PETScmultPointer)();
 
    Operator<Real>& A_;
 
    int iterations_;
 
    Mat PETScA_;
    bool destroyPETScA_;
 
    SLES solver_;
 
    Vec x_, y_;
  };
 
  // ************************************************************* PETScMat **
 
  class PETScMat : public Matrix<Real> {
  public:
    PETScMat(const Space<Real>& spcX, const Space<Real>& spcY);
 
    PETScMat(const Space<Real>& spc, BilinearForm<Real,Real>& bf);
 
    PETScMat(const SparseMatrix<Real>& matrix);
 
    ~PETScMat();
 
    virtual void operator()(const Function<double>& fncY,
                            Function<double>& fncX);
    virtual void operator()(const Function<std::complex<double> >& fncY,
                            Function<std::complex<double> >& fncX);
    void operator()(const Vector<double>& fncY,
                    Vector<double>& fncX);
    virtual void transpMult(const Vector<double>& fncY,
                            Vector<double>& fncX);
    virtual void transpMult(const Vector<std::complex<double> >& fncY,
                            Vector<std::complex<double> >& fncX);
 
    virtual Real operator()(const uint i, const uint j) const;
    virtual Real& operator()(const uint i, const uint j);
 
    virtual const Space<double>& spaceX() const { return spcX_; }
    virtual const Space<double>& spaceY() const { return spcY_; }
 
    operator Mat() { return A_; }
 
    void storeMatlab(const char* name) const;
  protected:
    std::ostream& info(std::ostream& os) const;
  private:
    const Space<Real>& spcX_;
 
    uint nX_;
 
    const Space<Real>& spcY_;
 
    uint nY_;
 
    Vec x_, y_;
 
    Mat A_;
  };
 
} // namespace concepts
 
#endif // op_PETSC_hh