Hello,
After I was not able to get rid of the zeros in the matrix by
>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES,
>> PETSC_TRUE);CHKERRQ(ierr);
I would like to create the matrix by the following routine:
ierr = MatCreate(PETSC_COMM_WORLD, &jacobian);CHKERRQ(ierr);
ierr = MatSetType(jacobian, MATMPIAIJ);CHKERRQ(ierr);
ierr = MatSetSizes(jacobian, PETSC_DECIDE, PETSC_DECIDE,
(PetscInt)(info.mx*info.my*4), (PetscInt)(info.mx*info.my*4));CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(jacobian, 11, PETSC_NULL, 18,
PETSC_NULL);CHKERRQ(ierr);
instead of using
ierr = DMGetMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr);
to save the memory. However, I found that the routine DMGetMatrix() was still
called inside DMMGSetSNESLocal(), therefore, the matrix I created was useless:
*********************************************
Breakpoint 2, main (argc=3, argv=0x7fff5fbff770) at twcartffxmhd.c:255
255 ierr = MatCreate(PETSC_COMM_WORLD, &jacobian);CHKERRQ(ierr);
(gdb) n
256 ierr = MatSetType(jacobian, MATMPIAIJ);CHKERRQ(ierr);
(gdb)
257 ierr = MatSetSizes(jacobian, PETSC_DECIDE, PETSC_DECIDE,
(PetscInt)(info.mx*info.my*4), (PetscInt)(info.mx*info.my*4));CHKERRQ(ierr);
(gdb)
258 ierr = MatMPIAIJSetPreallocation(jacobian, 11, PETSC_NULL, 18,
PETSC_NULL);CHKERRQ(ierr);
(gdb)
Breakpoint 3, main (argc=3, argv=0x7fff5fbff770) at twcartffxmhd.c:260
260 ierr = DMMGSetSNESLocal(dmmg, FormFunctionLocal,
FormJacobianLocal,0,0);CHKERRQ(ierr);
(gdb) s
DMMGSetSNESLocal_Private (dmmg=0x102004480, function=0x100016ba9
<FormFunctionLocal>, jacobian=0x100050e83 <FormJacobianLocal>, ad_function=0,
admf_function=0) at damgsnes.c:932
932 PetscInt i,nlevels = dmmg[0]->nlevels;
(gdb) n
934 PetscErrorCode
(*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*) = 0;
(gdb) s
937 PetscFunctionBegin;
(gdb) n
938 if (jacobian) computejacobian = DMMGComputeJacobian;
(gdb)
942 CHKMEMQ;
(gdb)
943 ierr = PetscObjectGetCookie((PetscObject)
dmmg[0]->dm,&cookie);CHKERRQ(ierr);
(gdb)
944 if (cookie == DM_COOKIE) {
(gdb)
948 ierr = PetscOptionsHasName(PETSC_NULL, "-dmmg_form_function_ghost",
&flag);CHKERRQ(ierr);
(gdb)
949 if (flag) {
(gdb)
952 ierr =
DMMGSetSNES(dmmg,DMMGFormFunction,computejacobian);CHKERRQ(ierr);
(gdb)
Matrix Object:
type=seqaij, rows=100, cols=100
total: nonzeros=5776, allocated nonzeros=5776
using I-node routines: found 25 nodes, limit used is 5
954 for (i=0; i<nlevels; i++) {
(gdb) b twcartffxmhd.c:280
Breakpoint 4 at 0x100004032: file twcartffxmhd.c, line 282.
(gdb) c
Continuing.
Matrix Object:
type=mpiaij, rows=100, cols=100
total: nonzeros=0, allocated nonzeros=2900
using I-node (on process 0) routines: found 20 nodes, limit used is 5
*********************************************
The final matrix after FormJacobianlLocal() call still have 5776 nonzeros:
% Size = 100 100
2 % Nonzeros = 5776
3 zzz = zeros(5776,3);
Is there a way that I can pass the matrix created via line 255-258 to
FormJacobianLocal()?
Thanks very much!
Cheers,
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/20120121/149cbb2d/attachment.htm>
