Actual source code: ex34.c
petsc-3.13.4 2020-08-01
2: static char help[] = "Tests solving linear system with KSPFGMRES + PCSOR (omega != 1) on a matrix obtained from MatTransposeMatMult.\n\n";
4: #include <petscksp.h>
6: int main(int argc,char **args)
7: {
8: Mat A,Ad,B;
10: PetscInt N = 10, M = 3;
11: PetscBool no_inodes=PETSC_TRUE,flg;
12: KSP ksp;
13: PC pc;
14: Vec x,y;
15: char mtype[256];
17: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
18: PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);
19: PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);
20: PetscOptionsGetString(NULL,NULL,"-mtype",mtype,256,&flg);
21: PetscOptionsGetBool(NULL,NULL,"-no_inodes",&no_inodes,NULL);
22: MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,M,N,NULL,&Ad);
23: MatSetRandom(Ad,NULL);
24: MatConvert(Ad,flg ? mtype : MATAIJ,MAT_INITIAL_MATRIX,&A);
25: MatProductCreate(A,A,NULL,&B);
26: MatProductSetType(B,MATPRODUCT_AtB);
27: MatProductSetAlgorithm(B,"default");
28: MatProductSetFill(B,PETSC_DEFAULT);
29: MatProductSetFromOptions(B);
30: MatProductSymbolic(B);
31: if (no_inodes) {
32: MatSetOption(B,MAT_USE_INODES,PETSC_FALSE);
33: }
34: MatProductNumeric(B);
35: MatTransposeMatMultEqual(A,A,B,10,&flg);
36: if (!flg) {
37: PetscPrintf(PETSC_COMM_WORLD,"Wrong MatTransposeMat");
38: }
39: KSPCreate(PETSC_COMM_WORLD,&ksp);
40: KSPSetOperators(ksp,B,B);
41: KSPGetPC(ksp,&pc);
42: PCSetType(pc,PCSOR);
43: PCSORSetOmega(pc,1.1);
44: KSPSetUp(ksp);
45: KSPView(ksp,NULL);
46: KSPGetPC(ksp,&pc);
47: MatCreateVecs(B,&y,&x);
48: VecSetRandom(x,NULL);
49: PCApply(pc,x,y);
50: KSPDestroy(&ksp);
51: VecDestroy(&x);
52: VecDestroy(&y);
53: MatDestroy(&B);
54: MatDestroy(&Ad);
55: MatDestroy(&A);
56: PetscFinalize();
57: return ierr;
58: }
60: /*TEST
62: test:
63: nsize: 1
64: suffix: 1
66: test:
67: nsize: 1
68: suffix: 1_mpiaij
69: args: -mtype mpiaij
71: test:
72: nsize: 3
73: suffix: 2
75: TEST*/