[petsc-users] How to simplify the nonzero structure of the jacobian matrix.
Matthew Knepley
knepley at gmail.com
Fri Jan 20 19:23:43 CST 2012
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/20120120/61e4f2e2/attachment.htm>
More information about the petsc-users
mailing list