1 /* file: del2mat.h */
2
3 #if !defined(DEL2MAT_H)
4 #define DEL2MAT_H
5
6 #include <petsc.h>
7 #include <petscvec.h>
8 #include <petscmat.h>
9
10 /* external Fortran 90 subroutine */
11 #define Del2Apply del2apply_
12 EXTERN_C_BEGIN
13 extern void Del2Apply(int*,double*,const double*,double*);
14 EXTERN_C_END
15
16 /* user data structure and routines
17 * defining the matrix-free operator */
18
19 typedef struct {
20 PetscInt N;
21 PetscScalar *F;
22 } Del2Mat;
23
24 /* y <- A * x */
Del2Mat_mult(Mat A,Vec x,Vec y)25 PetscErrorCode Del2Mat_mult(Mat A, Vec x, Vec y)
26 {
27 Del2Mat *ctx;
28 const PetscScalar *xx;
29 PetscScalar *yy;
30
31 PetscFunctionBegin;
32 PetscCall(MatShellGetContext(A,&ctx));
33 /* get raw vector arrays */
34 PetscCall(VecGetArrayRead(x,&xx));
35 PetscCall(VecGetArray(y,&yy));
36 /* call external Fortran subroutine */
37 Del2Apply(&ctx->N,ctx->F,xx,yy);
38 /* restore raw vector arrays */
39 PetscCall(VecRestoreArrayRead(x,&xx));
40 PetscCall(VecRestoreArray(y,&yy));
41 PetscFunctionReturn(PETSC_SUCCESS);
42 }
43
44 /*D_i <- A_ii */
Del2Mat_diag(Mat A,Vec D)45 PetscErrorCode Del2Mat_diag(Mat A, Vec D)
46 {
47 PetscFunctionBegin;
48 PetscCall(VecSet(D,6.0));
49 PetscFunctionReturn(PETSC_SUCCESS);
50 }
51
52 #endif
53