[petsc-users] How to simplify the nonzero structure of the jacobian matrix.

Matthew Knepley knepley at gmail.com
Fri Jan 20 19:01:15 CST 2012


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?

   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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120120/ff2915eb/attachment.htm>


More information about the petsc-users mailing list