Actual source code: ex80.c
1: static char help[] = "Test the Fischer-3 initial guess routine.\n\n";
3: #include <petscksp.h>
5: #define SIZE 3
7: int main(int argc,char **args)
8: {
9: PetscInt i;
10: {
11: Mat A;
12: PetscInt indices[SIZE] = {0,1,2};
13: PetscScalar values[SIZE] = {1.0,1.0,1.0};
14: Vec sol,rhs,newsol,newrhs;
16: PetscInitialize(&argc,&args,(char*)0,help);
18: /* common data structures */
19: MatCreateSeqDense(PETSC_COMM_SELF,SIZE,SIZE,NULL,&A);
20: for (i = 0; i < SIZE; ++i) {
21: MatSetValue(A,i,i,1.0,INSERT_VALUES);
22: }
23: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
24: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
26: VecCreateSeq(PETSC_COMM_SELF,SIZE,&sol);
27: VecDuplicate(sol,&rhs);
28: VecDuplicate(sol,&newrhs);
29: VecDuplicate(sol,&newsol);
31: VecSetValues(sol,SIZE,indices,values,INSERT_VALUES);
32: VecSetValues(rhs,SIZE - 1,indices,values,INSERT_VALUES);
33: VecSetValues(newrhs,SIZE - 2,indices,values,INSERT_VALUES);
34: VecAssemblyBegin(sol);
35: VecAssemblyBegin(rhs);
36: VecAssemblyBegin(newrhs);
37: VecAssemblyEnd(sol);
38: VecAssemblyEnd(rhs);
39: VecAssemblyEnd(newrhs);
41: /* Test one vector */
42: {
43: KSP ksp;
44: KSPGuess guess;
46: KSPCreate(PETSC_COMM_SELF,&ksp);
47: KSPSetOperators(ksp,A,A);
48: KSPSetFromOptions(ksp);
49: KSPGetGuess(ksp,&guess);
50: /* we aren't calling through the KSP so we call this ourselves */
51: KSPGuessSetUp(guess);
53: KSPGuessUpdate(guess,rhs,sol);
54: KSPGuessFormGuess(guess,newrhs,newsol);
55: VecView(newsol,PETSC_VIEWER_STDOUT_SELF);
57: KSPDestroy(&ksp);
58: }
60: /* Test a singular projection matrix */
61: {
62: KSP ksp;
63: KSPGuess guess;
65: KSPCreate(PETSC_COMM_SELF,&ksp);
66: KSPSetOperators(ksp,A,A);
67: KSPSetFromOptions(ksp);
68: KSPGetGuess(ksp,&guess);
69: KSPGuessSetUp(guess);
71: for (i = 0; i < 15; ++i) {
72: KSPGuessUpdate(guess,rhs,sol);
73: }
74: KSPGuessFormGuess(guess,newrhs,newsol);
75: VecView(newsol,PETSC_VIEWER_STDOUT_SELF);
77: KSPDestroy(&ksp);
78: }
79: VecDestroy(&newsol);
80: VecDestroy(&newrhs);
81: VecDestroy(&rhs);
82: VecDestroy(&sol);
84: MatDestroy(&A);
85: }
87: /* Test something triangular */
88: {
89: PetscInt triangle_size = 10;
90: Mat A;
92: MatCreateSeqDense(PETSC_COMM_SELF,triangle_size,triangle_size,NULL,&A);
93: for (i = 0; i < triangle_size; ++i) {
94: MatSetValue(A,i,i,1.0,INSERT_VALUES);
95: }
96: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
97: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
99: {
100: KSP ksp;
101: KSPGuess guess;
102: Vec sol,rhs;
103: PetscInt j,indices[] = {0,1,2,3,4};
104: PetscScalar values[] = {1.0,2.0,3.0,4.0,5.0};
106: KSPCreate(PETSC_COMM_SELF,&ksp);
107: KSPSetOperators(ksp,A,A);
108: KSPSetFromOptions(ksp);
109: KSPGetGuess(ksp,&guess);
110: KSPGuessSetUp(guess);
112: for (i = 0; i < 5; ++i) {
113: VecCreateSeq(PETSC_COMM_SELF,triangle_size,&sol);
114: VecCreateSeq(PETSC_COMM_SELF,triangle_size,&rhs);
115: for (j = 0; j < i; ++j) {
116: VecSetValue(sol,j,(PetscScalar)j,INSERT_VALUES);
117: VecSetValue(rhs,j,(PetscScalar)j,INSERT_VALUES);
118: }
119: VecAssemblyBegin(sol);
120: VecAssemblyBegin(rhs);
121: VecAssemblyEnd(sol);
122: VecAssemblyEnd(rhs);
124: KSPGuessUpdate(guess,rhs,sol);
126: VecDestroy(&rhs);
127: VecDestroy(&sol);
128: }
130: VecCreateSeq(PETSC_COMM_SELF,triangle_size,&sol);
131: VecCreateSeq(PETSC_COMM_SELF,triangle_size,&rhs);
132: VecSetValues(rhs,5,indices,values,INSERT_VALUES);
133: VecAssemblyBegin(sol);
134: VecAssemblyEnd(sol);
136: KSPGuessFormGuess(guess,rhs,sol);
137: VecView(sol,PETSC_VIEWER_STDOUT_SELF);
139: VecDestroy(&rhs);
140: VecDestroy(&sol);
141: KSPDestroy(&ksp);
142: }
143: MatDestroy(&A);
144: }
145: PetscFinalize();
146: return 0;
147: }
149: /* The relative tolerance here is strict enough to get rid of all the noise in both single and double precision: values as low as 5e-7 also work */
151: /*TEST
153: test:
154: args: -ksp_guess_type fischer -ksp_guess_fischer_model 3,10 -ksp_guess_fischer_monitor -ksp_guess_fischer_tol 1e-6
156: TEST*/