Actual source code: test25.c
slepc-3.17.2 2022-08-09
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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[] = "Solves a GNHEP problem with contour integral. "
12: "Based on ex7.\n"
13: "The command line options are:\n"
14: " -f1 <filename> -f2 <filename>, PETSc binary files containing A and B.\n\n";
16: #include <slepceps.h>
18: int main(int argc,char **argv)
19: {
20: EPS eps;
21: RG rg;
22: Mat A,B;
23: PetscBool flg;
24: EPSCISSExtraction extr;
25: EPSCISSQuadRule quad;
26: char filename[PETSC_MAX_PATH_LEN];
27: PetscViewer viewer;
29: SlepcInitialize(&argc,&argv,(char*)0,help);
30: PetscPrintf(PETSC_COMM_WORLD,"\nGNHEP problem with contour integral\n\n");
32: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
33: Load the matrices that define the eigensystem, Ax=kBx
34: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
36: PetscOptionsGetString(NULL,NULL,"-f1",filename,sizeof(filename),&flg);
39: #if defined(PETSC_USE_COMPLEX)
40: PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrices from binary files...\n");
41: #else
42: PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrices from binary files...\n");
43: #endif
44: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
45: MatCreate(PETSC_COMM_WORLD,&A);
46: MatSetFromOptions(A);
47: MatLoad(A,viewer);
48: PetscViewerDestroy(&viewer);
50: PetscOptionsGetString(NULL,NULL,"-f2",filename,sizeof(filename),&flg);
51: if (flg) {
52: PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
53: MatCreate(PETSC_COMM_WORLD,&B);
54: MatSetFromOptions(B);
55: MatLoad(B,viewer);
56: PetscViewerDestroy(&viewer);
57: } else {
58: PetscPrintf(PETSC_COMM_WORLD," Matrix B was not provided, setting B=I\n\n");
59: B = NULL;
60: }
62: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
63: Create the eigensolver and solve the problem
64: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
66: EPSCreate(PETSC_COMM_WORLD,&eps);
67: EPSSetOperators(eps,A,B);
68: EPSSetProblemType(eps,EPS_GNHEP);
69: EPSSetTolerances(eps,1e-9,PETSC_DEFAULT);
71: /* set CISS solver with various options */
72: EPSSetType(eps,EPSCISS);
73: EPSCISSSetExtraction(eps,EPS_CISS_EXTRACTION_HANKEL);
74: EPSCISSSetQuadRule(eps,EPS_CISS_QUADRULE_CHEBYSHEV);
75: EPSCISSSetUseST(eps,PETSC_TRUE);
76: EPSGetRG(eps,&rg);
77: RGSetType(rg,RGINTERVAL);
78: RGIntervalSetEndpoints(rg,-3000.0,0.0,0.0,0.0);
80: EPSSetFromOptions(eps);
82: EPSSolve(eps);
83: PetscObjectTypeCompare((PetscObject)eps,EPSCISS,&flg);
84: if (flg) {
85: EPSCISSGetExtraction(eps,&extr);
86: EPSCISSGetQuadRule(eps,&quad);
87: PetscPrintf(PETSC_COMM_WORLD," Solved with CISS using %s extraction with %s quadrature rule\n\n",EPSCISSExtractions[extr],EPSCISSQuadRules[quad]);
88: }
90: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
91: Display solution and clean up
92: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
94: EPSErrorView(eps,EPS_ERROR_BACKWARD,NULL);
95: EPSDestroy(&eps);
96: MatDestroy(&A);
97: MatDestroy(&B);
98: SlepcFinalize();
99: return 0;
100: }
102: /*TEST
104: testset:
105: args: -f1 ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62a.petsc -f2 ${SLEPC_DIR}/share/slepc/datafiles/matrices/bfw62b.petsc
106: output_file: output/test25_1.out
107: test:
108: suffix: 1
109: requires: double !complex !defined(PETSC_USE_64BIT_INDICES)
110: test:
111: suffix: 1_cuda
112: args: -mat_type aijcusparse -st_pc_factor_mat_solver_type cusparse
113: requires: cuda double !complex !defined(PETSC_USE_64BIT_INDICES)
115: testset:
116: args: -f1 ${DATAFILESPATH}/matrices/complex/qc324.petsc -rg_type ellipse -rg_ellipse_center 1-0.09i -rg_ellipse_radius 0.15 -rg_ellipse_vscale 0.1
117: output_file: output/test25_2.out
118: test:
119: suffix: 2
120: requires: double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
121: test:
122: suffix: 2_cuda
123: args: -mat_type aijcusparse -st_pc_factor_mat_solver_type cusparse
124: requires: cuda double complex datafilespath !defined(PETSC_USE_64BIT_INDICES)
126: TEST*/