Actual source code: test4.c
slepc-3.16.3 2022-04-11
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, 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: static char help[] = "Illustrates use of MFNSetBV().\n\n"
12: "The command line options are:\n"
13: " -t <sval>, where <sval> = scalar value that multiplies the argument.\n"
14: " -file <filename>, where <filename> = matrix file in PETSc binary form.\n\n";
16: #include <slepcmfn.h>
18: int main(int argc,char **argv)
19: {
20: Mat A; /* problem matrix */
21: MFN mfn;
22: FN f;
23: BV V;
24: PetscInt k=8;
25: PetscReal norm;
26: PetscScalar t=2.0;
27: Vec v,y;
28: PetscErrorCode ierr;
29: PetscViewer viewer;
30: PetscBool flg;
31: char filename[PETSC_MAX_PATH_LEN];
32: MFNConvergedReason reason;
34: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
36: PetscOptionsGetScalar(NULL,NULL,"-t",&t,NULL);
37: PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);
38: PetscPrintf(PETSC_COMM_WORLD,"\nMatrix exponential y=exp(t*A)*e, loaded from file\n\n");
40: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
41: Load matrix A from binary file
42: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
44: PetscOptionsGetString(NULL,NULL,"-file",filename,sizeof(filename),&flg);
45: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Must indicate a file name with the -file option");
47: #if defined(PETSC_USE_COMPLEX)
48: PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrix from a binary file...\n");
49: #else
50: PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrix from a binary file...\n");
51: #endif
52: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
53: MatCreate(PETSC_COMM_WORLD,&A);
54: MatSetFromOptions(A);
55: MatLoad(A,viewer);
56: PetscViewerDestroy(&viewer);
58: /* set v = ones(n,1) */
59: MatCreateVecs(A,NULL,&y);
60: MatCreateVecs(A,NULL,&v);
61: VecSet(v,1.0);
63: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
64: Create the BV object
65: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
67: BVCreate(PETSC_COMM_WORLD,&V);
68: PetscObjectSetName((PetscObject)V,"V");
69: BVSetSizesFromVec(V,v,k);
70: BVSetFromOptions(V);
72: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
73: Create the solver and set various options
74: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
76: MFNCreate(PETSC_COMM_WORLD,&mfn);
77: MFNSetOperator(mfn,A);
78: MFNSetBV(mfn,V);
79: MFNGetFN(mfn,&f);
80: FNSetType(f,FNEXP);
81: FNSetScale(f,t,1.0);
82: MFNSetDimensions(mfn,k);
83: MFNSetFromOptions(mfn);
85: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
86: Solve the problem, y=exp(t*A)*v
87: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
89: MFNSolve(mfn,v,y);
90: MFNGetConvergedReason(mfn,&reason);
91: if (reason<0) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Solver did not converge");
92: VecNorm(y,NORM_2,&norm);
93: PetscPrintf(PETSC_COMM_WORLD," Computed vector at time t=%.4g has norm %g\n\n",(double)PetscRealPart(t),(double)norm);
95: /*
96: Free work space
97: */
98: BVDestroy(&V);
99: MFNDestroy(&mfn);
100: MatDestroy(&A);
101: VecDestroy(&v);
102: VecDestroy(&y);
103: SlepcFinalize();
104: return ierr;
105: }
107: /*TEST
109: test:
110: suffix: 1
111: args: -file ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62b.petsc -k 12
112: requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
114: test:
115: suffix: 2
116: args: -file ${DATAFILESPATH}/matrices/complex/qc324.petsc -k 12
117: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
119: TEST*/