Actual source code: ex34.c
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;
9: PetscInt N = 10, M = 3;
10: PetscBool no_inodes=PETSC_TRUE,flg;
11: KSP ksp;
12: PC pc;
13: Vec x,y;
14: char mtype[256];
16: PetscInitialize(&argc,&args,(char*)0,help);
17: PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);
18: PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);
19: PetscOptionsGetString(NULL,NULL,"-mtype",mtype,sizeof(mtype),&flg);
20: PetscOptionsGetBool(NULL,NULL,"-no_inodes",&no_inodes,NULL);
21: MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,M,N,NULL,&Ad);
22: MatSetRandom(Ad,NULL);
23: MatConvert(Ad,flg ? mtype : MATAIJ,MAT_INITIAL_MATRIX,&A);
24: MatProductCreate(A,A,NULL,&B);
25: MatProductSetType(B,MATPRODUCT_AtB);
26: MatProductSetAlgorithm(B,"default");
27: MatProductSetFill(B,PETSC_DEFAULT);
28: MatProductSetFromOptions(B);
29: MatProductSymbolic(B);
30: if (no_inodes) {
31: MatSetOption(B,MAT_USE_INODES,PETSC_FALSE);
32: }
33: MatProductNumeric(B);
34: MatTransposeMatMultEqual(A,A,B,10,&flg);
35: if (!flg) {
36: PetscPrintf(PETSC_COMM_WORLD,"Wrong MatTransposeMat");
37: }
38: KSPCreate(PETSC_COMM_WORLD,&ksp);
39: KSPSetOperators(ksp,B,B);
40: KSPGetPC(ksp,&pc);
41: PCSetType(pc,PCSOR);
42: PCSORSetOmega(pc,1.1);
43: KSPSetUp(ksp);
44: KSPView(ksp,NULL);
45: KSPGetPC(ksp,&pc);
46: MatCreateVecs(B,&y,&x);
47: VecSetRandom(x,NULL);
48: PCApply(pc,x,y);
49: KSPDestroy(&ksp);
50: VecDestroy(&x);
51: VecDestroy(&y);
52: MatDestroy(&B);
53: MatDestroy(&Ad);
54: MatDestroy(&A);
55: PetscFinalize();
56: return 0;
57: }
59: /*TEST
61: test:
62: nsize: 1
63: suffix: 1
65: test:
66: nsize: 1
67: suffix: 1_mpiaij
68: args: -mtype mpiaij
70: test:
71: nsize: 3
72: suffix: 2
74: TEST*/