<html><head></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; ">Hello,<div><br></div><div>After I was not able to get rid of the zeros in the matrix by </div><div><br><div>
<span class="Apple-style-span" style="border-collapse: separate; color: rgb(0, 0, 0); font-family: Helvetica; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: 2; text-align: auto; text-indent: 0px; text-transform: none; white-space: normal; widows: 2; word-spacing: 0px; -webkit-border-horizontal-spacing: 0px; -webkit-border-vertical-spacing: 0px; -webkit-text-decorations-in-effect: none; -webkit-text-size-adjust: auto; -webkit-text-stroke-width: 0px; font-size: medium; "><span class="Apple-style-span" style="border-collapse: separate; color: rgb(0, 0, 0); font-family: Helvetica; font-style: normal; font-variant: normal; font-weight: normal; letter-spacing: normal; line-height: normal; orphans: 2; text-indent: 0px; text-transform: none; white-space: normal; widows: 2; word-spacing: 0px; -webkit-border-horizontal-spacing: 0px; -webkit-border-vertical-spacing: 0px; -webkit-text-decorations-in-effect: none; -webkit-text-size-adjust: auto; -webkit-text-stroke-width: 0px; font-size: medium; "><div style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; "><blockquote type="cite"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0.8ex; border-left-width: 1px; border-left-color: rgb(204, 204, 204); border-left-style: solid; padding-left: 1ex; position: static; z-index: auto; "><div style="word-wrap: break-word; "><div><div><div><blockquote type="cite"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0.8ex; border-left-width: 1px; border-left-color: rgb(204, 204, 204); border-left-style: solid; padding-left: 1ex; position: static; z-index: auto; "><div style="word-wrap: break-word; "><div><div>ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);CHKERRQ(ierr);</div></div></div></blockquote></div></blockquote></div></div></div></div></blockquote></div></blockquote><br><br></div></span></span>
</div><div>I would like to create the matrix by the following routine:</div><div><br></div><div><div><span class="Apple-tab-span" style="white-space:pre"> </span></div><div><span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatCreate(PETSC_COMM_WORLD, &jacobian);CHKERRQ(ierr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatSetType(jacobian, MATMPIAIJ);CHKERRQ(ierr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatSetSizes(jacobian, PETSC_DECIDE, PETSC_DECIDE, (PetscInt)(info.mx*info.my*4), (PetscInt)(info.mx*info.my*4));CHKERRQ(ierr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatMPIAIJSetPreallocation(jacobian, 11, PETSC_NULL, 18, PETSC_NULL);CHKERRQ(ierr);</div></div><div><br></div><div>instead of using </div><div><br></div><div><div>ierr = DMGetMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr);</div></div><div><br></div><div>to save the memory. However, I found that the routine DMGetMatrix() was still called inside DMMGSetSNESLocal(), therefore, the matrix I created was useless:</div><div><br></div><div><div>*********************************************</div><div><br></div><div>Breakpoint 2, main (argc=3, argv=0x7fff5fbff770) at twcartffxmhd.c:255</div><div>255<span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatCreate(PETSC_COMM_WORLD, &jacobian);CHKERRQ(ierr);</div><div>(gdb) n</div><div>256<span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatSetType(jacobian, MATMPIAIJ);CHKERRQ(ierr);</div><div>(gdb) </div><div>257<span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatSetSizes(jacobian, PETSC_DECIDE, PETSC_DECIDE, (PetscInt)(info.mx*info.my*4), (PetscInt)(info.mx*info.my*4));CHKERRQ(ierr);</div><div>(gdb) </div><div>258<span class="Apple-tab-span" style="white-space:pre"> </span>ierr = MatMPIAIJSetPreallocation(jacobian, 11, PETSC_NULL, 18, PETSC_NULL);CHKERRQ(ierr);</div></div><div><div>(gdb) </div><div><br></div><div>Breakpoint 3, main (argc=3, argv=0x7fff5fbff770) at twcartffxmhd.c:260</div><div>260<span class="Apple-tab-span" style="white-space:pre"> </span>ierr = DMMGSetSNESLocal(dmmg, FormFunctionLocal, FormJacobianLocal,0,0);CHKERRQ(ierr);</div><div>(gdb) s</div><div>DMMGSetSNESLocal_Private (dmmg=0x102004480, function=0x100016ba9 <FormFunctionLocal>, jacobian=0x100050e83 <FormJacobianLocal>, ad_function=0, admf_function=0) at damgsnes.c:932</div><div>932<span class="Apple-tab-span" style="white-space:pre"> </span> PetscInt i,nlevels = dmmg[0]->nlevels;</div><div>(gdb) n</div><div>934<span class="Apple-tab-span" style="white-space:pre"> </span> PetscErrorCode (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*) = 0;</div><div>(gdb) s</div><div>937<span class="Apple-tab-span" style="white-space:pre"> </span> PetscFunctionBegin;</div><div>(gdb) n</div><div>938<span class="Apple-tab-span" style="white-space:pre"> </span> if (jacobian) computejacobian = DMMGComputeJacobian;</div><div>(gdb) </div><div>942<span class="Apple-tab-span" style="white-space:pre"> </span> CHKMEMQ;</div><div>(gdb) </div><div>943<span class="Apple-tab-span" style="white-space:pre"> </span> ierr = PetscObjectGetCookie((PetscObject) dmmg[0]->dm,&cookie);CHKERRQ(ierr);</div><div>(gdb) </div><div>944<span class="Apple-tab-span" style="white-space:pre"> </span> if (cookie == DM_COOKIE) {</div><div>(gdb) </div></div><div><div>948<span class="Apple-tab-span" style="white-space:pre"> </span> ierr = PetscOptionsHasName(PETSC_NULL, "-dmmg_form_function_ghost", &flag);CHKERRQ(ierr);</div><div>(gdb) </div><div>949<span class="Apple-tab-span" style="white-space:pre"> </span> if (flag) {</div><div>(gdb) </div><div>952<span class="Apple-tab-span" style="white-space:pre"> </span> ierr = DMMGSetSNES(dmmg,DMMGFormFunction,computejacobian);CHKERRQ(ierr);</div><div>(gdb) </div><div>Matrix Object:</div><div> type=seqaij, rows=100, cols=100</div><div> total: nonzeros=5776, allocated nonzeros=5776</div><div> using I-node routines: found 25 nodes, limit used is 5</div><div>954<span class="Apple-tab-span" style="white-space:pre"> </span> for (i=0; i<nlevels; i++) {</div><div>(gdb) b twcartffxmhd.c:280</div><div>Breakpoint 4 at 0x100004032: file twcartffxmhd.c, line 282.</div><div>(gdb) c</div><div>Continuing.</div><div>Matrix Object:</div><div> type=mpiaij, rows=100, cols=100</div><div> total: nonzeros=0, allocated nonzeros=2900</div><div> using I-node (on process 0) routines: found 20 nodes, limit used is 5</div></div><div><br></div><div><div>*********************************************</div></div><div><br></div><div>The final matrix after FormJacobianlLocal() call still have 5776 nonzeros:</div><div><br></div><div><div>% Size = 100 100</div><div> 2 % Nonzeros = 5776</div><div> 3 zzz = zeros(5776,3);</div></div><div><br></div><div>Is there a way that I can pass the matrix created via line 255-258 to FormJacobianLocal()?</div><div><br></div><div>Thanks very much!</div><div><br></div><div>Cheers,</div><div><br></div><div>Rebecca</div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div>
<br><div><div>On Jan 20, 2012, at 5:23 PM, Matthew Knepley wrote:</div><br class="Apple-interchange-newline"><blockquote type="cite">On Fri, Jan 20, 2012 at 7:07 PM, Xuefei (Rebecca) Yuan <span dir="ltr"><<a href="mailto:xyuan@lbl.gov">xyuan@lbl.gov</a>></span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div style="word-wrap:break-word">Here is the output:<div><br></div><div></div></div><br><div style="word-wrap:break-word"><div></div><div><br><div><br></div><div><div><div style="word-wrap:break-word"><br></div></div><div>
<div>On Jan 20, 2012, at 5:01 PM, Matthew Knepley wrote:</div><br><blockquote type="cite">On Fri, Jan 20, 2012 at 6:55 PM, Xuefei (Rebecca) Yuan <span dir="ltr"><<a href="mailto:xyuan@lbl.gov" target="_blank">xyuan@lbl.gov</a>></span> wrote:<br>
<div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div style="word-wrap:break-word"><div>Hello Matt,</div><div><br></div><div>I tried several times for 3.1-p8 and dev version by putting MatSetOption </div></div></blockquote><div><br></div><div>Are you sure your entries are exactly 0.0?</div>
</div></blockquote></div></div></div></div></blockquote><div><br></div><div>Are you using ADD_VALUES?</div><div><br></div><div><a href="http://petsc.cs.iit.edu/petsc/petsc-dev/file/783e93230143/src/mat/impls/aij/seq/aij.c#l310">http://petsc.cs.iit.edu/petsc/petsc-dev/file/783e93230143/src/mat/impls/aij/seq/aij.c#l310</a></div>
<div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0.8ex; border-left-width: 1px; border-left-color: rgb(204, 204, 204); border-left-style: solid; padding-left: 1ex; position: static; z-index: auto; "><div style="word-wrap:break-word"><div><div><div><blockquote type="cite">
<div class="gmail_quote"><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0.8ex; border-left-width: 1px; border-left-color: rgb(204, 204, 204); border-left-style: solid; padding-left: 1ex; position: static; z-index: auto; "><div style="word-wrap:break-word"><div>1) right after creation of the matrix:</div>
<div><br></div><div><div>#ifdef petscDev</div><div> ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr);</div><div>#else</div><div> ierr = DAGetMatrix(DMMGGetDA(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr);</div>
<div>#endif</div><div> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);CHKERRQ(ierr);</div></div><div><br></div><div>2) at the beginning of the FormJacobianLocal() routine:</div><div><br></div><div>
<div> PetscFunctionBegin;</div><div> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);CHKERRQ(ierr);</div><div><br></div></div><div>3) before calling MatAssemblyBegin() in FormJacobianLocal() routine:</div>
<div><br></div><div><div> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);CHKERRQ(ierr);</div><div> ierr = MatAssemblyBegin(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div><div> ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div>
<div><br></div></div><div>None of those works. What is wrong here?</div><div><br></div><div>Thanks,</div><div><br></div><div>R</div><div><div style="word-wrap:break-word"><br></div>
</div>
<br><div><div>On Jan 20, 2012, at 4:32 PM, Matthew Knepley wrote:</div><br><blockquote type="cite">On Fri, Jan 20, 2012 at 6:28 PM, Xuefei (Rebecca) Yuan <span dir="ltr"><<a href="mailto:xyuan@lbl.gov" target="_blank">xyuan@lbl.gov</a>></span> wrote:<br>
<div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div style="word-wrap:break-word">Hello Matt,<div><br></div><div>I have changed the code as</div><div><br></div><div><div><span style="text-indent:0px"><span style="text-indent:0px"><div style="word-wrap:break-word"><div style="word-wrap:break-word">
ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);CHKERRQ(ierr);</div><div style="word-wrap:break-word"> ierr = MatAssemblyBegin(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div><div style="word-wrap:break-word">
ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div></div></span></span></div></div></div></blockquote><div><br></div><div>You have to set it before you start setting values, so we know to ignore them.</div>
<div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word"><div><div>
</div><div>but still get the same result as before, the matrix still has 5776 nonzeros:</div><div><br></div><div><div>% Size = 100 100</div><div> 2 % Nonzeros = 5776</div><div> 3 zzz = zeros(5776,3);</div><div><br></div>
</div><div>Then I switch the order as</div><div><br></div><div><div style="word-wrap:break-word"> ierr = MatAssemblyBegin(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div><div style="word-wrap:break-word"> ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div>
</div><div><div style="word-wrap:break-word"> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);CHKERRQ(ierr);</div></div><div><br></div><div>nothing changed.</div><div><br></div><div>The version is 3.1-p8.</div>
<div><br></div><div>Thanks very much!</div><div><br></div><div>Best regards,</div><div><br></div><div>Rebecca</div><div><br></div><div><br></div><div><br></div>
<br><div><div>On Jan 20, 2012, at 4:09 PM, Matthew Knepley wrote:</div><br><blockquote type="cite">On Fri, Jan 20, 2012 at 6:02 PM, Xuefei (Rebecca) Yuan <span dir="ltr"><<a href="mailto:xyuan@lbl.gov" target="_blank">xyuan@lbl.gov</a>></span> wrote:<br>
<div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Hello all,<br>
<br>
This is a test for np=1 case of the nonzero structure of the jacobian matrix. The jacobian matrix is created via<br>
<br>
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);<br>
<br>
ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr);<br>
<br>
After creation of the jacobian matrix,<br>
<br>
ierr = MatAssemblyBegin(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);<br>
ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);<br>
<br>
PetscViewer viewer;<br>
char fileName[120];<br>
sprintf(fileName, "jacobian_after_creation.m");CHKERRQ(ierr);<br>
<br>
FILE * fp;<br>
<br>
ierr = PetscViewerASCIIOpen(PETSC_COMM_WORLD,fileName,&viewer);CHKERRQ(ierr);<br>
ierr = PetscViewerSetFormat(viewer,PETSC_VIEWER_ASCII_MATLAB);CHKERRQ(ierr);<br>
ierr = MatView (jacobian, viewer); CHKERRQ (ierr);<br>
ierr = PetscFOpen(PETSC_COMM_WORLD,fileName,"a",&fp); CHKERRQ(ierr);<br>
ierr = PetscViewerASCIIPrintf(viewer,"spy((spconvert(zzz)));\n");CHKERRQ(ierr);<br>
ierr = PetscFClose(PETSC_COMM_WORLD,fp);CHKERRQ(ierr);<br>
PetscViewerDestroy(&viewer);<br>
<br>
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)<br>
<br>
<br><br>
<br>
Inside the FormJacobianLocal() function, I have selected the index to pass to the nonzero values to jacobian, for example,<br>
<br>
ierr = MatSetValuesStencil(jacobian, 1, &row, 6, col, v, INSERT_VALUES);CHKERRQ(ierr);<br>
<br>
where<br>
<br>
col[0].i = column[4].i;<br>
col[1].i = column[5].i;<br>
col[2].i = column[6].i;<br>
col[3].i = column[9].i;<br>
col[4].i = column[10].i;<br>
col[5].i = column[12].i;<br>
<br>
<br>
col[0].j = column[4].j;<br>
col[1].j = column[5].j;<br>
col[2].j = column[6].j;<br>
col[3].j = column[9].j;<br>
col[4].j = column[10].j;<br>
col[5].j = column[12].j;<br>
<br>
col[0].c = column[4].c;<br>
col[1].c = column[5].c;<br>
col[2].c = column[6].c;<br>
col[3].c = column[9].c;<br>
col[4].c = column[10].c;<br>
col[5].c = column[12].c;<br>
<br>
v[0] = value[4];<br>
v[1] = value[5];<br>
v[2] = value[6];<br>
v[3] = value[9];<br>
v[4] = value[10];<br>
v[5] = value[12];<br>
<br>
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<br>
<br>
<br><br>
<br>
for the true nonzero structures.<br>
<br>
But the ksp_view will give the nonzeros number as 5776, instead of 1075:<br>
<br>
linear system matrix = precond matrix:<br>
Matrix Object: Mat_0x84000000_1 1 MPI processes<br>
type: seqaij<br>
rows=100, cols=100<br>
total: nonzeros=5776, allocated nonzeros=5776<br>
<br>
It is a waste of memory to have all those values of zeros been stored in the jacobian.<br>
<br>
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.<br></blockquote><div><br></div><div>MatSetOption(<span style="font-size:medium;font-family:Times">MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE);</span></div>
<div><span style="font-size:medium;font-family:Times"><br></span></div><div><span style="font-size:medium;font-family:Times"> Matt</span></div><div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Thanks very much!<br>
<br>
Have a nice weekend!<br>
<br>
Cheers,<br>
<br>
Rebecca<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br></blockquote></div><br><br clear="all"><span><font color="#888888"><span><font color="#888888"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>
</font></span></font></span></blockquote></div><span><font color="#888888"><br></font></span></div></div></blockquote></div><span><font color="#888888"><br><br clear="all"><div><br></div>-- <br>
What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>
</font></span></blockquote></div><br></div></blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>
</blockquote></div><br></div></div></div><br></blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>
</blockquote></div><br></div></body></html>