Actual source code: test36.c
slepc-3.14.0 2020-09-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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 a HEP problem with Hermitian matrix.\n\n";
13: #include <slepceps.h>
15: int main(int argc,char **argv)
16: {
17: Mat A; /* matrix */
18: EPS eps; /* eigenproblem solver context */
19: PetscInt N,n=20,m,Istart,Iend,II,i,j;
20: PetscBool flag;
23: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
24: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
25: PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag);
26: if (!flag) m=n;
27: N = n*m;
28: PetscPrintf(PETSC_COMM_WORLD,"\nHermitian Eigenproblem, N=%D (%Dx%D grid)\n\n",N,n,m);
29: #if !defined(PETSC_USE_COMPLEX)
30: SETERRQ(PETSC_COMM_WORLD,1,"This example requires complex scalars!");
31: #endif
33: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
34: Compute the matrix that defines the eigensystem, Ax=kx
35: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
37: MatCreate(PETSC_COMM_WORLD,&A);
38: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N);
39: MatSetFromOptions(A);
40: MatSetUp(A);
42: MatGetOwnershipRange(A,&Istart,&Iend);
43: for (II=Istart;II<Iend;II++) {
44: i = II/n; j = II-i*n;
45: if (i>0) { MatSetValue(A,II,II-n,-1.0-0.1*PETSC_i,INSERT_VALUES); }
46: if (i<m-1) { MatSetValue(A,II,II+n,-1.0+0.1*PETSC_i,INSERT_VALUES); }
47: if (j>0) { MatSetValue(A,II,II-1,-1.0-0.1*PETSC_i,INSERT_VALUES); }
48: if (j<n-1) { MatSetValue(A,II,II+1,-1.0+0.1*PETSC_i,INSERT_VALUES); }
49: MatSetValue(A,II,II,4.0,INSERT_VALUES);
50: }
51: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
52: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
54: MatSetOption(A,MAT_HERMITIAN,PETSC_TRUE);
56: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
57: Create the eigensolver and solve the problem
58: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
60: EPSCreate(PETSC_COMM_WORLD,&eps);
61: EPSSetOperators(eps,A,NULL);
62: EPSSetProblemType(eps,EPS_HEP);
63: EPSSetFromOptions(eps);
64: EPSSolve(eps);
65: EPSErrorView(eps,EPS_ERROR_BACKWARD,NULL);
67: EPSDestroy(&eps);
68: MatDestroy(&A);
69: SlepcFinalize();
70: return ierr;
71: }
73: /*TEST
75: build:
76: requires: complex
78: testset:
79: args: -m 18 -n 19 -eps_nev 4 -eps_max_it 1000
80: requires: !single complex
81: output_file: output/test36_1.out
82: test:
83: suffix: 1
84: args: -eps_type {{krylovschur subspace arnoldi gd jd lapack}}
85: test:
86: suffix: 1_elemental
87: args: -eps_type elemental
88: requires: elemental
90: test:
91: suffix: 2
92: args: -eps_nev 4 -eps_smallest_real -eps_type {{lobpcg rqcg lapack}}
93: requires: !single complex
95: TEST*/