conjugategrad.h

Go to the documentation of this file.

/******************************************************************************
 *
 * MantaFlow fluid solver framework
 * Copyright 2011 Tobias Pfaff, Nils Thuerey 
 *
 * This program is free software, distributed under the terms of the
 * GNU General Public License (GPL) 
 * http://www.gnu.org/licenses
 *
 * Conjugate gradient solver
 *
 ******************************************************************************/

#ifndef _CONJUGATEGRADIENT_H
#define _CONJUGATEGRADIENT_H

#include "vectorbase.h"
#include "grid.h"
#include "kernel.h"
#include "multigrid.h"

namespace Manta { 

static const bool CG_DEBUG = false;

class GridCgInterface {
    public:
        enum PreconditionType { PC_None=0, PC_ICP, PC_mICP, PC_MGP };
        
        GridCgInterface() : mUseL2Norm(true) {};
        virtual ~GridCgInterface() {};

        // solving functions
        virtual bool iterate() = 0;
        virtual void solve(int maxIter) = 0;

        // precond
        virtual void setICPreconditioner(PreconditionType method, Grid<Real> *A0, Grid<Real> *Ai, Grid<Real> *Aj, Grid<Real> *Ak) = 0;
        virtual void setMGPreconditioner(PreconditionType method, GridMg* MG) = 0;

        // access
        virtual Real getSigma() const = 0;
        virtual Real getIterations() const = 0;
        virtual Real getResNorm() const = 0;
        virtual void setAccuracy(Real set) = 0;
        virtual Real getAccuracy() const = 0;

        virtual void forceReinit() = 0;

        void setUseL2Norm(bool set) { mUseL2Norm = set; }

    protected:

        // use l2 norm of residualfor threshold? (otherwise uses max norm)
        bool mUseL2Norm; 
};



template<class APPLYMAT>
class GridCg : public GridCgInterface {
    public:
        GridCg(Grid<Real>& dst, Grid<Real>& rhs, Grid<Real>& residual, Grid<Real>& search, FlagGrid& flags, Grid<Real>& tmp, 
                Grid<Real>* A0, Grid<Real>* pAi, Grid<Real>* pAj, Grid<Real>* pAk);
        ~GridCg() {}
        
        void doInit();
        bool iterate();
        void solve(int maxIter);
        void setICPreconditioner(PreconditionType method, Grid<Real> *A0, Grid<Real> *Ai, Grid<Real> *Aj, Grid<Real> *Ak);
        void setMGPreconditioner(PreconditionType method, GridMg* MG);
        void forceReinit() { mInited = false; }
        
        // Accessors        
        Real getSigma() const { return mSigma; }
        Real getIterations() const { return mIterations; }

        Real getResNorm() const { return mResNorm; }

        void setAccuracy(Real set) { mAccuracy=set; }
        Real getAccuracy() const { return mAccuracy; }

    protected:
        bool mInited;
        int mIterations;
        // grids
        Grid<Real>& mDst;
        Grid<Real>& mRhs;
        Grid<Real>& mResidual;
        Grid<Real>& mSearch;
        FlagGrid& mFlags;
        Grid<Real>& mTmp;

        Grid<Real> *mpA0, *mpAi, *mpAj, *mpAk;

        PreconditionType mPcMethod;
        Grid<Real> *mpPCA0, *mpPCAi, *mpPCAj, *mpPCAk;
        GridMg* mMG;

        Real mSigma;
        Real mAccuracy;
        Real mResNorm;
}; // GridCg


void ApplyMatrix (FlagGrid& flags, Grid<Real>& dst, Grid<Real>& src, Grid<Real>& A0, Grid<Real>& Ai, Grid<Real>& Aj, Grid<Real>& Ak) {}


void ApplyMatrix2D (FlagGrid& flags, Grid<Real>& dst, Grid<Real>& src, Grid<Real>& A0, Grid<Real>& Ai, Grid<Real>& Aj, Grid<Real>& Ak) {}


void MakeLaplaceMatrix(FlagGrid& flags, Grid<Real>& A0, Grid<Real>& Ai, Grid<Real>& Aj, Grid<Real>& Ak, MACGrid* fractions = 0) {}





} // namespace

#endif