Actual source code: lmvmmat.h

  1: #ifndef __LMVMMAT_H


  5: #include "private/matimpl.h"
  6: #include "tao_sys.h"

  8: #define MatLMVM_Scale_None                0
  9: #define MatLMVM_Scale_Scalar                1
 10: #define MatLMVM_Scale_Broyden                2
 11: #define MatLMVM_Scale_Types             3

 13: #define MatLMVM_Rescale_None                0
 14: #define MatLMVM_Rescale_Scalar                1
 15: #define MatLMVM_Rescale_GL                2
 16: #define MatLMVM_Rescale_Types                  3

 18: #define MatLMVM_Limit_None                0
 19: #define MatLMVM_Limit_Average                1
 20: #define MatLMVM_Limit_Relative                2
 21: #define MatLMVM_Limit_Absolute                3
 22: #define MatLMVM_Limit_Types                4

 24: #define TAO_ZERO_SAFEGUARD        1e-8
 25: #define TAO_INF_SAFEGUARD        1e+8

 27: typedef struct{
 28:     PetscBool allocated;
 29:     PetscInt lm;
 30:     PetscReal eps;
 31:     PetscInt limitType;
 32:     PetscInt scaleType;
 33:     PetscInt rScaleType;
 34:     
 35:     PetscReal s_alpha;        /*  Factor for scalar scaling */
 36:     PetscReal r_alpha;        /*  Factor on scalar for rescaling diagonal matrix */
 37:     PetscReal r_beta;        /*  Factor on diagonal for rescaling diagonal matrix */
 38:     PetscReal mu;                /*  Factor for using historical information */
 39:     PetscReal nu;                /*  Factor for using historical information */
 40:     PetscReal phi;                /*  Factor for Broyden scaling */

 42:   PetscInt scalar_history;        /*  Amount of history to keep for scalar scaling */
 43:   PetscReal *yy_history;        /*  Past information for scalar scaling */
 44:   PetscReal *ys_history;        /*  Past information for scalar scaling */
 45:   PetscReal *ss_history;        /*  Past information for scalar scaling */

 47:   PetscInt rescale_history;  /*  Amount of history to keep for rescaling diagonal */
 48:   PetscReal *yy_rhistory;        /*  Past information for scalar rescaling */
 49:   PetscReal *ys_rhistory;        /*  Past information for scalar rescaling */
 50:   PetscReal *ss_rhistory;        /*  Past information for scalar rescaling */

 52:   PetscReal delta_max;        /*  Maximum value for delta */
 53:   PetscReal delta_min;        /*  Minimum value for delta */

 55:   PetscInt lmnow;
 56:   PetscInt iter;
 57:   PetscInt nupdates;
 58:   PetscInt nrejects;

 60:   Vec *S;
 61:   Vec *Y;
 62:   Vec Gprev;
 63:   Vec Xprev;

 65:   Vec D;
 66:   Vec U;
 67:   Vec V;
 68:   Vec W;
 69:   Vec P;
 70:   Vec Q;

 72:   PetscReal delta;
 73:   PetscReal sigma;
 74:   PetscReal *rho;
 75:   PetscReal *beta;

 77:   PetscBool useDefaultH0;
 78:   Mat H0;

 80:   PetscBool useScale;
 81:   Vec scale;
 82:      

 84: } MatLMVMCtx;




 90: /* PETSc Mat overrides */

 94: /*
 95: int MatMultTranspose_LMVM(Mat,Vec,Vec);
 96: int MatDiagonalShift_LMVM(Vec,Mat);
 97: int MatDestroy_LMVM(Mat);
 98: int MatShift_LMVM(Mat,PetscReal);
 99: int MatDuplicate_LMVM(Mat,MatDuplicateOption,Mat*);
100: int MatEqual_LMVM(Mat,Mat,PetscBool*);
101: int MatScale_LMVM(Mat,PetscReal);
102: int MatGetSubMatrix_LMVM(Mat,IS,IS,int,MatReuse,Mat *);
103: int MatGetSubMatrices_LMVM(Mat,int,IS*,IS*,MatReuse,Mat**);
104: int MatTranspose_LMVM(Mat,Mat*);
105: int MatGetDiagonal_LMVM(Mat,Vec);
106: int MatGetColumnVector_LMVM(Mat,Vec, int);
107: int MatNorm_LMVM(Mat,NormType,PetscReal *);
108: */

110: /* Functions used by TAO */
111: PetscErrorCode MatLMVMReset(Mat);
112: PetscErrorCode MatLMVMUpdate(Mat,Vec, Vec);
113: PetscErrorCode MatLMVMSetDelta(Mat,PetscReal);
114: PetscErrorCode MatLMVMSetScale(Mat,Vec);
115: PetscErrorCode MatLMVMGetRejects(Mat,PetscInt*);
116: PetscErrorCode MatLMVMSetH0(Mat,Mat);
117: PetscErrorCode MatLMVMSetPrev(Mat,Vec,Vec);
118: PetscErrorCode MatLMVMGetX0(Mat,Vec);
119: PetscErrorCode MatLMVMRefine(Mat, Mat, Vec, Vec);
120: PetscErrorCode MatLMVMAllocateVectors(Mat m, Vec v);
121: PetscErrorCode MatLMVMSolve(Mat, Vec, Vec);


124: #endif