-dm_preallocate_only http://www.mcs.anl.gov/petsc/petsc-dev/docs/manualpages/DM/DMSetMatrixPreallocateOnly.html
On Mon, Jan 23, 2012 at 15:52, Xuefei (Rebecca) Yuan <xyuan at lbl.gov> wrote: > Hello all, > > There is one example that has used the > > ierr = MatSetOption(*B,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE);CHKERRQ(ierr); > > in /src/mat/examples/tutorials/ex12.c, > > and this option works for this case because the matrix B is created via > > ierr = MatCreateSeqAIJ(PETSC_COMM_SELF,n+1,n+1,0,cnt,B);CHKERRQ(ierr); > > and the values are inserted by > > ierr = MatSetValues(*B,1,&i,aij->i[i+1] - aij->i[i],aij->j + > aij->i[i],aij->a + aij->i[i],INSERT_VALUES);CHKERRQ(ierr); > > > However, the option MAT_IGNORE_ZERO_ENTRIES does not work for the matrix > if it is generated via > > ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr); > > no matter if MatSetValues() or MatSetValuesStencil() is used. > > In the other word, if the matrix is generated based on the stencil > information, the zero entries cannot be ignored even if the option is > called before inserting/adding values. > > Is there anyway that can get rid of those zeros in the matrix generated > based on the stencil information? > > > Thanks very much! > > Best regards, > > Rebecca > > > > > > On Jan 20, 2012, at 5:23 PM, Matthew Knepley wrote: > > On Fri, Jan 20, 2012 at 7:07 PM, Xuefei (Rebecca) Yuan <xyuan at > lbl.gov>wrote: > >> Here is the output: >> >> >> >> >> >> On Jan 20, 2012, at 5:01 PM, Matthew Knepley wrote: >> >> On Fri, Jan 20, 2012 at 6:55 PM, Xuefei (Rebecca) Yuan <xyuan at >> lbl.gov>wrote: >> >>> Hello Matt, >>> >>> I tried several times for 3.1-p8 and dev version by putting MatSetOption >>> >> >> Are you sure your entries are exactly 0.0? >> >> > Are you using ADD_VALUES? > > > http://petsc.cs.iit.edu/petsc/petsc-dev/file/783e93230143/src/mat/impls/aij/seq/aij.c#l310 > > Matt > > >> Matt >> >> >>> 1) right after creation of the matrix: >>> >>> #ifdef petscDev >>> ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, >>> &jacobian);CHKERRQ(ierr); >>> #else >>> ierr = DAGetMatrix(DMMGGetDA(dmmg), MATAIJ, >>> &jacobian);CHKERRQ(ierr); >>> #endif >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); >>> >>> 2) at the beginning of the FormJacobianLocal() routine: >>> >>> PetscFunctionBegin; >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); >>> >>> 3) before calling MatAssemblyBegin() in FormJacobianLocal() routine: >>> >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); >>> ierr = MatAssemblyBegin(jacobian, >>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>> ierr = MatAssemblyEnd(jacobian, >>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>> >>> None of those works. What is wrong here? >>> >>> Thanks, >>> >>> R >>> >>> >>> On Jan 20, 2012, at 4:32 PM, Matthew Knepley wrote: >>> >>> On Fri, Jan 20, 2012 at 6:28 PM, Xuefei (Rebecca) Yuan <xyuan at >>> lbl.gov>wrote: >>> >>>> Hello Matt, >>>> >>>> I have changed the code as >>>> >>>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>>> PETSC_TRUE);CHKERRQ(ierr); >>>> ierr = MatAssemblyBegin(jacobian, >>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> ierr = MatAssemblyEnd(jacobian, >>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> >>> >>> You have to set it before you start setting values, so we know to ignore >>> them. >>> >>> Matt >>> >>> >>>> but still get the same result as before, the matrix still has 5776 >>>> nonzeros: >>>> >>>> % Size = 100 100 >>>> 2 % Nonzeros = 5776 >>>> 3 zzz = zeros(5776,3); >>>> >>>> Then I switch the order as >>>> >>>> ierr = MatAssemblyBegin(jacobian, >>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> ierr = MatAssemblyEnd(jacobian, >>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>>> PETSC_TRUE);CHKERRQ(ierr); >>>> >>>> nothing changed. >>>> >>>> The version is 3.1-p8. >>>> >>>> Thanks very much! >>>> >>>> Best regards, >>>> >>>> Rebecca >>>> >>>> >>>> >>>> >>>> On Jan 20, 2012, at 4:09 PM, Matthew Knepley wrote: >>>> >>>> On Fri, Jan 20, 2012 at 6:02 PM, Xuefei (Rebecca) Yuan <xyuan at >>>> lbl.gov>wrote: >>>> >>>>> Hello all, >>>>> >>>>> This is a test for np=1 case of the nonzero structure of the jacobian >>>>> matrix. The jacobian matrix is created via >>>>> >>>>> ierr = >>>>> DMDACreate2d(comm,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_STENCIL_BOX, >>>>> parameters.mxfield, parameters.myfield, PETSC_DECIDE, PETSC_DECIDE, 4, 2, >>>>> 0, 0, &da);CHKERRQ(ierr); >>>>> >>>>> ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, >>>>> &jacobian);CHKERRQ(ierr); >>>>> >>>>> After creation of the jacobian matrix, >>>>> >>>>> ierr = MatAssemblyBegin(jacobian, >>>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>>> ierr = MatAssemblyEnd(jacobian, >>>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>>> >>>>> PetscViewer viewer; >>>>> char fileName[120]; >>>>> sprintf(fileName, >>>>> "jacobian_after_creation.m");CHKERRQ(ierr); >>>>> >>>>> FILE * fp; >>>>> >>>>> ierr = >>>>> PetscViewerASCIIOpen(PETSC_COMM_WORLD,fileName,&viewer);CHKERRQ(ierr); >>>>> ierr = >>>>> PetscViewerSetFormat(viewer,PETSC_VIEWER_ASCII_MATLAB);CHKERRQ(ierr); >>>>> ierr = MatView (jacobian, viewer); CHKERRQ (ierr); >>>>> ierr = PetscFOpen(PETSC_COMM_WORLD,fileName,"a",&fp); >>>>> CHKERRQ(ierr); >>>>> ierr = >>>>> PetscViewerASCIIPrintf(viewer,"spy((spconvert(zzz)));\n");CHKERRQ(ierr); >>>>> ierr = PetscFClose(PETSC_COMM_WORLD,fp);CHKERRQ(ierr); >>>>> PetscViewerDestroy(&viewer); >>>>> >>>>> I took a look at the structure of the jacobian by storing it in the >>>>> matlab format, the matrix has 5776 nonzeros entries, however, those values >>>>> are all zeros at the moment as I have not insert or add any values into it >>>>> yet, the structure shows: (the following figure shows a global replacement >>>>> of 0.0 by 1.0 for those 5776 numbers) >>>>> >>>>> >>>>> >>>>> >>>>> Inside the FormJacobianLocal() function, I have selected the index to >>>>> pass to the nonzero values to jacobian, for example, >>>>> >>>>> ierr = MatSetValuesStencil(jacobian, 1, &row, 6, col, v, >>>>> INSERT_VALUES);CHKERRQ(ierr); >>>>> >>>>> where >>>>> >>>>> col[0].i = column[4].i; >>>>> col[1].i = column[5].i; >>>>> col[2].i = column[6].i; >>>>> col[3].i = column[9].i; >>>>> col[4].i = column[10].i; >>>>> col[5].i = column[12].i; >>>>> >>>>> >>>>> col[0].j = column[4].j; >>>>> col[1].j = column[5].j; >>>>> col[2].j = column[6].j; >>>>> col[3].j = column[9].j; >>>>> col[4].j = column[10].j; >>>>> col[5].j = column[12].j; >>>>> >>>>> col[0].c = column[4].c; >>>>> col[1].c = column[5].c; >>>>> col[2].c = column[6].c; >>>>> col[3].c = column[9].c; >>>>> col[4].c = column[10].c; >>>>> col[5].c = column[12].c; >>>>> >>>>> v[0] = value[4]; >>>>> v[1] = value[5]; >>>>> v[2] = value[6]; >>>>> v[3] = value[9]; >>>>> v[4] = value[10]; >>>>> v[5] = value[12]; >>>>> >>>>> and did not pass the zero entries into the jacobian matrix. However, >>>>> after inserting or adding all values to the matrix, by the same routine >>>>> above to take a look at the jacobian matrix in matlab format, the matrix >>>>> still has 5776 nonzeros, in which 1075 numbers are nonzeros, and the other >>>>> 4701 numbers are all zeros. The spy() gives >>>>> >>>>> >>>>> >>>>> >>>>> for the true nonzero structures. >>>>> >>>>> But the ksp_view will give the nonzeros number as 5776, instead of >>>>> 1075: >>>>> >>>>> linear system matrix = precond matrix: >>>>> Matrix Object: Mat_0x84000000_1 1 MPI processes >>>>> type: seqaij >>>>> rows=100, cols=100 >>>>> total: nonzeros=5776, allocated nonzeros=5776 >>>>> >>>>> It is a waste of memory to have all those values of zeros been stored >>>>> in the jacobian. >>>>> >>>>> Is there anyway to get rid of those zero values in jacobian and has >>>>> the only nonzero numbers stored in jacobian? In such a case, the ksp_view >>>>> will tell that total: nonzeros=1075. >>>>> >>>> >>>> MatSetOption(MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE); >>>> >>>> Matt >>>> >>>> >>>>> Thanks very much! >>>>> >>>>> Have a nice weekend! >>>>> >>>>> Cheers, >>>>> >>>>> Rebecca >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>> >>>> >>>> -- >>>> What most experimenters take for granted before they begin their >>>> experiments is infinitely more interesting than any results to which their >>>> experiments lead. >>>> -- Norbert Wiener >>>> >>>> >>>> >>> >>> >>> -- >>> What most experimenters take for granted before they begin their >>> experiments is infinitely more interesting than any results to which their >>> experiments lead. >>> -- Norbert Wiener >>> >>> >>> >> >> >> -- >> What most experimenters take for granted before they begin their >> experiments is infinitely more interesting than any results to which their >> experiments lead. >> -- Norbert Wiener >> >> >> >> > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > > -------------- next part -------------- An HTML attachment was scrubbed... URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120123/786b85f9/attachment.htm>
