<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><div class=""><br class=""></div> This is becoming a nightmare. <div class=""><br class=""></div><div class=""> You pointed out something I had not thought about. Not only is the nonzero structure different in your case but so is the local to global mapping. Because certain of your ghost points are connected to different global points than for the standard case, thus the calls to MatSetValuesStencil() will not work correctly for the nonstandard ghost points.</div><div class=""><br class=""></div><div class="">The mapping is constructed in DMSetUp_DA_3D() and it is very tedious code. In theory, the idx[] entries could be patched just before the call to ISLocalToGlobalMappingCreate() but it may be tedious to figure out exactly which local points need to have their global point changed. In theory it is easy but the code is extremely tedious to get to exactly the place that needs change.</div><div class=""><br class=""></div><div class=""> I'm sorry, we may have been attempting to push the DMDA construct too far from its simple, save home territory and thus ending up in a morass. </div><div class=""><br class=""></div><div class=""> I suggest looking at the DMSetUp_DA_3D() to see if you can figure out how to get the specific ghost points properly mapped but I am not sure if that is a reasonable task. Otherwise, I have no good suggestions. If the local to global mapping cannot be fixed then the entire model of using MatSetValuesStencil falls apart.</div><div class=""><br class=""></div><div class="">Barry</div><div class=""><br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Sep 13, 2022, at 12:03 AM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" class="">Jiannan_Tu@uml.edu</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><meta charset="UTF-8" class=""><div class="WordSection1" style="page: WordSection1; caret-color: rgb(0, 0, 0); font-family: Helvetica; font-size: 18px; font-style: normal; font-variant-caps: normal; font-weight: 400; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none;"><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Indeed, it is not a standard grid structure that DMDA handles by default. Normally, a grid (k, j, i) has connectivity at grid (k, j-1, i). I understand if j > 0, DMDA correctly allocates memory for it. When j = 0, j-1 = -1, this (k, -1, i) grid is actually the grid across the pole with k-index = (k+a3/2) % a3, where a3 is the number of grids for k-index. Attached figure illustrates a pair of such grids. Assume red grid is at (k = 10, j = 0, i), its corresponding grid *blue dot) cross the pole is at (k = 55, j =0, i) if a3 = 90.<o:p class=""></o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""><o:p class=""> </o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">The question is how to make DMDA allocate memory at such local column stencils like (col.k = 55, col.j = 0, col.i = i, col.c = 0), instead of the column stencil (col.k = 10, col.j = -1, col.i = i, col.c = 0).<o:p class=""></o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""><o:p class=""> </o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I tried insert at stencils like (col.k = 10, col.j = -1, col.i = i, col.c = 0). The matrix assemble has no issues but TSSolve failed. I guess probably due to wrong column locations of those Jacobian elements.<o:p class=""></o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""><o:p class=""> </o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you,<o:p class=""></o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan<o:p class=""></o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""><o:p class=""> </o:p></span></div><div class=""><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="Apple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><span class="Apple-converted-space"> </span><br class=""><b class="">Sent:</b><span class="Apple-converted-space"> </span>Monday, September 12, 2022 4:18 PM<br class=""><b class="">To:</b><span class="Apple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Cc:</b><span class="Apple-converted-space"> </span>petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a>><br class=""><b class="">Subject:</b><span class="Apple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><o:p class=""> </o:p></div><div style="border: 1pt hidden rgb(255, 202, 202); padding: 2pt;" class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span style="font-size: 12pt;" class=""><span class="Apple-converted-space"> </span>This email was sent from outside the UMass Lowell network.<o:p class=""></o:p></span></div></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><o:p class=""> </o:p></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><o:p class=""> </o:p></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> I understand in principle, I just cannot wrap my head around a sphere :-)<o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><o:p class=""> </o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><br class=""><br class=""><o:p class=""></o:p></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 12, 2022, at 3:45 PM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><o:p class=""> </o:p></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Simulation domain is in spherical coordinates. The grid along co-latitude direction is that j = 0 on one side of the pole and j = -1 (in fact index 0) on the other side of the pole. The longitude index k for j = 0, say is 0 degree then for j = -1 is 180 degrees. kcm accounts for this change in k-index when going from one side of the pole to another side.<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I'll see if I can solve the problem with the method you suggested.<span class="apple-converted-space"> </span><o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you very much.<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan<o:p class=""></o:p></span></div></div><div class="MsoNormal" align="center" style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; text-align: center;"><hr size="2" width="977" align="center" style="width: 732.8pt;" class=""></div><div id="divRplyFwdMsg" class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="apple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><br class=""><b class="">Sent:</b><span class="apple-converted-space"> </span>Monday, September 12, 2022 9:31 AM<br class=""><b class="">To:</b><span class="apple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Cc:</b><span class="apple-converted-space"> </span>petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a>><br class=""><b class="">Subject:</b><span class="apple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""></o:p></span></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> <o:p class=""></o:p></span></div></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="apple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><span style="font-size: 12pt;" class=""><o:p class=""></o:p></span></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> I am not able to visualize the pattern you refer to but I am pretty sure it is not a type of connectivity that DMDA handles by default since DMDA is for relatively simply structured meshes. Here is what I suggest.<span class="Apple-converted-space"> </span><o:p class=""></o:p></span></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> The routine that does the preallocation for DMDA that you are using is in src/dm/impls/da/fdda.c called DMCreateMatrix_DA_3d_MPIAIJ_Fill(). On each MPI rank it determines the connectivity for each local grid vertex (calledthe variable row in the routine) with all the other grid vertices (these are stored temporarily in the array cols[]).<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">The call to MatPreallocateSetLocal() takes this coupling information and stores it temporarily. The calls to <o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> PetscCall(MatSeqAIJSetPreallocation(J, 0, dnz));<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> PetscCall(MatMPIAIJSetPreallocation(J, 0, dnz, 0, onz));<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> MatPreallocateEnd(dnz, onz);<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> PetscCall(MatSetLocalToGlobalMapping(J, ltog, ltog));<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">at the end of the loop over the local grid vertices then provide the information to the matrix so it now has the correct preallocation and correct local to global mapping needed by MatSetValuesStencil()<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">The same loop structure is then called again where it inserts zero values into the matrix at all the correctly preallocated locations.<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">You can copy this routine and modify it slightly to include the "special coupling" that your problem requires. I do not understand your formula for kcm so cannot tell you exactly what you need to change but it may correspond to remapping the kk values at the extremes of the loop for (kk = kstart; kk < kend + 1; kk++) { so that the<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">values nc * (slot + ii + gnx * jj + gnx * gny * kk) correspond to the correct coupling location. I recommend drawing a picture of a case with a very small size for nx,ny, and nz so you can see exactly how the mapping takes place and trace through your modified code that it does what you need.<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">Barry<o:p class=""></o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><br class=""><br class=""><o:p class=""></o:p></span></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">On Sep 11, 2022, at 11:01 PM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></span></div></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><o:p class=""> </o:p></span></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Barry,</span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Creating DNDA with periodic boundary BC and using km=k-1; kp=k+1 even at the BCs solve the “augment out of range” problem for inserting elements at theae boundaries of k-index. Now the problem is due to wrapping around the north pole when j=0 and jm = -1 is reset as j=0 but on the other side of the pole with k-index replacing by kcm = [kmax+1)/2\ % (kmax+1). In my case, the “Augment out of range” happens at global row/column (26, 5265026), exactly at wrapping jm=-1 to j=0 to the other side.</span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I am thinking about how to properly use ghost grid j = -1 by setting appropriate BC there and inserting at that location without wrapping.</span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan</span><o:p class=""></o:p></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div><div class=""><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="xapple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><span class="xapple-converted-space"> </span><br class=""><b class="">Sent:</b><span class="xapple-converted-space"> </span>Sunday, September 11, 2022 12:10 PM<br class=""><b class="">To:</b><span class="xapple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Cc:</b><span class="xapple-converted-space"> </span>petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a>><br class=""><b class="">Subject:</b><span class="xapple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><o:p class=""> </o:p></p></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 11, 2022, at 9:56 AM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Barry,</span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you.</span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Yes, the DMDA is created with periodic BC. But it might not be the periodic BC since in azimuthal direction, the ghost grid -1 is actually grid Nmax, and Nmax+1 is grid 0.</span><o:p class=""></o:p></div></div></div></div></blockquote><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Consider the following, one has n points from 0 to n-1 where n = 4, so I think with your definition of Nmax, Nmax is 3, Nmax + 1 is 4,<o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> -1 0 1 2 3 4<o:p class=""></o:p></div></div></div><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> + * * * * +<o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Thus there is one grid point at each end. The * represent regular grid points and the + ghost locations. Is this your situation? <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">This is what we call periodic with a stencil width of 1 with DMDA and the DMCreateMatrix() should automatically provide the correct preallocation and nonzero structure.<o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Can you please do a MatView() on the result from DMCreateMatrix() for a very small grid before you set any values and you should see the locations for the coupling across the periodic boundary.<o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">It sounds like you have a more complicated cross periodic coupling, we need to understand it exactly so we can determine how to proceed.<o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Barry<o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"> <o:p class=""></o:p></p></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Still the problem is how to properly pre-allocate matrix. Otherwise either “Augument out of range” error or code stuck at inserting values if MatSetOption() is called.</span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan</span><o:p class=""></o:p></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div><div class=""><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="xapple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><span class="xapple-converted-space"> </span><br class=""><b class="">Sent:</b><span class="xapple-converted-space"> </span>Sunday, September 11, 2022 12:10 AM<br class=""><b class="">To:</b><span class="xapple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Cc:</b><span class="xapple-converted-space"> </span>petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a>><br class=""><b class="">Subject:</b><span class="xapple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Have you created the DMDA with the periodic boundary conditions arguments? Or are your boundary conditions not simply periodic? <o:p class=""></o:p></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Compare DMCreateMatrix_DA_3d_MPIAIJ() and DMCreateMatrix_DA_3d_MPIAIJ_Fill() in src/dm/impls/da/fdda.c I think they both attempt to handle the periodic case.<o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Barry<o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><o:p class=""></o:p></p></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 10, 2022, at 10:14 PM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Barry,</span><o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I think I understand how to set values for the intra-grid coupling. But I am not clear about the way to set values for inter-grid coupling. When a component couples to itself at<span class="xapple-converted-space"> </span></span>above, below, right, left, front and back grids, diagonal is set as 1 for the second array of DMDASetBlockFills(). How about if involving more grids? How does DMDA know where to allocate memory? In my case, a component couples to itself at the other end of the grid because of periodic boundary condition. There is an error message “Augment out of range . Insert new nonzero at global row/column (26, 520526) into matrix” because of this far away coupling. The column index 520526 is related to the other end of the grid in azimuthal direction.<o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Thanks,<o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Jiannan<o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div><div class=""><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="xapple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><span class="xapple-converted-space"> </span><br class=""><b class="">Sent:</b><span class="xapple-converted-space"> </span>Saturday, September 10, 2022 1:10 PM<br class=""><b class="">To:</b><span class="xapple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Cc:</b><span class="xapple-converted-space"> </span>petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a>><br class=""><b class="">Subject:</b><span class="xapple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> If you are using DMCreateMatrix() then <span style="font-size: 12pt;" class="">MatMPIAIJSetPreallocation is not needed, it is automatically handled by the DM.</span><o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> It sounds like that the values you set in DMDASetBlockFills() did not include all the coupling, hence the preallocation was wrong, hence the time to insert values was too long.</span><o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> Barry</span><o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 10, 2022, at 8:16 AM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Just found MatMPIAIJSetPreallocation is needed. Now Jacobin insertion and assembly is fast.</span><o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div><div class=""><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="xapple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><span class="xapple-converted-space"> </span><br class=""><b class="">Sent:</b><span class="xapple-converted-space"> </span>Thursday, September 8, 2022 9:53 PM<br class=""><b class="">To:</b><span class="xapple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Subject:</b><span class="xapple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 8, 2022, at 8:34 PM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Barry,</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you.</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I set up these two arrays according to the non-zero pattern of the Jacobian. The DMDASetBlockFills() is then called after the call to DMCreateMatrix(da, &A).</span><o:p class=""></o:p></div></div></div></div></div></div></div></blockquote><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Absolutely not! It must be called before DMCreateMatrix() so that the DM knows how to construct exactly the matrix you need.<o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Nevertheless, the program execution is still killed with exit code 9, which I assume is due to overuse of the memory. My desktop machine has 32 GB RAM.</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I ran the code with a smaller mesh, 100x45x90, the result is the same. With that number of grids, 26 unknowns, and up to 10 non-zero at each matrix row, the estimated memory to store 100x45x90x26x10 elements should be less than 1 GB (double precision). I am wondering why the matrix still takes too much memory. Maybe I have to use matrix-free?</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you,</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan</span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div><div class=""><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="xapple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><span class="xapple-converted-space"> </span><br class=""><b class="">Sent:</b><span class="xapple-converted-space"> </span>Thursday, September 8, 2022 4:19 PM<br class=""><b class="">To:</b><span class="xapple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Subject:</b><span class="xapple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 8, 2022, at 3:39 PM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class="">Barry,</span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class="">Thank you very much for the detailed description.</span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class="">So the first array, dfill, is for intra grids and the second one, ofill is for the inter-grid coupling?</span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class="">In the case for inter-grid coupling, a component is only coupled to itself, say x in your example is coupled to itself on above, below, right, left, front, and back, how to set values for the second array?</span></span><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div></blockquote><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;">Then the off-diagonal one would just have values on its diagonal.<br class=""><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Jiannan<o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:<span class="xapple-converted-space"> </span></b><a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">Barry Smith</a><br class=""><b class="">Sent:<span class="xapple-converted-space"> </span></b>Thursday, September 8, 2022 2:12 PM<br class=""><b class="">To:<span class="xapple-converted-space"> </span></b><a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Tu, Jiannan</a><br class=""><b class="">Cc:<span class="xapple-converted-space"> </span></b><a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a><br class=""><b class="">Subject:<span class="xapple-converted-space"> </span></b>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 8, 2022, at 1:47 PM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class="">Barry,</span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class="">Thank you very much.</span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">DMDASetBlockFills() needs two input arrays specifying sparsity pattern. In my case, I have 26 components at each grid. Is it correct that I need to have two arrays of 26x26 elements with value 1 for coupling and 0 for no coupling?<o:p class=""></o:p></div></div></div></div></div></div></div></div></div></blockquote><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Exactly<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Also how to determine when the value should be 1 or 0?<span class="xapple-converted-space"> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></blockquote><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> That is problem dependent. Given paper-and-pencil you know for your stencil which components are coupled with other components at each grid point and which components are coupled with components at neighboring grid points. <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Note the second array contains information about all the neighboring grid points. Above, below, right, left, front back. Normally the coupling is symmetric so the nonzero structures of those 6 blocks are identical. If, for you particular problem, the nonzero structures of the six are not identical you use the union of the nonzero structures of the six of them.<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> 2d example say you have x,y,z at each grid point and the inter-point coupling is x is connected to y and y is connect to x then the first array is<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">[ 1 1 0 ; 1 1; 0 ; 0 0 1]; the zeros are because x is not coupled to z and y is not coupled to z.<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">For intra grid points, say x is coupled to z and itself, z is coupled to x and itself and y is not coupled to itself then<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">[1 0 1; 0 0 0; 1 0 1]; the middle row is zero because y is not coupled to anything and the zeros in the first and last are because x and z are not coupled to y.<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"> <br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Each component involves 7 stencils (i-1, j, k), (I, j, k), (i+1, j, k), (I, j-1, k), (I, j+1, k), (I, j, k-1), and (I, j, k+1), and couples with several other components.<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Jiannan<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">Sent from<span class="xapple-converted-space"> </span><a href="https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgo.microsoft.com%2Ffwlink%2F%3FLinkId%3D550986&data=05%7C01%7CJiannan_Tu%40uml.edu%7Cc61340312eb04774dc8708da94fbdca2%7C4c25b8a617f746f983f054734ab81fb1%7C0%7C0%7C637986106759336186%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=7yOcSuUJIzB88palFhM7LvgDsXFS0zDeJQUXCH8QztA%3D&reserved=0" style="color: blue; text-decoration: underline;" class="">Mail</a><span class="xapple-converted-space"> </span>for Windows<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div style="border-style: solid none none; border-top-width: 1pt; border-top-color: rgb(225, 225, 225); padding: 3pt 0in 0in;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:<span class="xapple-converted-space"> </span></b><a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">Barry Smith</a><br class=""><b class="">Sent:<span class="xapple-converted-space"> </span></b>Wednesday, September 7, 2022 11:53 AM<br class=""><b class="">To:<span class="xapple-converted-space"> </span></b><a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Tu, Jiannan</a><br class=""><b class="">Cc:<span class="xapple-converted-space"> </span></b><a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a><br class=""><b class="">Subject:<span class="xapple-converted-space"> </span></b>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span class="xdefaultfonthxmailstyle0"><span style="font-size: 12pt;" class=""> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">DMDASetBlockFills() or DMDASetBlockFillsSparse() (they are just different ways of providing the same information) will help you here enormously, your type of case is exactly what they were designed for.<o:p class=""></o:p></div></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> Barry<o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></blockquote></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><o:p class=""></o:p></p></div></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class="">On Sep 7, 2022, at 10:47 AM, Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>> wrote:<o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""> <o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Barry and Hong,</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you.<span class="xapple-converted-space"> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">There are 26 components at each grid and there are not fully coupled in terms of stiff functions. Mutual coupling is among about 6 components.<span class="xapple-converted-space"> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I would prefer not using matrix-free since the Jacobina is not difficult to calculate and only up to 10 non-zeros at each row. I'll try DMDASetBlockFills() or DMDASetBlockFillsSparse() and see how they can reduce the memory usage.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><E6340AADF8494DB7861EA4659E7713CE.png><o:p class=""></o:p></div></div></div></div></div></div></div></div><div id="x_divRplyFwdMsg" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><b class="">From:</b><span class="xapple-converted-space"> </span>Barry Smith <<a href="mailto:bsmith@petsc.dev" style="color: blue; text-decoration: underline;" class="">bsmith@petsc.dev</a>><br class=""><b class="">Sent:</b><span class="xapple-converted-space"> </span>Tuesday, September 6, 2022 11:33 PM<br class=""><b class="">To:</b><span class="xapple-converted-space"> </span>Tu, Jiannan <<a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class="">Jiannan_Tu@uml.edu</a>><br class=""><b class="">Cc:</b><span class="xapple-converted-space"> </span><a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a><span class="xapple-converted-space"> </span><<a href="mailto:petsc-users@mcs.anl.gov" style="color: blue; text-decoration: underline;" class="">petsc-users@mcs.anl.gov</a>><br class=""><b class="">Subject:</b><span class="xapple-converted-space"> </span>Re: [petsc-users] Using matrix-free with KSP<o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif; line-height: 16pt; background-color: rgb(255, 202, 202);" class=""><b class=""><span style="font-size: 12pt;" class="">CAUTION:</span></b><span class="xapple-converted-space"><span style="font-size: 12pt;" class=""> </span></span><span style="font-size: 12pt;" class="">This email was sent from outside the UMass Lowell network.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""></span><o:p class=""></o:p></p></div></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">On Sep 6, 2022, at 11:00 PM, Tu, Jiannan <</span><a href="mailto:Jiannan_Tu@uml.edu" style="color: blue; text-decoration: underline;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">Jiannan_Tu@uml.edu</span></a><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class="">> wrote:</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I am using TS IMEX to solve a large DAE system. The DAE is obtained by applying finite FV method to 3-D multi-specie ion/neutral fluids equations with magnetic induction equation. The Jacobian for stiff part is formed by using MatSetValuesStencil(). The Jacobian matrix is very sparse, no more than 10 non-zeros on each row. MatSetValuesStencil requires local to global mapping by calling ISLocaltoGlobalMapping(). Could you please instruct me how to use local to global mapping?<span class="xxapple-converted-space"> </span></span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></div></div></div></blockquote><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> DMDA automatically sets up the </span><span style="font-size: 12pt;" class="">ISLocaltoGlobalMapping() so you should not need to.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">I also tried using DMCreateMatrix() to create the Jacobian. While local to global mapping is not necessary, the matrix takes too much memory and requires 64-bit indices. I would prefer to take the advantage of sparsity of the Jacobian, pre-allocate the matrix to use as less as possible memory so that the code can be run on a multi-core desktop.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></div></div></div></blockquote><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> If using matrix-free does not work for you because the linear solves do not converge or converge too slowly. Then you might be able to decrease the memory used in the matrix. The DMCreateMatrix() does take advantage of sparsity and tries to preallocate only what it needs. Often what it preallocates is the best possible, but for multicomponent problems it has to assume there is full coupling within all the degrees of freedom that live at each at each grid point. How many components live at each grid point in your model and are they not all coupled to each other in the equations? If they are not fully coupled you can use either DMDASetBlockFills() or DMDASetBlockFillsSparse() to indicate the reduced coupling that actually exists for you model.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""> Barry</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><p class="MsoNormal" style="margin: 0in 0in 12pt; font-size: 11pt; font-family: Calibri, sans-serif;"><span style="font-size: 13.5pt; font-family: Helvetica, sans-serif;" class=""><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""><br class=""></span><o:p class=""></o:p></p></div></div></div></div></div></div></div><blockquote style="margin-top: 5pt; margin-bottom: 5pt;" class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Thank you very much for your advice.</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class=""> </span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div class=""><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><span style="font-size: 12pt;" class="">Jiannan</span><o:p class=""></o:p></div></div></div></div></div></div></div></div></div></div></div></div></blockquote></div></div></div></div></blockquote></div></div></div></blockquote></div></div></div></blockquote></div></div></div></blockquote></div></div></div></blockquote></div></div></div></div></blockquote></div></div></div></div></blockquote></div></div></div></div></blockquote></div></div></div></div></div></blockquote></div><div style="margin: 0in; font-size: 11pt; font-family: Calibri, sans-serif;" class=""><o:p class=""> </o:p></div></div></div><span id="cid:F5CE0FDF-74E1-452E-9793-641CFA2BA366"><fig2.png></span></div></blockquote></div><br class=""></div></body></html>