Actual source code: test28.c
slepc-3.18.3 2023-03-24
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[] = "Tests multiple calls to EPSSolve with different matrix of different size.\n\n"
12: "The command line options are:\n"
13: " -n <n>, where <n> = number of grid subdivisions in x dimension.\n"
14: " -m <m>, where <m> = number of grid subdivisions in y dimension.\n\n";
16: #include <slepceps.h>
18: int main(int argc,char **argv)
19: {
20: Mat A,B;
21: EPS eps;
22: PetscInt N,n=10,m=11,Istart,Iend,II,nev=3,i,j;
23: PetscBool flag,terse;
26: SlepcInitialize(&argc,&argv,(char*)0,help);
27: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
28: PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag);
29: N = n*m;
30: PetscPrintf(PETSC_COMM_WORLD,"\n2-D Laplacian Eigenproblem, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,m);
32: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
33: Create the 2-D Laplacian with coarse mesh
34: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
36: MatCreate(PETSC_COMM_WORLD,&A);
37: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N);
38: MatSetFromOptions(A);
39: MatSetUp(A);
40: MatGetOwnershipRange(A,&Istart,&Iend);
41: for (II=Istart;II<Iend;II++) {
42: i = II/n; j = II-i*n;
43: if (i>0) MatSetValue(A,II,II-n,-1.0,INSERT_VALUES);
44: if (i<m-1) MatSetValue(A,II,II+n,-1.0,INSERT_VALUES);
45: if (j>0) MatSetValue(A,II,II-1,-1.0,INSERT_VALUES);
46: if (j<n-1) MatSetValue(A,II,II+1,-1.0,INSERT_VALUES);
47: MatSetValue(A,II,II,4.0,INSERT_VALUES);
48: }
49: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
50: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
52: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
53: Create the eigensolver, set options and solve the eigensystem
54: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
56: EPSCreate(PETSC_COMM_WORLD,&eps);
57: EPSSetOperators(eps,A,NULL);
58: EPSSetProblemType(eps,EPS_HEP);
59: EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);
60: EPSSetDimensions(eps,nev,PETSC_DEFAULT,PETSC_DEFAULT);
61: EPSSetFromOptions(eps);
63: EPSSolve(eps);
65: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
66: Display solution of first solve
67: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
69: PetscOptionsHasName(NULL,NULL,"-terse",&terse);
70: if (terse) EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
71: else {
72: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL);
73: EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD);
74: EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD);
75: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
76: }
78: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
79: Create the 2-D Laplacian with finer mesh
80: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
82: n *= 2;
83: m *= 2;
84: N = n*m;
85: PetscPrintf(PETSC_COMM_WORLD,"\n2-D Laplacian Eigenproblem, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,m);
87: MatCreate(PETSC_COMM_WORLD,&B);
88: MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,N,N);
89: MatSetFromOptions(B);
90: MatSetUp(B);
91: MatGetOwnershipRange(B,&Istart,&Iend);
92: for (II=Istart;II<Iend;II++) {
93: i = II/n; j = II-i*n;
94: if (i>0) MatSetValue(B,II,II-n,-1.0,INSERT_VALUES);
95: if (i<m-1) MatSetValue(B,II,II+n,-1.0,INSERT_VALUES);
96: if (j>0) MatSetValue(B,II,II-1,-1.0,INSERT_VALUES);
97: if (j<n-1) MatSetValue(B,II,II+1,-1.0,INSERT_VALUES);
98: MatSetValue(B,II,II,4.0,INSERT_VALUES);
99: }
100: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
101: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
103: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
104: Solve again, calling EPSReset() since matrix size has changed
105: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
107: /* EPSReset(eps); */ /* not required, will be called in EPSSetOperators() */
108: EPSSetOperators(eps,B,NULL);
109: EPSSolve(eps);
111: if (terse) EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
112: else {
113: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL);
114: EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD);
115: EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD);
116: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
117: }
119: EPSDestroy(&eps);
120: MatDestroy(&A);
121: MatDestroy(&B);
122: SlepcFinalize();
123: return 0;
124: }
126: /*TEST
128: testset:
129: requires: !single
130: output_file: output/test28_1.out
131: test:
132: suffix: 1
133: args: -eps_type {{krylovschur arnoldi gd jd rqcg lobpcg lapack}} -terse
134: test:
135: suffix: 1_lanczos
136: args: -eps_type lanczos -eps_lanczos_reorthog local -terse
138: test:
139: suffix: 2
140: args: -eps_type {{power subspace}} -eps_target 8 -st_type sinvert -terse
142: testset:
143: args: -eps_interval 0.5,0.67 -terse
144: output_file: output/test28_3.out
145: test:
146: suffix: 3
147: args: -st_type sinvert -st_pc_type cholesky
148: requires: !single
149: test:
150: suffix: 3_evsl
151: nsize: {{1 2}}
152: args: -eps_type evsl -eps_evsl_dos_method lanczos
153: requires: evsl
155: TEST*/