Actual source code: test4.c
slepc-3.21.1 2024-04-26
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;
25: PetscFunctionBeginUser;
26: PetscCall(SlepcInitialize(&argc,&argv,(char*)0,help));
27: PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
28: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nTest ST with four matrices, n=%" PetscInt_FMT "\n\n",n));
29: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
30: Compute the operator matrices
31: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
33: PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
34: PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n));
35: PetscCall(MatSetFromOptions(A));
37: PetscCall(MatCreate(PETSC_COMM_WORLD,&B));
38: PetscCall(MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,n,n));
39: PetscCall(MatSetFromOptions(B));
41: PetscCall(MatCreate(PETSC_COMM_WORLD,&C));
42: PetscCall(MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,n,n));
43: PetscCall(MatSetFromOptions(C));
45: PetscCall(MatCreate(PETSC_COMM_WORLD,&D));
46: PetscCall(MatSetSizes(D,PETSC_DECIDE,PETSC_DECIDE,n,n));
47: PetscCall(MatSetFromOptions(D));
49: PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
50: for (i=Istart;i<Iend;i++) {
51: PetscCall(MatSetValue(A,i,i,2.0,INSERT_VALUES));
52: if (i>0) {
53: PetscCall(MatSetValue(A,i,i-1,-1.0,INSERT_VALUES));
54: PetscCall(MatSetValue(B,i,i,(PetscScalar)i,INSERT_VALUES));
55: } else PetscCall(MatSetValue(B,i,i,-1.0,INSERT_VALUES));
56: if (i<n-1) PetscCall(MatSetValue(A,i,i+1,-1.0,INSERT_VALUES));
57: PetscCall(MatSetValue(C,i,n-i-1,1.0,INSERT_VALUES));
58: PetscCall(MatSetValue(D,i,i,i*.1,INSERT_VALUES));
59: if (i==0) PetscCall(MatSetValue(D,0,n-1,1.0,INSERT_VALUES));
60: if (i==n-1) PetscCall(MatSetValue(D,n-1,0,1.0,INSERT_VALUES));
61: }
63: PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
64: PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
65: PetscCall(MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY));
66: PetscCall(MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY));
67: PetscCall(MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY));
68: PetscCall(MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY));
69: PetscCall(MatAssemblyBegin(D,MAT_FINAL_ASSEMBLY));
70: PetscCall(MatAssemblyEnd(D,MAT_FINAL_ASSEMBLY));
71: PetscCall(MatCreateVecs(A,&v,&w));
72: PetscCall(VecSet(v,1.0));
74: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
75: Create the spectral transformation object
76: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
78: PetscCall(STCreate(PETSC_COMM_WORLD,&st));
79: mat[0] = A;
80: mat[1] = B;
81: mat[2] = C;
82: mat[3] = D;
83: PetscCall(STSetMatrices(st,4,mat));
84: PetscCall(STGetKSP(st,&ksp));
85: PetscCall(KSPSetTolerances(ksp,100*PETSC_MACHINE_EPSILON,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT));
86: PetscCall(STSetFromOptions(st));
88: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
89: Apply the transformed operator for several ST's
90: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
92: /* shift, sigma=0.0 */
93: PetscCall(STSetUp(st));
94: PetscCall(STGetType(st,&type));
95: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"ST type %s\n",type));
96: for (i=0;i<4;i++) {
97: PetscCall(STMatMult(st,i,v,w));
98: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i));
99: PetscCall(VecView(w,NULL));
100: }
101: PetscCall(STMatSolve(st,v,w));
102: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"solve\n"));
103: PetscCall(VecView(w,NULL));
105: /* shift, sigma=0.1 */
106: sigma = 0.1;
107: PetscCall(STSetShift(st,sigma));
108: PetscCall(STGetShift(st,&sigma));
109: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"With shift=%g\n",(double)PetscRealPart(sigma)));
110: for (i=0;i<4;i++) {
111: PetscCall(STMatMult(st,i,v,w));
112: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i));
113: PetscCall(VecView(w,NULL));
114: }
115: PetscCall(STMatSolve(st,v,w));
116: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"solve\n"));
117: PetscCall(VecView(w,NULL));
119: /* sinvert, sigma=0.1 */
120: PetscCall(STPostSolve(st));
121: PetscCall(STSetType(st,STSINVERT));
122: PetscCall(STGetType(st,&type));
123: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"ST type %s\n",type));
124: PetscCall(STGetShift(st,&sigma));
125: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"With shift=%g\n",(double)PetscRealPart(sigma)));
126: for (i=0;i<4;i++) {
127: PetscCall(STMatMult(st,i,v,w));
128: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i));
129: PetscCall(VecView(w,NULL));
130: }
131: PetscCall(STMatSolve(st,v,w));
132: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"solve\n"));
133: PetscCall(VecView(w,NULL));
135: /* sinvert, sigma=-0.5 */
136: sigma = -0.5;
137: PetscCall(STSetShift(st,sigma));
138: PetscCall(STGetShift(st,&sigma));
139: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"With shift=%g\n",(double)PetscRealPart(sigma)));
140: for (i=0;i<4;i++) {
141: PetscCall(STMatMult(st,i,v,w));
142: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"k= %" PetscInt_FMT "\n",i));
143: PetscCall(VecView(w,NULL));
144: }
145: PetscCall(STMatSolve(st,v,w));
146: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"solve\n"));
147: PetscCall(VecView(w,NULL));
149: PetscCall(STDestroy(&st));
150: PetscCall(MatDestroy(&A));
151: PetscCall(MatDestroy(&B));
152: PetscCall(MatDestroy(&C));
153: PetscCall(MatDestroy(&D));
154: PetscCall(VecDestroy(&v));
155: PetscCall(VecDestroy(&w));
156: PetscCall(SlepcFinalize());
157: return 0;
158: }
160: /*TEST
162: test:
163: suffix: 1
164: args: -st_transform -st_matmode {{copy shell}}
165: requires: !single
167: test:
168: suffix: 2
169: args: -st_matmode {{copy shell}}
171: TEST*/