<div dir="ltr">Hi Barry,<div><br></div><div>Thanks for your message. When I use DMCreateInterpolation, I got the error message like "[0]PETSC ERROR: Invalid argument</div><div>[0]PETSC ERROR: Coarsening factor in x must be 2".  The error will be gone if I change refine_x = 3 to refine_x=2. Why must the coarsening factor be a factor of 2?</div><div><br></div><div>Thank you.</div><div><br></div><div>Xiangdong</div><div><br></div><div>Below is my test code:</div><div><br></div><div><div>int main(int argc, char **argv)</div><div>{</div><div>  PetscErrorCode ierr;</div><div><br></div><div>  PetscInitialize(&argc, &argv, PETSC_NULL, PETSC_NULL);</div><div><br></div><div>  int Nx=6, Ny=4, Nz=1;</div><div>    </div><div>  DM da;</div><div>  </div><div>  ierr = DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE, DMDA_STENCIL_STAR,Nx,Ny,PETSC_DECIDE, PETSC_DECIDE,1,1,NULL,NULL,&da);CHKERRQ(ierr);</div><div><br></div><div>  int refine_x=3, refine_y=2, refine_z=1;</div><div>  ierr = DMDASetRefinementFactor(da,refine_x,refine_y,refine_z);</div><div>    </div><div>  ierr = DMDASetInterpolationType(da,DMDA_Q0);</div><div><br></div><div>  ierr = DMSetFromOptions(da); CHKERRQ(ierr);</div><div>  ierr = DMSetUp(da); CHKERRQ(ierr);</div><div><br></div><div>  DM dac;</div><div>  ierr = DMCoarsen(da,NULL,&dac);</div><div> </div><div>  Mat M;</div><div>  Vec V;</div><div>  </div><div>  ierr = DMCreateInterpolation(dac,da,&M,&V);</div><div>  </div><div>  MatView(M,PETSC_VIEWER_STDOUT_WORLD);</div><div><br></div><div><br></div><div>  ierr = PetscFinalize(); CHKERRQ(ierr);</div><div>  PetscFunctionReturn(0);</div><div>}</div></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Mon, May 7, 2018 at 6:16 PM, Smith, Barry F. <span dir="ltr"><<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><span class=""><br>
<br>
<br>
> On May 7, 2018, at 3:25 PM, Xiangdong <<a href="mailto:epscodes@gmail.com">epscodes@gmail.com</a>> wrote:<br>
> <br>
> Hello everyone,<br>
> <br>
> I have a 2D DMDA with Nx=6 and Ny=4. If I set refine_x=3 and refine_y=2, I can get a coarse DMDA by DMCoarsen(). However, when I tried to call the DMCreateRestriction between these two DMDAs, I got error about " No support for this operation for this object type. DMCreateRestriction not implemented for this type."<br>
> <br>
> Can you help me get around this?<br>
<br>
</span>   Us DMCreateInterpolation() and then use MatMultTranspose() for example to apply the restriction.<br>
<span class=""><br>
> <br>
> Another simple question, where can I find the explanation/meaning of DMDA_Q0 and DMDA_Q1?<br>
<br>
</span>    This is finite element terminology, Q0 means piecewise constant basis functions (and hence piecewise constant interpolation) while Q1 means linear interpolation.<br>
For Q0 you should think of the unknowns as living on the centers of each element while for Q1 the unknowns live on the vertices of the elements.<br>
<span class="HOEnZb"><font color="#888888"><br>
   Barry<br>
</font></span><div class="HOEnZb"><div class="h5"><br>
> <br>
> Below is my test example.<br>
> <br>
> Thank you.<br>
> <br>
> Best,<br>
> Xiangdong<br>
> <br>
> #include <petsc.h><br>
> <br>
> int main(int argc, char **argv)<br>
> {<br>
>   PetscErrorCode ierr;<br>
> <br>
>   PetscInitialize(&argc, &argv, PETSC_NULL, PETSC_NULL);<br>
> <br>
>   int Nx=6, Ny=4, Nz=1;<br>
>     <br>
>   DM da;<br>
>   <br>
>   ierr = DMDACreate2d(PETSC_COMM_WORLD,<wbr>DM_BOUNDARY_NONE,DM_BOUNDARY_<wbr>NONE, DMDA_STENCIL_STAR,Nx,Ny,PETSC_<wbr>DECIDE, PETSC_DECIDE,1,1,NULL,NULL,&<wbr>da);CHKERRQ(ierr);<br>
> <br>
>   int refine_x=3, refine_y=2, refine_z=1;<br>
>   ierr = DMDASetRefinementFactor(da,<wbr>refine_x,refine_y,refine_z);<br>
>     <br>
>   ierr = DMDASetInterpolationType(da,<wbr>DMDA_Q0);<br>
> <br>
>   ierr = DMSetFromOptions(da); CHKERRQ(ierr);<br>
>   ierr = DMSetUp(da); CHKERRQ(ierr);<br>
> <br>
>   DM dac;<br>
>   ierr = DMCoarsen(da,NULL,&dac);<br>
>  <br>
>   Mat M;<br>
>   Vec V;<br>
>   <br>
>   ierr = DMCreateRestriction(dac,da,&M)<wbr>;<br>
>   <br>
>   MatView(M,PETSC_VIEWER_STDOUT_<wbr>WORLD);<br>
> <br>
> <br>
>   ierr = PetscFinalize(); CHKERRQ(ierr);<br>
>   PetscFunctionReturn(0);<br>
> }<br>
> <br>
<br>
</div></div></blockquote></div><br></div>