Actual source code: lmbasis.h

  1: #pragma once
  2: #include <petscmat.h>

  4: PETSC_INTERN PetscLogEvent LMBASIS_GEMM;
  5: PETSC_INTERN PetscLogEvent LMBASIS_GEMV;
  6: PETSC_INTERN PetscLogEvent LMBASIS_GEMVH;

  8: typedef struct _n_VecLink *VecLink;

 10: struct _n_VecLink {
 11:   Vec     vec;
 12:   VecLink next;
 13: };

 15: // Limited Memory Basis
 16: typedef struct _n_LMBasis *LMBasis;
 17: struct _n_LMBasis {
 18:   PetscInt         m;                   // Number of vectors in the limited memory window
 19:   PetscInt         k;                   // Index of the history-order next vector to be inserted
 20:   Mat              vecs;                // Dense matrix backing storage of vectors
 21:   PetscObjectId    operator_id;         // If these vecs include the output of an operator (like B0 * S), the id of the operator
 22:   PetscObjectState operator_state;      // The state of the operator when vectors in S were computed, to determine when basis vectors are stale because B0 has changed
 23:   Vec              cached_product;      // Some methods will cache v <- S^T g during MatLMVMUpdate(B, x, g) for use in MatSolve(B, g, p), this is that v
 24:   PetscObjectId    cached_vec_id;       // The id of g, to help determine when an input vector is g that was used to compute v
 25:   PetscObjectState cached_vec_state;    // The state of g when v was computed, to determine when v is stale because g has changed
 26:   VecLink          work_vecs_available; // work vectors the layout of column vectors of S
 27:   VecLink          work_vecs_in_use;
 28:   VecLink          work_rows_available; // work vectors the layout of row vectors of S
 29:   VecLink          work_rows_in_use;
 30: };

 32: PETSC_INTERN PetscErrorCode LMBasisCreate(Vec, PetscInt, LMBasis *);
 33: PETSC_INTERN PetscErrorCode LMBasisDestroy(LMBasis *);
 34: PETSC_INTERN PetscErrorCode LMBasisReset(LMBasis);
 35: PETSC_INTERN PetscErrorCode LMBasisGetNextVec(LMBasis, Vec *);
 36: PETSC_INTERN PetscErrorCode LMBasisRestoreNextVec(LMBasis, Vec *);
 37: PETSC_INTERN PetscErrorCode LMBasisSetNextVec(LMBasis, Vec);
 38: PETSC_INTERN PetscErrorCode LMBasisGetVecRead(LMBasis, PetscInt, Vec *);
 39: PETSC_INTERN PetscErrorCode LMBasisRestoreVecRead(LMBasis, PetscInt, Vec *);
 40: PETSC_INTERN PetscErrorCode LMBasisCopy(LMBasis, LMBasis);
 41: PETSC_INTERN PetscErrorCode LMBasisCreateRow(LMBasis, Vec *);
 42: PETSC_INTERN PetscErrorCode LMBasisGetWorkRow(LMBasis, Vec *);
 43: PETSC_INTERN PetscErrorCode LMBasisRestoreWorkRow(LMBasis, Vec *);
 44: PETSC_INTERN PetscErrorCode LMBasisGetWorkVec(LMBasis, Vec *);
 45: PETSC_INTERN PetscErrorCode LMBasisRestoreWorkVec(LMBasis, Vec *);
 46: PETSC_INTERN PetscErrorCode LMBasisGetRange(LMBasis, PetscInt *, PetscInt *);
 47: PETSC_INTERN PetscErrorCode LMBasisGEMV(LMBasis, PetscInt, PetscInt, PetscScalar, Vec, PetscScalar, Vec);
 48: PETSC_INTERN PetscErrorCode LMBasisGEMVH(LMBasis, PetscInt, PetscInt, PetscScalar, Vec, PetscScalar, Vec);
 49: PETSC_INTERN PetscErrorCode LMBasisGEMMH(LMBasis, PetscInt, PetscInt, LMBasis, PetscInt, PetscInt, PetscScalar, PetscScalar, Mat);
 50: PETSC_INTERN PetscErrorCode LMBasisSetCachedProduct(LMBasis, Vec, Vec);