16679dcc1SBarry Smith static char help[] = "Creates MatSeqBAIJ matrix of given BS for timing tests of MatMult().\n"; 26679dcc1SBarry Smith 36679dcc1SBarry Smith #include <petscmat.h> 46679dcc1SBarry Smith 56679dcc1SBarry Smith int main(int argc,char **args) 66679dcc1SBarry Smith { 76679dcc1SBarry Smith Mat A; 86679dcc1SBarry Smith Vec x,y; 96679dcc1SBarry Smith PetscErrorCode ierr; 106679dcc1SBarry Smith PetscInt m=50000,bs=12,i,j,k,l,row,col,M; 116679dcc1SBarry Smith PetscScalar rval,*vals; 126679dcc1SBarry Smith PetscRandom rdm; 136679dcc1SBarry Smith 146679dcc1SBarry Smith ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 156679dcc1SBarry Smith ierr = PetscOptionsGetInt(NULL,NULL,"-mat_block_size",&bs,NULL);CHKERRQ(ierr); 166679dcc1SBarry Smith ierr = PetscOptionsGetInt(NULL,NULL,"-mat_size",&m,NULL);CHKERRQ(ierr); 176679dcc1SBarry Smith M = m*bs; 186679dcc1SBarry Smith ierr = MatCreateSeqBAIJ(PETSC_COMM_SELF,bs,M,M,27,NULL,&A);CHKERRQ(ierr); 19e269983cSBarry Smith ierr = MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);CHKERRQ(ierr); 206679dcc1SBarry Smith 216679dcc1SBarry Smith ierr = PetscRandomCreate(PETSC_COMM_SELF,&rdm);CHKERRQ(ierr); 226679dcc1SBarry Smith ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr); 236679dcc1SBarry Smith ierr = VecCreateSeq(PETSC_COMM_SELF,M,&x);CHKERRQ(ierr); 246679dcc1SBarry Smith ierr = VecDuplicate(x,&y);CHKERRQ(ierr); 256679dcc1SBarry Smith 26e269983cSBarry Smith /* For each block row insert at most 27 blocks */ 276679dcc1SBarry Smith ierr = PetscMalloc1(bs*bs,&vals);CHKERRQ(ierr); 286679dcc1SBarry Smith for (i=0; i<m; i++) { 296679dcc1SBarry Smith row = i; 306679dcc1SBarry Smith for (j=0; j<27; j++) { 316679dcc1SBarry Smith ierr = PetscRandomGetValue(rdm,&rval);CHKERRQ(ierr); 326679dcc1SBarry Smith col = (PetscInt)(PetscRealPart(rval)*m); 336679dcc1SBarry Smith for (k=0; k<bs; k++) { 346679dcc1SBarry Smith for (l=0; l<bs; l++) { 356679dcc1SBarry Smith ierr = PetscRandomGetValue(rdm,&rval);CHKERRQ(ierr); 366679dcc1SBarry Smith vals[k*bs + l] = rval; 376679dcc1SBarry Smith } 386679dcc1SBarry Smith } 396679dcc1SBarry Smith ierr = MatSetValuesBlocked(A,1,&row,1,&col,vals,INSERT_VALUES);CHKERRQ(ierr); 406679dcc1SBarry Smith } 416679dcc1SBarry Smith } 426679dcc1SBarry Smith ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 436679dcc1SBarry Smith ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 446679dcc1SBarry Smith ierr = PetscFree(vals);CHKERRQ(ierr); 456679dcc1SBarry Smith 466679dcc1SBarry Smith /* Time MatMult(), MatMultAdd() */ 476679dcc1SBarry Smith for (i=0; i<25; i++) { 486679dcc1SBarry Smith ierr = VecSetRandom(x,rdm);CHKERRQ(ierr); 496679dcc1SBarry Smith ierr = MatMult(A,x,y);CHKERRQ(ierr); 506679dcc1SBarry Smith ierr = VecSetRandom(x,rdm);CHKERRQ(ierr); 516679dcc1SBarry Smith ierr = VecSetRandom(y,rdm);CHKERRQ(ierr); 526679dcc1SBarry Smith ierr = MatMultAdd(A,x,y,y);CHKERRQ(ierr); 536679dcc1SBarry Smith } 546679dcc1SBarry Smith 556679dcc1SBarry Smith ierr = MatDestroy(&A);CHKERRQ(ierr); 566679dcc1SBarry Smith ierr = VecDestroy(&x);CHKERRQ(ierr); 576679dcc1SBarry Smith ierr = VecDestroy(&y);CHKERRQ(ierr); 586679dcc1SBarry Smith ierr = PetscRandomDestroy(&rdm);CHKERRQ(ierr); 596679dcc1SBarry Smith ierr = PetscFinalize(); 606679dcc1SBarry Smith return ierr; 616679dcc1SBarry Smith } 626679dcc1SBarry Smith 636679dcc1SBarry Smith 646679dcc1SBarry Smith /*TEST 656679dcc1SBarry Smith 66*cc384ab3SSatish Balay testset: 67041d8c91SPierre Jolivet requires: define(PETSC_USING_64BIT_PTR) 68*cc384ab3SSatish Balay output_file: output/ex238_1.out 69*cc384ab3SSatish Balay test: 70*cc384ab3SSatish Balay suffix: 1 71*cc384ab3SSatish Balay args: -mat_block_size 1 72*cc384ab3SSatish Balay test: 73*cc384ab3SSatish Balay suffix: 2 74*cc384ab3SSatish Balay args: -mat_block_size 2 75*cc384ab3SSatish Balay test: 76*cc384ab3SSatish Balay suffix: 4 77*cc384ab3SSatish Balay args: -mat_block_size 4 78*cc384ab3SSatish Balay test: 79*cc384ab3SSatish Balay suffix: 5 80*cc384ab3SSatish Balay args: -mat_block_size 5 81*cc384ab3SSatish Balay test: 82*cc384ab3SSatish Balay suffix: 6 83*cc384ab3SSatish Balay args: -mat_block_size 6 84*cc384ab3SSatish Balay test: 85*cc384ab3SSatish Balay suffix: 8 86*cc384ab3SSatish Balay args: -mat_block_size 8 87*cc384ab3SSatish Balay test: 88*cc384ab3SSatish Balay suffix: 12 89*cc384ab3SSatish Balay args: -mat_block_size 12 90*cc384ab3SSatish Balay test: 91*cc384ab3SSatish Balay suffix: 15 92*cc384ab3SSatish Balay args: -mat_block_size 15 936679dcc1SBarry Smith 946679dcc1SBarry Smith TEST*/ 95