Actual source code: fnimpl.h
slepc-3.9.0 2018-04-12
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2018, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #if !defined(_FNIMPL)
12: #define _FNIMPL
14: #include <slepcfn.h>
15: #include <slepc/private/slepcimpl.h>
17: PETSC_EXTERN PetscBool FNRegisterAllCalled;
18: PETSC_EXTERN PetscErrorCode FNRegisterAll(void);
19: PETSC_EXTERN PetscLogEvent FN_Evaluate;
21: typedef struct _FNOps *FNOps;
23: struct _FNOps {
24: PetscErrorCode (*evaluatefunction)(FN,PetscScalar,PetscScalar*);
25: PetscErrorCode (*evaluatederivative)(FN,PetscScalar,PetscScalar*);
26: PetscErrorCode (*evaluatefunctionmat[FN_MAX_SOLVE])(FN,Mat,Mat);
27: PetscErrorCode (*evaluatefunctionmatvec[FN_MAX_SOLVE])(FN,Mat,Vec);
28: PetscErrorCode (*setfromoptions)(PetscOptionItems*,FN);
29: PetscErrorCode (*view)(FN,PetscViewer);
30: PetscErrorCode (*duplicate)(FN,MPI_Comm,FN*);
31: PetscErrorCode (*destroy)(FN);
32: };
34: #define FN_MAX_W 6
36: struct _p_FN {
37: PETSCHEADER(struct _FNOps);
38: /*------------------------- User parameters --------------------------*/
39: PetscScalar alpha; /* inner scaling (argument) */
40: PetscScalar beta; /* outer scaling (result) */
41: PetscInt method; /* the method to compute matrix functions */
42: FNParallelType pmode; /* parallel mode (redundant or synchronized) */
44: /*---------------------- Cached data and workspace -------------------*/
45: Mat W[FN_MAX_W]; /* workspace matrices */
46: PetscInt nw; /* number of allocated W matrices */
47: PetscInt cw; /* current W matrix */
48: void *data;
49: };
51: /*
52: FN_AllocateWorkMat - Allocate a work Mat of the same dimension of A and copy
53: its contents. The work matrix is returned in M and should be freed with
54: FN_FreeWorkMat().
55: */
56: PETSC_STATIC_INLINE PetscErrorCode FN_AllocateWorkMat(FN fn,Mat A,Mat *M)
57: {
59: PetscInt n,na;
60: PetscBool create=PETSC_FALSE;
63: *M = NULL;
64: if (fn->cw==FN_MAX_W) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Too many requested work matrices %D",fn->cw);
65: if (fn->nw<=fn->cw) {
66: create=PETSC_TRUE;
67: fn->nw++;
68: } else {
69: MatGetSize(fn->W[fn->cw],&n,NULL);
70: MatGetSize(A,&na,NULL);
71: if (n!=na) {
72: MatDestroy(&fn->W[fn->cw]);
73: create=PETSC_TRUE;
74: }
75: }
76: if (create) {
77: MatDuplicate(A,MAT_COPY_VALUES,&fn->W[fn->cw]);
78: PetscLogObjectParent((PetscObject)fn,(PetscObject)fn->W[fn->cw]);
79: } else {
80: MatCopy(A,fn->W[fn->cw],SAME_NONZERO_PATTERN);
81: }
82: *M = fn->W[fn->cw];
83: fn->cw++;
84: return(0);
85: }
87: /*
88: FN_FreeWorkMat - Release a work matrix created with FN_AllocateWorkMat().
89: */
90: PETSC_STATIC_INLINE PetscErrorCode FN_FreeWorkMat(FN fn,Mat *M)
91: {
93: if (!fn->cw) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"There are no work matrices");
94: fn->cw--;
95: if (fn->W[fn->cw]!=*M) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Work matrices must be freed in the reverse order of their creation");
96: *M = NULL;
97: return(0);
98: }
100: PETSC_INTERN PetscErrorCode SlepcMatDenseSqrt(PetscBLASInt,PetscScalar*,PetscBLASInt);
101: PETSC_INTERN PetscErrorCode SlepcSqrtmSchur(PetscBLASInt,PetscScalar*,PetscBLASInt,PetscBool);
102: PETSC_INTERN PetscErrorCode SlepcSqrtmDenmanBeavers(PetscBLASInt,PetscScalar*,PetscBLASInt,PetscBool);
103: PETSC_INTERN PetscErrorCode SlepcSqrtmNewtonSchulz(PetscBLASInt,PetscScalar*,PetscBLASInt,PetscBool);
104: PETSC_INTERN PetscErrorCode FNEvaluateFunctionMat_Private(FN,Mat,Mat,PetscBool);
105: PETSC_INTERN PetscErrorCode FNEvaluateFunctionMatVec_Private(FN,Mat,Vec,PetscBool);
107: #endif