Actual source code: test3.c
slepc-3.20.1 2023-11-27
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.\n\n";
13: #include <slepceps.h>
15: int main(int argc,char **argv)
16: {
17: Mat A1,A2; /* problem matrices */
18: EPS eps; /* eigenproblem solver context */
19: PetscReal tol=PETSC_SMALL,v;
20: Vec d;
21: PetscInt n=30,i,Istart,Iend;
22: PetscRandom myrand;
24: PetscFunctionBeginUser;
25: PetscCall(SlepcInitialize(&argc,&argv,(char*)0,help));
27: PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
28: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal with random diagonal, n=%" PetscInt_FMT "\n\n",n));
30: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
31: Create matrix tridiag([-1 0 -1])
32: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
33: PetscCall(MatCreate(PETSC_COMM_WORLD,&A1));
34: PetscCall(MatSetSizes(A1,PETSC_DECIDE,PETSC_DECIDE,n,n));
35: PetscCall(MatSetFromOptions(A1));
36: PetscCall(MatSetUp(A1));
38: PetscCall(MatGetOwnershipRange(A1,&Istart,&Iend));
39: for (i=Istart;i<Iend;i++) {
40: if (i>0) PetscCall(MatSetValue(A1,i,i-1,-1.0,INSERT_VALUES));
41: if (i<n-1) PetscCall(MatSetValue(A1,i,i+1,-1.0,INSERT_VALUES));
42: }
43: PetscCall(MatAssemblyBegin(A1,MAT_FINAL_ASSEMBLY));
44: PetscCall(MatAssemblyEnd(A1,MAT_FINAL_ASSEMBLY));
46: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: Create two matrices by filling the diagonal with rand values
48: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
49: PetscCall(MatDuplicate(A1,MAT_COPY_VALUES,&A2));
50: PetscCall(MatCreateVecs(A1,NULL,&d));
51: PetscCall(PetscRandomCreate(PETSC_COMM_WORLD,&myrand));
52: PetscCall(PetscRandomSetFromOptions(myrand));
53: PetscCall(PetscRandomSetInterval(myrand,0.0,1.0));
54: for (i=Istart;i<Iend;i++) {
55: PetscCall(PetscRandomGetValueReal(myrand,&v));
56: PetscCall(VecSetValue(d,i,v,INSERT_VALUES));
57: }
58: PetscCall(VecAssemblyBegin(d));
59: PetscCall(VecAssemblyEnd(d));
60: PetscCall(MatDiagonalSet(A1,d,INSERT_VALUES));
61: for (i=Istart;i<Iend;i++) {
62: PetscCall(PetscRandomGetValueReal(myrand,&v));
63: PetscCall(VecSetValue(d,i,v,INSERT_VALUES));
64: }
65: PetscCall(VecAssemblyBegin(d));
66: PetscCall(VecAssemblyEnd(d));
67: PetscCall(MatDiagonalSet(A2,d,INSERT_VALUES));
68: PetscCall(VecDestroy(&d));
69: PetscCall(PetscRandomDestroy(&myrand));
71: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: Create the eigensolver
73: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
74: PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
75: PetscCall(EPSSetProblemType(eps,EPS_HEP));
76: PetscCall(EPSSetTolerances(eps,tol,PETSC_DEFAULT));
77: PetscCall(EPSSetOperators(eps,A1,NULL));
78: PetscCall(EPSSetFromOptions(eps));
80: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
81: Solve first eigenproblem
82: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
83: PetscCall(EPSSolve(eps));
84: PetscCall(PetscPrintf(PETSC_COMM_WORLD," - - - First matrix - - -\n"));
85: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
87: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
88: Solve second eigenproblem
89: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
90: PetscCall(EPSSetOperators(eps,A2,NULL));
91: PetscCall(EPSSolve(eps));
92: PetscCall(PetscPrintf(PETSC_COMM_WORLD," - - - Second matrix - - -\n"));
93: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
95: PetscCall(EPSDestroy(&eps));
96: PetscCall(MatDestroy(&A1));
97: PetscCall(MatDestroy(&A2));
98: PetscCall(SlepcFinalize());
99: return 0;
100: }
102: /*TEST
104: testset:
105: args: -eps_nev 4
106: requires: !single
107: output_file: output/test3_1.out
108: test:
109: suffix: 1
110: args: -eps_type {{krylovschur subspace arnoldi lapack}}
111: test:
112: suffix: 1_lanczos
113: args: -eps_type lanczos -eps_lanczos_reorthog local
114: test:
115: suffix: 1_power
116: args: -eps_type power -eps_max_it 20000
117: test:
118: suffix: 1_jd
119: args: -eps_type jd -eps_jd_initial_size 7
120: test:
121: suffix: 1_gd
122: args: -eps_type gd -eps_gd_initial_size 7
123: test:
124: suffix: 1_gd2
125: args: -eps_type gd -eps_gd_double_expansion
126: test:
127: suffix: 1_arpack
128: args: -eps_type arpack
129: requires: arpack
130: test:
131: suffix: 1_primme
132: args: -eps_type primme -eps_conv_abs -eps_primme_blocksize 4
133: requires: primme
134: test:
135: suffix: 1_trlan
136: args: -eps_type trlan
137: requires: trlan
138: test:
139: suffix: 1_scalapack
140: args: -eps_type scalapack
141: requires: scalapack
142: test:
143: suffix: 1_elpa
144: args: -eps_type elpa
145: requires: elpa
146: test:
147: suffix: 1_elemental
148: args: -eps_type elemental
149: requires: elemental
151: testset:
152: args: -eps_nev 4 -eps_smallest_real -eps_max_it 500
153: output_file: output/test3_2.out
154: test:
155: suffix: 2_rqcg
156: args: -eps_type rqcg -eps_rqcg_reset 5 -eps_ncv 32
157: test:
158: suffix: 2_lobpcg
159: args: -eps_type lobpcg -eps_lobpcg_blocksize 5 -st_pc_type none
160: test:
161: suffix: 2_lanczos
162: args: -eps_type lanczos -eps_lanczos_reorthog local
163: requires: !single
164: test:
165: suffix: 2_lanczos_delayed
166: args: -eps_type lanczos -eps_lanczos_reorthog delayed -eps_tol 1e-8
167: requires: !single
168: test:
169: suffix: 2_trlan
170: args: -eps_type trlan
171: requires: trlan
172: test:
173: suffix: 2_blopex
174: args: -eps_type blopex -eps_conv_abs -st_shift -2
175: requires: blopex
177: TEST*/