Actual source code: ex141.c
2: static char help[] = "Tests converting a SBAIJ matrix to BAIJ format with MatConvert. Modified from ex55.c\n\n";
4: #include <petscmat.h>
5: /* Usage: ./ex141 -mat_view */
7: int main(int argc,char **args)
8: {
9: Mat C,B;
10: PetscInt i,bs=2,mbs,m,block,d_nz=6,col[3];
11: PetscMPIInt size;
12: char file[PETSC_MAX_PATH_LEN];
13: PetscViewer fd;
14: PetscBool equal,loadmat;
15: PetscScalar value[4];
16: PetscInt ridx[2],cidx[2];
18: PetscInitialize(&argc,&args,(char*)0,help);
19: PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&loadmat);
20: MPI_Comm_size(PETSC_COMM_WORLD,&size);
23: /* input matrix C */
24: if (loadmat) {
25: /* Open binary file. Load a sbaij matrix, then destroy the viewer. */
26: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
27: MatCreate(PETSC_COMM_WORLD,&C);
28: MatSetType(C,MATSEQSBAIJ);
29: MatLoad(C,fd);
30: PetscViewerDestroy(&fd);
31: } else { /* Create a sbaij mat with bs>1 */
32: mbs =8;
33: PetscOptionsGetInt(NULL,NULL,"-mbs",&mbs,NULL);
34: m = mbs*bs;
35: MatCreate(PETSC_COMM_WORLD,&C);
36: MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,m,m);
37: MatSetType(C,MATSBAIJ);
38: MatSetFromOptions(C);
39: MatSeqSBAIJSetPreallocation(C,bs,d_nz,NULL);
40: MatSetUp(C);
41: MatSetOption(C,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE);
43: for (block=0; block<mbs; block++) {
44: /* diagonal blocks */
45: value[0] = -1.0; value[1] = 4.0; value[2] = -1.0;
46: for (i=1+block*bs; i<bs-1+block*bs; i++) {
47: col[0] = i-1; col[1] = i; col[2] = i+1;
48: MatSetValues(C,1,&i,3,col,value,INSERT_VALUES);
49: }
50: i = bs - 1+block*bs; col[0] = bs - 2+block*bs; col[1] = bs - 1+block*bs;
52: value[0]=-1.0; value[1]=4.0;
53: MatSetValues(C,1,&i,2,col,value,INSERT_VALUES);
55: i = 0+block*bs; col[0] = 0+block*bs; col[1] = 1+block*bs;
57: value[0]=4.0; value[1] = -1.0;
58: MatSetValues(C,1,&i,2,col,value,INSERT_VALUES);
59: }
60: /* off-diagonal blocks */
61: value[0]=-1.0; value[1] = -0.1; value[2] = 0.0; value[3] = -1.0; /* row-oriented */
62: for (block=0; block<mbs-1; block++) {
63: for (i=0; i<bs; i++) {
64: ridx[i] = block*bs+i; cidx[i] = (block+1)*bs+i;
65: }
66: MatSetValues(C,bs,ridx,bs,cidx,value,INSERT_VALUES);
67: }
68: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
69: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
70: }
72: /* convert C to BAIJ format */
73: MatConvert(C,MATSEQBAIJ,MAT_INITIAL_MATRIX,&B);
74: MatMultEqual(B,C,10,&equal);
77: MatDestroy(&B);
78: MatDestroy(&C);
79: PetscFinalize();
80: return 0;
81: }
83: /*TEST
85: test:
86: output_file: output/ex141.out
88: TEST*/