On Wed, Nov 30, 2011 at 9:59 AM, Treue, Frederik <span dir="ltr"><<a href="mailto:frtr@risoe.dtu.dk">frtr@risoe.dtu.dk</a>></span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div lang="DA" link="blue" vlink="purple"><div><p class="MsoNormal">Hi everyone,<u></u><u></u></p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"><span lang="EN-US">Caveat: I have just started using petsc, so the answer to my question may very well be fairly trivial.</span></p>
</div></div></blockquote><div><br></div><div>See <a href="http://www.mcs.anl.gov/petsc/petsc-dev/src/snes/examples/tutorials/ex5.c.html">SNES ex5</a> for the right way to interact with the DMDA. We will preallocate the matrix for you and allow</div>
<div>you to set values using a stencil.</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 lang="DA" link="blue" vlink="purple">
<div><p class="MsoNormal"> </p><p class="MsoNormal"><span lang="EN-US">I’m trying to run the following bits of code:<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"><u></u> <u></u></span></p><p class="MsoNormal">
<span lang="EN-US">DMDACreate2d(PETSC_COMM_WORLD, DMDA_BOUNDARY_GHOSTED, DMDA_BOUNDARY_GHOSTED, DMDA_STENCIL_BOX,10,10,PETSC_DECIDE,PETSC_DECIDE,1,1,PETSC_NULL,PETSC_NULL,&da);<u></u><u></u></span></p><p class="MsoNormal">
<span lang="EN-US">[snip]<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US">MatCreate(PETSC_COMM_WORLD,&((*FD).ddx));<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> MatSetSizes((*FD).ddx,PETSC_DECIDE,PETSC_DECIDE,100,100);<u></u><u></u></span></p>
<p class="MsoNormal"><span lang="EN-US"> </span>MatSetFromOptions((*FD).ddx);<u></u><u></u></p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"> for (i=0;i<10;i++) {<u></u><u></u></p><p class="MsoNormal">
<span lang="EN-US">col[0]=i*10;col[1]=i*10+1; row[0]=i*10;<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> val[0]=1;val[1]=1;<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> MatSetValues((*FD).ddx,1,row,2,col,val,INSERT_VALUES);<u></u><u></u></span></p>
<p class="MsoNormal"> for (j=1;j<10-1;j++) {<u></u><u></u></p><p class="MsoNormal"> col[0]=i*10+j-1;col[1]=i*10+j+1; row[0]=i*10+j;<u></u><u></u></p><p class="MsoNormal"> <span lang="EN-US">val[0]=-1;val[1]=1;<u></u><u></u></span></p>
<p class="MsoNormal"><span lang="EN-US"> MatSetValues((*FD).ddx,1,row,2,col,val,INSERT_VALUES);<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> }<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> col[0]=i*10+10-2;col[1]=i*10+10-1; row[0]=i*10+10-1;<u></u><u></u></span></p>
<p class="MsoNormal"><span lang="EN-US"> val[0]=-1;val[1]=-1;<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> MatSetValues((*FD).ddx,1,row,2,col,val,INSERT_VALUES);<u></u><u></u></span></p><p class="MsoNormal">
<span lang="EN-US"> }<u></u><u></u></span></p><p class="MsoNormal"><span lang="EN-US"> MatAssemblyBegin((*FD).ddx,MAT_FINAL_ASSEMBLY);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US">MatAssemblyEnd((*FD).ddx,MAT_FINAL_ASSEMBLY);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"><u></u> <u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt">MatScale((*FD).ddx,1/(2*(1/9)));<u></u><u></u></p><p class="MsoNormal" style="text-indent:4.5pt">
[snip]<u></u><u></u></p><p class="MsoNormal" style="text-indent:4.5pt"> DMCreateGlobalVector(da,&tmpvec2);<u></u><u></u></p><p class="MsoNormal" style="text-indent:4.5pt"> <span lang="EN-US">VecSet(tmpvec2,1.0);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> VecAssemblyBegin(tmpvec2);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> VecAssemblyEnd(tmpvec2);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> </span>DMCreateGlobalVector(da,&tmpvec3);<u></u><u></u></p><p class="MsoNormal" style="text-indent:4.5pt"> VecSet(tmpvec3,1.0);<u></u><u></u></p><p class="MsoNormal" style="text-indent:4.5pt">
<span lang="EN-US">VecAssemblyBegin(tmpvec3);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> VecAssemblyEnd(tmpvec3);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt">
<span lang="EN-US"> MatView((*FD).ddx,PETSC_VIEWER_STDOUT_WORLD);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> VecView(tmpvec2,PETSC_VIEWER_STDOUT_WORLD);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> MatMult((*FD).ddx,tmpvec2,tmpvec3);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> VecView(tmpvec3,PETSC_VIEWER_STDOUT_WORLD);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> int tid,first,last;<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> MPI_Comm_rank(PETSC_COMM_WORLD, &tid);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> sleep(1);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> MatGetOwnershipRange((*FD).ddx,&first,&last);<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"> printf("rank: %d,first: %d,last: %d\n",tid,first,last);<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US"><u></u> <u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US">When running it on a single processor, everything works as expected, see attached file seqRes<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt">
<span lang="EN-US">However when running with 4 processors (mpirun –np 4 ./progname) I get the output in mpiRes. Notice that there really is a difference, its not just a surprising division of points between the processes – I checked this with PETSC_VIEWER_DRAW_WORLD. How come? I notice that although in the end each process postulates that it has 25 rows, the result of matview is<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US">Matrix Object: 1 MPI processes<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US">type: mpiaij<u></u><u></u></span></p>
<p class="MsoNormal" style="text-indent:4.5pt"><span lang="EN-US">Is this OK? And if not, what am I doing wrong, presumably in the matrix allocation code?<u></u><u></u></span></p><p class="MsoNormal" style="text-indent:4.5pt">
<span lang="EN-US"><u></u> <u></u></span></p><p class="MsoNormal"><span lang="EN-US"><u></u> <u></u></span></p><p class="MsoNormal">---<u></u><u></u></p><p class="MsoNormal">yours sincerily<u></u><u></u></p><p class="MsoNormal">
Frederik Treue<u></u><u></u></p><p class="MsoNormal"><u></u> <u></u></p></div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>