Actual source code: test39.c
slepc-3.16.2 2022-02-01
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[] = "Tests multiple calls to EPSSolve with matrices of different local 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: /*
19: Create 2-D Laplacian matrix
20: */
21: PetscErrorCode Laplacian(MPI_Comm comm,PetscInt n,PetscInt m,PetscInt shift,Mat *A)
22: {
24: PetscInt N = n*m,i,j,II,Istart,Iend,nloc;
25: PetscMPIInt rank;
28: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
29: nloc = PETSC_DECIDE;
30: PetscSplitOwnership(comm,&nloc,&N);
31: if (rank==0) nloc += shift;
32: else if (rank==1) nloc -= shift;
34: MatCreate(comm,A);
35: MatSetSizes(*A,nloc,nloc,N,N);
36: MatSetFromOptions(*A);
37: MatSetUp(*A);
38: MatGetOwnershipRange(*A,&Istart,&Iend);
39: for (II=Istart;II<Iend;II++) {
40: i = II/n; j = II-i*n;
41: if (i>0) { MatSetValue(*A,II,II-n,-1.0,INSERT_VALUES); }
42: if (i<m-1) { MatSetValue(*A,II,II+n,-1.0,INSERT_VALUES); }
43: if (j>0) { MatSetValue(*A,II,II-1,-1.0,INSERT_VALUES); }
44: if (j<n-1) { MatSetValue(*A,II,II+1,-1.0,INSERT_VALUES); }
45: MatSetValue(*A,II,II,4.0,INSERT_VALUES);
46: }
47: MatAssemblyBegin(*A,MAT_FINAL_ASSEMBLY);
48: MatAssemblyEnd(*A,MAT_FINAL_ASSEMBLY);
49: return(0);
50: }
52: int main(int argc,char **argv)
53: {
54: Mat A,B;
55: EPS eps;
56: PetscInt N,n=10,m=11,nev=3;
57: PetscMPIInt size;
58: PetscBool flag,terse;
61: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
62: MPI_Comm_size(PETSC_COMM_WORLD,&size);
63: if (size==1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This example requires at least two processes");
64: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
65: PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag);
66: N = n*m;
67: PetscPrintf(PETSC_COMM_WORLD,"\n2-D Laplacian Eigenproblem, N=%D (%Dx%D grid)\n\n",N,n,m);
69: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
70: Create 2-D Laplacian matrices
71: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
73: Laplacian(PETSC_COMM_WORLD,n,m,1,&A);
74: Laplacian(PETSC_COMM_WORLD,n,m,-1,&B);
76: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
77: Create the eigensolver, set options and solve the eigensystem
78: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
80: PetscPrintf(PETSC_COMM_WORLD,"First solve:\n\n");
81: EPSCreate(PETSC_COMM_WORLD,&eps);
82: EPSSetOperators(eps,A,NULL);
83: EPSSetProblemType(eps,EPS_HEP);
84: EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);
85: EPSSetDimensions(eps,nev,PETSC_DEFAULT,PETSC_DEFAULT);
86: EPSSetFromOptions(eps);
88: EPSSolve(eps);
90: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
91: Display solution of first solve
92: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
94: PetscOptionsHasName(NULL,NULL,"-terse",&terse);
95: if (terse) {
96: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
97: } else {
98: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL);
99: EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD);
100: EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD);
101: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
102: }
104: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
105: Solve with second matrix
106: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
108: PetscPrintf(PETSC_COMM_WORLD,"\nSecond solve:\n\n");
109: /*EPSReset(eps);*/ /* not required, will be called in EPSSetOperators() */
110: EPSSetOperators(eps,B,NULL);
111: EPSSolve(eps);
113: if (terse) {
114: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
115: } else {
116: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL);
117: EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD);
118: EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD);
119: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
120: }
122: EPSDestroy(&eps);
123: MatDestroy(&A);
124: MatDestroy(&B);
125: SlepcFinalize();
126: return ierr;
127: }
129: /*TEST
131: testset:
132: nsize: 2
133: requires: !single
134: output_file: output/test39_1.out
135: test:
136: suffix: 1
137: args: -eps_type {{krylovschur arnoldi lobpcg lapack}} -terse
138: test:
139: suffix: 1_lanczos
140: args: -eps_type lanczos -eps_lanczos_reorthog local -terse
142: TEST*/