[petsc-users] How to simplify the nonzero structure of the jacobian matrix.
Matthew Knepley
knepley at gmail.com
Fri Jan 20 18:09:02 CST 2012
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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120120/2dcaefea/attachment.htm>
More information about the petsc-users
mailing list