Actual source code: ex83.c
1: static char help[] = "Test the Fischer-1 initial guess routine with VECNEST.\n\n";
3: #include <petscksp.h>
5: int main(int argc, char **args)
6: {
7: /* Test that exceeding the number of stored vectors works correctly - this used to not work with VecNest */
8: PetscInt triangle_size = 10;
9: Mat A, A_nest;
10: KSP ksp;
11: KSPGuess guess;
12: Vec sol, rhs, sol_nest, rhs_nest;
13: PetscInt i, j, indices[] = {0, 1, 2, 3, 4};
14: PetscScalar values[] = {1.0, 2.0, 3.0, 4.0, 5.0};
16: PetscInitialize(&argc, &args, (char *)0, help);
18: MatCreateSeqDense(PETSC_COMM_SELF, triangle_size, triangle_size, NULL, &A);
19: for (i = 0; i < triangle_size; ++i) { MatSetValue(A, i, i, 1.0, INSERT_VALUES); }
20: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
21: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
22: MatCreateNest(PETSC_COMM_SELF, 1, NULL, 1, NULL, &A, &A_nest);
23: MatNestSetVecType(A_nest, VECNEST);
25: KSPCreate(PETSC_COMM_SELF, &ksp);
26: KSPSetOperators(ksp, A_nest, A_nest);
27: KSPSetFromOptions(ksp);
28: KSPGetGuess(ksp, &guess);
29: KSPGuessSetUp(guess);
31: for (i = 0; i < 5; ++i) {
32: VecCreateSeq(PETSC_COMM_SELF, triangle_size, &sol);
33: VecCreateSeq(PETSC_COMM_SELF, triangle_size, &rhs);
34: for (j = 0; j < i; ++j) {
35: VecSetValue(sol, j, (PetscScalar)j, INSERT_VALUES);
36: VecSetValue(rhs, j, (PetscScalar)j, INSERT_VALUES);
37: }
38: VecAssemblyBegin(sol);
39: VecAssemblyBegin(rhs);
40: VecAssemblyEnd(sol);
41: VecAssemblyEnd(rhs);
43: VecCreateNest(PETSC_COMM_SELF, 1, NULL, &sol, &sol_nest);
44: VecCreateNest(PETSC_COMM_SELF, 1, NULL, &rhs, &rhs_nest);
46: KSPGuessUpdate(guess, rhs_nest, sol_nest);
48: VecDestroy(&rhs_nest);
49: VecDestroy(&sol_nest);
50: VecDestroy(&rhs);
51: VecDestroy(&sol);
52: }
54: VecCreateSeq(PETSC_COMM_SELF, triangle_size, &sol);
55: VecCreateSeq(PETSC_COMM_SELF, triangle_size, &rhs);
56: VecSetValues(rhs, 5, indices, values, INSERT_VALUES);
57: VecAssemblyBegin(rhs);
58: VecAssemblyEnd(rhs);
60: VecCreateNest(PETSC_COMM_SELF, 1, NULL, &sol, &sol_nest);
61: VecCreateNest(PETSC_COMM_SELF, 1, NULL, &rhs, &rhs_nest);
63: KSPGuessFormGuess(guess, rhs_nest, sol_nest);
64: VecView(sol_nest, PETSC_VIEWER_STDOUT_SELF);
66: VecDestroy(&rhs_nest);
67: VecDestroy(&sol_nest);
68: VecDestroy(&rhs);
69: VecDestroy(&sol);
71: KSPDestroy(&ksp);
73: MatDestroy(&A_nest);
74: MatDestroy(&A);
76: PetscFinalize();
77: return 0;
78: }
80: /*TEST
82: test:
83: args: -ksp_guess_type fischer -ksp_guess_fischer_model 1,3 -ksp_guess_fischer_monitor
85: TEST*/