Actual source code: test5.c
slepc-3.23.0 2025-03-29
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[] = "Test various ST interface functions.\n\n";
13: #include <slepceps.h>
15: int main (int argc,char **argv)
16: {
17: ST st;
18: KSP ksp;
19: PC pc;
20: Mat A,mat[1],Op;
21: Vec v,w;
22: PetscInt N,n=4,i,j,II,Istart,Iend;
23: PetscScalar d,val1,val2;
24: PetscBool test_compl=PETSC_FALSE;
26: PetscFunctionBeginUser;
27: PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
28: PetscCall(PetscOptionsGetBool(NULL,NULL,"-complex",&test_compl,NULL));
29: PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
30: N = n*n;
32: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
33: Create non-symmetric matrix
34: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
36: PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
37: PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
38: PetscCall(MatSetFromOptions(A));
40: PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
41: #if defined(PETSC_USE_COMPLEX)
42: val1 = test_compl? PetscCMPLX(1.0,-0.2): 1.0;
43: val2 = test_compl? PetscCMPLX(-1.0,0.4): -1.0;
44: #else
45: val1 = 1.0;
46: val2 = -1.0;
47: #endif
48: for (II=Istart;II<Iend;II++) {
49: i = II/n; j = II-i*n;
50: d = 0.0;
51: if (i>0) { PetscCall(MatSetValue(A,II,II-n,val1,INSERT_VALUES)); d=d+1.0; }
52: if (i<n-1) { PetscCall(MatSetValue(A,II,II+n,val2,INSERT_VALUES)); d=d+1.0; }
53: if (j>0) { PetscCall(MatSetValue(A,II,II-1,val1,INSERT_VALUES)); d=d+1.0; }
54: if (j<n-1) { PetscCall(MatSetValue(A,II,II+1,val2,INSERT_VALUES)); d=d+1.0; }
55: PetscCall(MatSetValue(A,II,II,d,INSERT_VALUES));
56: }
57: PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
58: PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
60: #if defined(PETSC_USE_COMPLEX)
61: val1 = test_compl? PetscCMPLX(-0.5,0.01): -0.5;
62: #else
63: val1 = -0.5;
64: #endif
65: PetscCall(MatCreateVecs(A,&v,&w));
66: PetscCall(VecSetValue(v,0,val1,INSERT_VALUES));
67: PetscCall(VecSetValue(v,1,1.5,INSERT_VALUES));
68: PetscCall(VecSetValue(v,2,2,INSERT_VALUES));
69: PetscCall(VecAssemblyBegin(v));
70: PetscCall(VecAssemblyEnd(v));
71: PetscCall(VecView(v,NULL));
73: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: Create the spectral transformation object
75: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
76: PetscCall(STCreate(PETSC_COMM_WORLD,&st));
77: mat[0] = A;
78: PetscCall(STSetMatrices(st,1,mat));
79: PetscCall(STSetType(st,STCAYLEY));
80: PetscCall(STSetShift(st,2.0));
81: PetscCall(STCayleySetAntishift(st,1.0));
82: PetscCall(STSetTransform(st,PETSC_TRUE));
84: PetscCall(KSPCreate(PETSC_COMM_WORLD,&ksp));
85: PetscCall(KSPSetType(ksp,KSPPREONLY));
86: PetscCall(KSPGetPC(ksp,&pc));
87: PetscCall(PCSetType(pc,PCLU));
88: PetscCall(KSPSetTolerances(ksp,100*PETSC_MACHINE_EPSILON,PETSC_CURRENT,PETSC_CURRENT,PETSC_CURRENT));
89: PetscCall(KSPSetFromOptions(ksp));
90: PetscCall(STSetKSP(st,ksp));
91: PetscCall(KSPDestroy(&ksp));
93: PetscCall(STSetFromOptions(st));
95: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
96: Apply the operator
97: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
98: PetscCall(STApply(st,v,w));
99: PetscCall(VecView(w,NULL));
101: PetscCall(STApplyTranspose(st,v,w));
102: PetscCall(VecView(w,NULL));
104: if (test_compl) {
105: PetscCall(STApplyHermitianTranspose(st,v,w));
106: PetscCall(VecView(w,NULL));
107: }
109: PetscCall(STMatMult(st,1,v,w));
110: PetscCall(VecView(w,NULL));
112: PetscCall(STMatMultTranspose(st,1,v,w));
113: PetscCall(VecView(w,NULL));
115: if (test_compl) {
116: PetscCall(STMatMultHermitianTranspose(st,1,v,w));
117: PetscCall(VecView(w,NULL));
118: }
120: PetscCall(STMatSolve(st,v,w));
121: PetscCall(VecView(w,NULL));
123: PetscCall(STMatSolveTranspose(st,v,w));
124: PetscCall(VecView(w,NULL));
126: if (test_compl) {
127: PetscCall(STMatSolveHermitianTranspose(st,v,w));
128: PetscCall(VecView(w,NULL));
129: }
131: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
132: Get the operator matrix
133: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
134: PetscCall(STGetOperator(st,&Op));
135: PetscCall(MatMult(Op,v,w));
136: PetscCall(VecView(w,NULL));
137: PetscCall(MatMultTranspose(Op,v,w));
138: PetscCall(VecView(w,NULL));
139: if (test_compl) {
140: PetscCall(MatMultHermitianTranspose(Op,v,w));
141: PetscCall(VecView(w,NULL));
142: }
143: PetscCall(STRestoreOperator(st,&Op));
145: PetscCall(STDestroy(&st));
146: PetscCall(MatDestroy(&A));
147: PetscCall(VecDestroy(&v));
148: PetscCall(VecDestroy(&w));
149: PetscCall(SlepcFinalize());
150: return 0;
151: }
153: /*TEST
155: testset:
156: output_file: output/test5_1.out
157: requires: !single
158: test:
159: suffix: 1
160: args: -st_matmode {{copy inplace}}
161: test:
162: suffix: 1_shell
163: args: -st_matmode shell -ksp_type bcgs -pc_type jacobi
165: testset:
166: output_file: output/test5_2.out
167: requires: complex !single
168: args: -complex
169: test:
170: suffix: 2
171: args: -st_matmode {{copy inplace}}
172: test:
173: suffix: 2_shell
174: args: -st_matmode shell -ksp_type bcgs -pc_type jacobi
176: TEST*/