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: PetscCall(PetscInitialize(&argc, &args, NULL, help));
18: PetscCall(MatCreateSeqDense(PETSC_COMM_SELF, triangle_size, triangle_size, NULL, &A));
19: for (i = 0; i < triangle_size; ++i) PetscCall(MatSetValue(A, i, i, 1.0, INSERT_VALUES));
20: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
21: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
22: PetscCall(MatCreateNest(PETSC_COMM_SELF, 1, NULL, 1, NULL, &A, &A_nest));
23: PetscCall(MatNestSetVecType(A_nest, VECNEST));
25: PetscCall(KSPCreate(PETSC_COMM_SELF, &ksp));
26: PetscCall(KSPSetOperators(ksp, A_nest, A_nest));
27: PetscCall(KSPSetFromOptions(ksp));
28: PetscCall(KSPGetGuess(ksp, &guess));
29: PetscCall(KSPGuessSetUp(guess));
31: for (i = 0; i < 5; ++i) {
32: PetscCall(VecCreateSeq(PETSC_COMM_SELF, triangle_size, &sol));
33: PetscCall(VecCreateSeq(PETSC_COMM_SELF, triangle_size, &rhs));
34: for (j = 0; j < i; ++j) {
35: PetscCall(VecSetValue(sol, j, (PetscScalar)j, INSERT_VALUES));
36: PetscCall(VecSetValue(rhs, j, (PetscScalar)j, INSERT_VALUES));
37: }
38: PetscCall(VecAssemblyBegin(sol));
39: PetscCall(VecAssemblyBegin(rhs));
40: PetscCall(VecAssemblyEnd(sol));
41: PetscCall(VecAssemblyEnd(rhs));
43: PetscCall(VecCreateNest(PETSC_COMM_SELF, 1, NULL, &sol, &sol_nest));
44: PetscCall(VecCreateNest(PETSC_COMM_SELF, 1, NULL, &rhs, &rhs_nest));
46: PetscCall(KSPGuessUpdate(guess, rhs_nest, sol_nest));
48: PetscCall(VecDestroy(&rhs_nest));
49: PetscCall(VecDestroy(&sol_nest));
50: PetscCall(VecDestroy(&rhs));
51: PetscCall(VecDestroy(&sol));
52: }
54: PetscCall(VecCreateSeq(PETSC_COMM_SELF, triangle_size, &sol));
55: PetscCall(VecCreateSeq(PETSC_COMM_SELF, triangle_size, &rhs));
56: PetscCall(VecSetValues(rhs, 5, indices, values, INSERT_VALUES));
57: PetscCall(VecAssemblyBegin(rhs));
58: PetscCall(VecAssemblyEnd(rhs));
60: PetscCall(VecCreateNest(PETSC_COMM_SELF, 1, NULL, &sol, &sol_nest));
61: PetscCall(VecCreateNest(PETSC_COMM_SELF, 1, NULL, &rhs, &rhs_nest));
63: PetscCall(KSPGuessFormGuess(guess, rhs_nest, sol_nest));
64: PetscCall(VecView(sol_nest, PETSC_VIEWER_STDOUT_SELF));
66: PetscCall(VecDestroy(&rhs_nest));
67: PetscCall(VecDestroy(&sol_nest));
68: PetscCall(VecDestroy(&rhs));
69: PetscCall(VecDestroy(&sol));
71: PetscCall(KSPDestroy(&ksp));
73: PetscCall(MatDestroy(&A_nest));
74: PetscCall(MatDestroy(&A));
76: PetscCall(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*/