Actual source code: ex185.c
1: static char help[] = "Tests MatCreateConstantDiagonal().\n"
2: "\n";
4: #include <petscmat.h>
6: /*T
7: Concepts: Mat
8: T*/
10: int main(int argc, char **args)
11: {
12: Vec X, Y;
13: Mat A,B,Af;
14: PetscBool flg;
15: PetscReal xnorm,ynorm,anorm;
17: PetscInitialize(&argc,&args,(char*)0,help);
19: MatCreateConstantDiagonal(PETSC_COMM_WORLD,PETSC_DETERMINE,PETSC_DETERMINE,20,20,3.0,&A);
20: MatCreateVecs(A,&X,&Y);
21: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
23: VecSetRandom(X,NULL);
24: VecNorm(X,NORM_2,&xnorm);
25: MatMult(A,X,Y);
26: VecNorm(Y,NORM_2,&ynorm);
28: MatShift(A,5.0);
29: MatScale(A,.5);
30: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
31: MatNorm(A,NORM_FROBENIUS,&anorm);
34: /* Convert to AIJ (exercises MatGetRow/MatRestoreRow) */
35: MatConvert(A,MATAIJ,MAT_INITIAL_MATRIX,&B);
36: MatMultEqual(A,B,10,&flg);
37: if (!flg) PetscPrintf(PETSC_COMM_WORLD,"Error MatMult\n");
38: MatMultAddEqual(A,B,10,&flg);
39: if (!flg) PetscPrintf(PETSC_COMM_WORLD,"Error MatMultAdd\n");
40: MatMultTransposeEqual(A,B,10,&flg);
41: if (!flg) PetscPrintf(PETSC_COMM_WORLD,"Error MatMultTranspose\n");
42: MatMultTransposeAddEqual(A,B,10,&flg);
43: if (!flg) PetscPrintf(PETSC_COMM_WORLD,"Error MatMultTransposeAdd\n");
45: MatGetDiagonal(A,Y);
46: MatGetFactor(A,MATSOLVERPETSC,MAT_FACTOR_LU,&Af);
47: MatLUFactorSymbolic(Af,A,NULL,NULL,NULL);
48: MatLUFactorNumeric(Af,A,NULL);
49: MatSolve(Af,X,Y);
50: VecNorm(Y,NORM_2,&ynorm);
53: MatDestroy(&A);
54: MatDestroy(&B);
55: MatDestroy(&Af);
56: VecDestroy(&X);
57: VecDestroy(&Y);
59: PetscFinalize();
60: return 0;
61: }
63: /*TEST
65: test:
66: nsize: 2
68: TEST*/