Actual source code: test1.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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[] = "Computes exp(t*A)*v for a matrix loaded from file.\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: PetscReal norm;
24: PetscScalar t=2.0;
25: Vec v,y;
26: PetscErrorCode ierr;
27: PetscViewer viewer;
28: PetscBool flg;
29: char filename[PETSC_MAX_PATH_LEN];
30: MFNConvergedReason reason;
32: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
34: PetscOptionsGetScalar(NULL,NULL,"-t",&t,NULL);
35: PetscPrintf(PETSC_COMM_WORLD,"\nMatrix exponential y=exp(t*A)*e, loaded from file\n\n");
37: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
38: Load matrix A from binary file
39: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
41: PetscOptionsGetString(NULL,NULL,"-file",filename,PETSC_MAX_PATH_LEN,&flg);
42: if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate a file name with the -file option");
44: #if defined(PETSC_USE_COMPLEX)
45: PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrix from a binary file...\n");
46: #else
47: PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrix from a binary file...\n");
48: #endif
49: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
50: MatCreate(PETSC_COMM_WORLD,&A);
51: MatSetFromOptions(A);
52: MatLoad(A,viewer);
53: PetscViewerDestroy(&viewer);
55: /* set v = ones(n,1) */
56: MatCreateVecs(A,NULL,&y);
57: MatCreateVecs(A,NULL,&v);
58: VecSet(v,1.0);
60: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
61: Create the solver and set various options
62: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
64: MFNCreate(PETSC_COMM_WORLD,&mfn);
65: MFNSetOperator(mfn,A);
66: MFNGetFN(mfn,&f);
67: FNSetType(f,FNEXP);
68: FNSetScale(f,t,1.0);
69: MFNSetFromOptions(mfn);
71: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: Solve the problem, y=exp(t*A)*v
73: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
75: MFNSolve(mfn,v,y);
76: MFNGetConvergedReason(mfn,&reason);
77: if (reason<0) SETERRQ(PETSC_COMM_WORLD,1,"Solver did not converge");
78: VecNorm(y,NORM_2,&norm);
79: PetscPrintf(PETSC_COMM_WORLD," Computed vector at time t=%.4g has norm %g\n\n",(double)PetscRealPart(t),(double)norm);
81: /*
82: Free work space
83: */
84: MFNDestroy(&mfn);
85: MatDestroy(&A);
86: VecDestroy(&v);
87: VecDestroy(&y);
88: SlepcFinalize();
89: return ierr;
90: }
92: /*TEST
94: test:
95: suffix: 1
96: args: -file ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62b.petsc -mfn_type {{krylov expokit}}
97: requires: double !complex !define(PETSC_USE_64BIT_INDICES)
98: output_file: output/test1.out
100: test:
101: suffix: 2
102: args: -file ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62b.petsc -mfn_type {{krylov expokit}} -mat_type aijcusparse
103: requires: cuda double !complex !define(PETSC_USE_64BIT_INDICES)
104: output_file: output/test1.out
106: TEST*/