On Fri, Jun 25, 2010 at 7:25 AM, Li, Zhisong (lizs) <span dir="ltr">&lt;<a href="mailto:lizs@mail.uc.edu">lizs@mail.uc.edu</a>&gt;</span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">





<div>
<div style="direction:ltr;font-family:Tahoma;color:rgb(0, 0, 0);font-size:13px">
<div><br>
Hi, Petsc Team, <br>
<br>
Recently I encounter a weird problem for segmentation violation. I wrote a simple test code to describe it. Here the line  &quot; pp = sk[j+1][i].p; &quot; causes segmentation violation trouble when I try to invoke values of ghost points in j direction. If I change it
 into &quot;pp = sk[j][i+1].p;&quot; invoking ghost point values in i diection, then it works smoothly. I check previous archives about segmentation violation, but cannot find any clue for this. Can you point out where is wrong here or is it a bug?<br>
</div></div></div></blockquote><div><br></div><div>When you call DAVecGetArray() on a GLOBAL vector, there are no ghost values. You need</div><div>a LOCAL vector. The i+1 is still wrong, but does not SEGV since you look into memory you own.</div>
<div><br></div><div>   Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;"><div><div style="direction:ltr;font-family:Tahoma;color:rgb(0, 0, 0);font-size:13px">
<div>
Thank you.<br>
<br>
<br>
Zhisong Li<br>
<br>
<br>
<br>
static char help[] = &quot;test&quot;;<br>
#include &quot;petscda.h&quot;<br>
<br>
typedef struct { PetscScalar p; } Field;<br>
<br>
int main(int argc, char **args)<br>
{ Vec  xx;  <br>
  PetscInt   dof = 1, m = 24, n= 32, i, j, xs, ys, xm, ym, yints, yinte, xints, xinte;<br>
  PetscScalar pp;<br>
  DA   da;<br>
  Field  **sk;<br>
  PetscInitialize(&amp;argc, &amp;args, (char *)0, help) ;<br>
    <br>
  DACreate2d(PETSC_COMM_WORLD, DA_NONPERIODIC, DA_STENCIL_STAR, m, n, PETSC_DECIDE, PETSC_DECIDE, dof, 1, PETSC_NULL, PETSC_NULL, &amp;da);<br>
  DACreateGlobalVector(da, &amp;xx); <br>
  DAGetCorners(da, &amp;xs, &amp;ys, 0, &amp;xm, &amp;ym, 0);<br>
      xints = xs;   xinte = xs+xm;   yints = ys;   yinte = ys+ym;    <br>
   <br>
  VecSet(xx,1.0);<br>
  <br>
  DAVecGetArray(da, xx, &amp;sk);    <br>
  if (xints == 0){ xints = xints + 1; }  <br>
  if (yints == 0){ yints = yints + 1; }                                                                                                      <br>
  if (xinte == m){ xinte = xinte - 1; }             <br>
  if (yinte == n){ yinte = yinte - 1; }  <br>
<br>
    for (j=yints; j&lt;yinte; j++){ <br>
      for (i=xints; i&lt;xinte; i++)   { pp = sk[j+1][i].p; }<br>
                               }      <br>
                          <br>
  DAVecRestoreArray(da, xx, &amp;sk);                                                     <br>
  VecDestroy(xx);<br>
  DADestroy(da);         <br>
  PetscFinalize();<br>
  PetscFunctionReturn(0);<br>
}<br>
<br>
</div>
</div>
</div>

</blockquote></div><br><br clear="all"><br>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener<br>