Actual source code: test4.c
slepc-3.18.0 2022-10-01
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 ST with four matrices.\n\n";
13: #include <slepcst.h>
15: int main(int argc,char **argv)
16: {
17: Mat A,B,C,D,mat[4];
18: ST st;
19: KSP ksp;
20: Vec v,w;
21: STType type;
22: PetscScalar sigma;
23: PetscInt n=10,i,Istart,Iend;
26: SlepcInitialize(&argc,&argv,(char*)0,help);
27: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
28: PetscPrintf(PETSC_COMM_WORLD,"\nTest ST with four matrices, n=%" PetscInt_FMT "\n\n",n);
29: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
30: Compute the operator matrices
31: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
33: MatCreate(PETSC_COMM_WORLD,&A);
34: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
35: MatSetFromOptions(A);
36: MatSetUp(A);
38: MatCreate(PETSC_COMM_WORLD,&B);
39: MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,n,n);
40: MatSetFromOptions(B);
41: MatSetUp(B);
43: MatCreate(PETSC_COMM_WORLD,&C);
44: MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,n,n);
45: MatSetFromOptions(C);
46: MatSetUp(C);
48: MatCreate(PETSC_COMM_WORLD,&D);
49: MatSetSizes(D,PETSC_DECIDE,PETSC_DECIDE,n,n);
50: MatSetFromOptions(D);
51: MatSetUp(D);
53: MatGetOwnershipRange(A,&Istart,&Iend);
54: for (i=Istart;i<Iend;i++) {
55: MatSetValue(A,i,i,2.0,INSERT_VALUES);
56: if (i>0) {
57: MatSetValue(A,i,i-1,-1.0,INSERT_VALUES);
58: MatSetValue(B,i,i,(PetscScalar)i,INSERT_VALUES);
59: } else MatSetValue(B,i,i,-1.0,INSERT_VALUES);
60: if (i<n-1) MatSetValue(A,i,i+1,-1.0,INSERT_VALUES);
61: MatSetValue(C,i,n-i-1,1.0,INSERT_VALUES);
62: MatSetValue(D,i,i,i*.1,INSERT_VALUES);
63: if (i==0) MatSetValue(D,0,n-1,1.0,INSERT_VALUES);
64: if (i==n-1) MatSetValue(D,n-1,0,1.0,INSERT_VALUES);
65: }
67: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
68: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
69: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
70: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
71: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
72: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
73: MatAssemblyBegin(D,MAT_FINAL_ASSEMBLY);
74: MatAssemblyEnd(D,MAT_FINAL_ASSEMBLY);
75: MatCreateVecs(A,&v,&w);
76: VecSet(v,1.0);
78: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
79: Create the spectral transformation object
80: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
82: STCreate(PETSC_COMM_WORLD,&st);
83: mat[0] = A;
84: mat[1] = B;
85: mat[2] = C;
86: mat[3] = D;
87: STSetMatrices(st,4,mat);
88: STGetKSP(st,&ksp);
89: KSPSetTolerances(ksp,100*PETSC_MACHINE_EPSILON,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
90: STSetFromOptions(st);
92: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
93: Apply the transformed operator for several ST's
94: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
96: /* shift, sigma=0.0 */
97: STSetUp(st);
98: STGetType(st,&type);
99: PetscPrintf(PETSC_COMM_WORLD,"ST type %s\n",type);
100: for (i=0;i<4;i++) {
101: STMatMult(st,i,v,w);
102: PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i);
103: VecView(w,NULL);
104: }
105: STMatSolve(st,v,w);
106: PetscPrintf(PETSC_COMM_WORLD,"solve\n");
107: VecView(w,NULL);
109: /* shift, sigma=0.1 */
110: sigma = 0.1;
111: STSetShift(st,sigma);
112: STGetShift(st,&sigma);
113: PetscPrintf(PETSC_COMM_WORLD,"With shift=%g\n",(double)PetscRealPart(sigma));
114: for (i=0;i<4;i++) {
115: STMatMult(st,i,v,w);
116: PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i);
117: VecView(w,NULL);
118: }
119: STMatSolve(st,v,w);
120: PetscPrintf(PETSC_COMM_WORLD,"solve\n");
121: VecView(w,NULL);
123: /* sinvert, sigma=0.1 */
124: STPostSolve(st);
125: STSetType(st,STSINVERT);
126: STGetType(st,&type);
127: PetscPrintf(PETSC_COMM_WORLD,"ST type %s\n",type);
128: STGetShift(st,&sigma);
129: PetscPrintf(PETSC_COMM_WORLD,"With shift=%g\n",(double)PetscRealPart(sigma));
130: for (i=0;i<4;i++) {
131: STMatMult(st,i,v,w);
132: PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i);
133: VecView(w,NULL);
134: }
135: STMatSolve(st,v,w);
136: PetscPrintf(PETSC_COMM_WORLD,"solve\n");
137: VecView(w,NULL);
139: /* sinvert, sigma=-0.5 */
140: sigma = -0.5;
141: STSetShift(st,sigma);
142: STGetShift(st,&sigma);
143: PetscPrintf(PETSC_COMM_WORLD,"With shift=%g\n",(double)PetscRealPart(sigma));
144: for (i=0;i<4;i++) {
145: STMatMult(st,i,v,w);
146: PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i);
147: VecView(w,NULL);
148: }
149: STMatSolve(st,v,w);
150: PetscPrintf(PETSC_COMM_WORLD,"solve\n");
151: VecView(w,NULL);
153: STDestroy(&st);
154: MatDestroy(&A);
155: MatDestroy(&B);
156: MatDestroy(&C);
157: MatDestroy(&D);
158: VecDestroy(&v);
159: VecDestroy(&w);
160: SlepcFinalize();
161: return 0;
162: }
164: /*TEST
166: test:
167: suffix: 1
168: args: -st_transform -st_matmode {{copy shell}}
169: output_file: output/test4_1.out
170: requires: !single
172: test:
173: suffix: 2
174: args: -st_matmode {{copy shell}}
175: output_file: output/test4_2.out
177: TEST*/