<p>Ex48 uses periodic in y and z and it works.</p>
<div class="gmail_quote">On Oct 9, 2011 8:08 PM, "Barry Smith" <<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>> wrote:<br type="attribution"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
Ok 3.1-p8 also fails with the same problems. Likely periodic boundary conditions have never worked in 3d? Yikes this looks bad<br>
<br>
Barry<br>
<br>
On Oct 9, 2011, at 3:28 PM, Barry Smith wrote:<br>
<br>
><br>
> I have created a new test example in petsc-3.2 for 3D global to local with periodic boundary conditions and it fails rather completely so whomever introduced all the IXs, IXe stuff into DMSetUp_DA_3d() likely introduced some bugs, or the bugs have always been there and no one has tested it properly.<br>
><br>
> For example (with periodicity indicated only in the x direction) one gets:<br>
><br>
><br>
> barry-smiths-macbook-pro:tests barrysmith$ petscmpiexec -n 2 ./ex38 -da_view -stencil_width 1 | more<br>
> Processor [0] M 2 N 3 P 4 m 1 n 1 p 2 w 1 s 1<br>
> X range of indices: 0 2, Y range of indices: 0 3, Z range of indices: 0 2<br>
> Processor [1] M 2 N 3 P 4 m 1 n 1 p 2 w 1 s 1<br>
> X range of indices: 0 2, Y range of indices: 0 3, Z range of indices: 2 4<br>
> Vector Object: 1 MPI processes<br>
> type: seq<br>
> 1<br>
> 0<br>
> 1<br>
> 0<br>
> 3<br>
> 2<br>
> 3<br>
> 2<br>
> 5<br>
> 4<br>
> 5<br>
> 4<br>
> 7<br>
> 6<br>
> 7<br>
> 6<br>
> 9<br>
> 8<br>
> 9<br>
> 8<br>
> 11<br>
> 10<br>
> 11<br>
> 10<br>
> -1<br>
> 112<br>
> 113<br>
> -1<br>
> -1<br>
> 114<br>
> 115<br>
> -1<br>
> -1<br>
> 116<br>
> 117<br>
> -1<br>
><br>
> so it is not properly handling the periodicity on the "ghost plane" brought over from the rank 1 process to the rank 0.<br>
><br>
> Note this happens also with a stencil width of 1, so it is not just a question of not handling the stencil width correctly.<br>
><br>
> Barry<br>
><br>
><br>
> On Oct 9, 2011, at 1:32 PM, Richard Katz wrote:<br>
><br>
>><br>
>> Here's how I set up the problem:<br>
>><br>
>> ierr = DMDACreate3d(user->comm,DMDA_BOUNDARY_PERIODIC,DMDA_BOUNDARY_PERIODIC,DMDA_BOUNDARY_NONE,<br>
>> DMDA_STENCIL_BOX,p->Ni,p->Nj,p->Nk,1,1,PETSC_DECIDE,<br>
>> (PetscInt)(sizeof(SolVel3D)/sizeof(PetscReal)),2,<br>
>> PETSC_NULL,PETSC_NULL,PETSC_NULL,&user->dav);CHKERRQ(ierr);<br>
>> ierr = DMCreateGlobalVector(user->dav,&user->v);CHKERRQ(ierr);<br>
>> ierr = DMCreateGlobalVector(user->daa,&user->an);CHKERRQ(ierr);<br>
>> ierr = VecDuplicate(user->v,&user->vo);CHKERRQ(ierr);<br>
>> ierr = VecDuplicate(user->v,&user->vm);CHKERRQ(ierr);<br>
>> ierr = VecDuplicate(user->v,&user->rv);CHKERRQ(ierr);<br>
>> ierr = VecDuplicate(user->an,&user->ano);CHKERRQ(ierr);<br>
>> ierr = SNESCreate(user->comm,&user->snes);CHKERRQ(ierr);<br>
>> if (D2(p)) { ierr = SNESSetFunction(user->snes,user->rv,FormResidual_v2d,user);CHKERRQ(ierr); }<br>
>> else { ierr = SNESSetFunction(user->snes,user->rv,FormResidual_v3d,user);CHKERRQ(ierr); }<br>
>> ierr = DMGetColoring(user->dav,IS_COLORING_GLOBAL,MATAIJ,&iscoloring);<br>
>> ierr = DMGetMatrix(user->dav,MATAIJ,&user->Jv);CHKERRQ(ierr);<br>
>> ierr = MatFDColoringCreate(user->Jv,iscoloring,&user->mfdcv);<br>
>> ierr = MatFDColoringSetFromOptions(user->mfdcv);CHKERRQ(ierr);<br>
>> if (D2(p)) { ierr = MatFDColoringSetFunction(user->mfdcv,(PetscErrorCode(*)(void))FormResidual_v2d,user); }<br>
>> else { ierr = MatFDColoringSetFunction(user->mfdcv,(PetscErrorCode(*)(void))FormResidual_v3d,user); }<br>
>> ierr = SNESSetJacobian(user->snes,user->Jv,user->Jv,SNESDefaultComputeJacobianColor,user->mfdcv);<br>
>> ierr = ISColoringDestroy(&iscoloring);<br>
>> ierr = SNESSetOptionsPrefix(user->snes,"PV_");CHKERRQ(ierr);<br>
>> ierr = SNESSetFromOptions(user->snes);CHKERRQ(ierr);<br>
>><br>
>><br>
>> D2(p) is a macro that returns 1 if the problem is 2D, and 0 if 3D.<br>
>><br>
>> Rich<br>
>><br>
>><br>
>> On 9 Oct 2011, at 19:30, Barry Smith wrote:<br>
>><br>
>>><br>
>>> box stencil or star stencil?<br>
>>><br>
>>><br>
>>> On Oct 9, 2011, at 1:26 PM, Richard Katz wrote:<br>
>>><br>
>>>><br>
>>>> On 9 Oct 2011, at 19:24, Barry Smith wrote:<br>
>>>><br>
>>>>><br>
>>>>> To clarify: With this build valgrind generates no messages at all?<br>
>>>><br>
>>>> That's correct.<br>
>>>><br>
>>>><br>
>>>> [magmox:test]$ $PETSC_DIR/bin/petscmpiexec -valgrind -n 2 ../porbandsim -nd 3 -alpha 2 -theta 45 -output_file test2<br>
>>>> --94613-- run: /usr/bin/dsymutil "../porbandsim"<br>
>>>> --94612-- run: /usr/bin/dsymutil "../porbandsim"<br>
>>>> --------------porbandsim: Sun Oct 9 17:24:24 2011---------------<br>
>>>> Processes: 2, DOFs/Process: 2800<br>
>>>> #PETSc Option Table entries:<br>
>>>> -PV_snes_converged_reason<br>
>>>> -PV_snes_monitor<br>
>>>> -alpha 2<br>
>>>> -ksp_gmres_restart 150<br>
>>>> -nd 3<br>
>>>> -output_file test2<br>
>>>> -snes_converged_reason<br>
>>>> -snes_monitor<br>
>>>> -theta 45<br>
>>>> #End of PETSc Option Table entries<br>
>>>> -----------------------------------------------------------------<br>
>>>> PetscBag Object: ParameterBag params for porosity band simulation<br>
>>>> takei = TRUE; Use Takei & Holztman (2009) viscosity (T/F)<br>
>>>> hperm = FALSE; Use harmonic mean for computing permeability (T/F)<br>
>>>> lambda = 27; Porosity-weakening exponent in shear viscosity, dim'less<br>
>>>> R = 0.57735; Ratio compaction length to domain height, dim'less<br>
>>>> sigma_sat = 1; Saturation stress-difference for anisotropy, dim'less<br>
>>>> n = 3; Exponent in porosity--permeability relationship, dim'less<br>
>>>> phi0 = 0.1; Initial porosity, dim'less<br>
>>>> phi1 = 0.01; Amplitude of porosity perturbation, dim'less<br>
>>>> theta = 45; Fixed anisotropy angle, degrees (-1234567 indicates dynamic)<br>
>>>> alpha = 2; Fixed anisotropy magnitude, 0<=alpha<=2 (-1234567 indicates dynamic)<br>
>>>> viscratio = 1.66667; Ratio of bulk to shear viscosity at phi=phi0, dim'less<br>
>>>> nsig = 1; Exponent in stress--strain rate relationship, dim'less<br>
>>>> init_harmonic = FALSE; Use harmonic perturbation as IC (T/F)<br>
>>>> input_file = null; Name base for input porosity file<br>
>>>> test = FALSE; Just load parameters and then exit<br>
>>>> output_file = test2; Name base for output files<br>
>>>> output_all = FALSE; Make output after each time-step (T/F)<br>
>>>> output_err = TRUE; Make output after failed time-step (T/F)<br>
>>>> output_aux = TRUE; Output auxilliary fields (T/F)<br>
>>>> converged = TRUE; Indicates whether a time-step converged (NO SET)<br>
>>>> nb = 2; Number of buffer rows at top and bottom of domain<br>
>>>> ns = 1000000; Maximum number of time-steps<br>
>>>> nits = 1; Number of Picard iterations<br>
>>>> nplot = 0; Output counter<br>
>>>> advect = 1; Advection scheme: 0=FV, 1=Fromm, 2=TVD<br>
>>>> time = 0; Time at current step, dim'less<br>
>>>> cfl = 0.5; Courant-Friedrichs-Lewy number, dim'less<br>
>>>> tmax = 20; Maximimum simulated time, dim'less<br>
>>>> tout = 0.05; Output interval, dim'less<br>
>>>> nd = 3; Number of spatial dimensions, 2 or 3<br>
>>>> ni = 20; Grid points in x-direction (shear-direction, periodic)<br>
>>>> nj = 5; Grid points in y-direction (cross-direction, periodic)<br>
>>>> nk = 14; Grid points in z-direction, incl buffer points<br>
>>>> -------------------------END PARAM REPORT------------------------<br>
>>>> Continuation in lambda for initial condition<br>
>>>> 0 SNES Function norm 2.512586736188e+01<br>
>>>> 1 SNES Function norm 2.149396609776e-03<br>
>>>> 2 SNES Function norm 1.863405745330e-08<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Solved with lambda = 1.0 in 205.35 seconds<br>
>>>> 0 SNES Function norm 1.874271649438e+00<br>
>>>> 1 SNES Function norm 8.642606950729e-06<br>
>>>> 2 SNES Function norm 7.636114324949e-11<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Solved with lambda = 2.0 in 206.47 seconds<br>
>>>> 0 SNES Function norm 3.748572573117e+00<br>
>>>> 1 SNES Function norm 1.823985563954e-05<br>
>>>> 2 SNES Function norm 1.446658479082e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Solved with lambda = 4.0 in 207.02 seconds<br>
>>>> 0 SNES Function norm 7.498313656718e+00<br>
>>>> 1 SNES Function norm 3.949293531885e-05<br>
>>>> 2 SNES Function norm 3.092503935042e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Solved with lambda = 8.0 in 207.44 seconds<br>
>>>> 0 SNES Function norm 1.500973207738e+01<br>
>>>> 1 SNES Function norm 8.256012403101e-05<br>
>>>> 2 SNES Function norm 6.728525270322e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Solved with lambda = 16.0 in 207.38 seconds<br>
>>>> 0 SNES Function norm 2.066862203903e+01<br>
>>>> 1 SNES Function norm 6.201694512543e-05<br>
>>>> 2 SNES Function norm 5.537777484572e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Solved with lambda = 27.0 in 207.32 seconds<br>
>>>> Generating output file: "test2_0000"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.914413986104e+01<br>
>>>> 1 SNES Function norm 7.835942994856e-05<br>
>>>> 2 SNES Function norm 2.520454050827e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 2.038885321297e+01<br>
>>>> 1 SNES Function norm 1.727646593582e-04<br>
>>>> 2 SNES Function norm 1.325429214752e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.050 in 213.4 seconds<br>
>>>> Generating output file: "test2_0001"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.592319862889e+01<br>
>>>> 1 SNES Function norm 5.602653646480e-05<br>
>>>> 2 SNES Function norm 2.453571869974e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 1.321599310609e+01<br>
>>>> 1 SNES Function norm 1.261023473773e-04<br>
>>>> 2 SNES Function norm 1.163959209038e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.100 in 212.6 seconds<br>
>>>> Generating output file: "test2_0002"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.422948512067e+01<br>
>>>> 1 SNES Function norm 4.242156538283e-05<br>
>>>> 2 SNES Function norm 1.625303993154e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 9.675556760630e+00<br>
>>>> 1 SNES Function norm 1.362209365650e-04<br>
>>>> 2 SNES Function norm 1.223382062021e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.150 in 212.9 seconds<br>
>>>> Generating output file: "test2_0003"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.307890574732e+01<br>
>>>> 1 SNES Function norm 3.214708350215e-05<br>
>>>> 2 SNES Function norm 1.093103013350e-10<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 7.636762082099e+00<br>
>>>> 1 SNES Function norm 1.308928801753e-04<br>
>>>> 2 SNES Function norm 1.143217709967e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.200 in 212.8 seconds<br>
>>>> Generating output file: "test2_0004"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.216229064046e+01<br>
>>>> 1 SNES Function norm 2.441284335468e-05<br>
>>>> 2 SNES Function norm 7.517551623210e-11<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 6.371107649815e+00<br>
>>>> 1 SNES Function norm 1.185031269819e-04<br>
>>>> 2 SNES Function norm 1.095566218658e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.250 in 213.0 seconds<br>
>>>> Generating output file: "test2_0005"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.136564194223e+01<br>
>>>> 1 SNES Function norm 1.865453423897e-05<br>
>>>> 2 SNES Function norm 5.239501925142e-11<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 5.520895480607e+00<br>
>>>> 1 SNES Function norm 1.025932275448e-04<br>
>>>> 2 SNES Function norm 1.071648118984e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.300 in 213.0 seconds<br>
>>>> Generating output file: "test2_0006"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 1.064481900239e+01<br>
>>>> 1 SNES Function norm 1.438574893573e-05<br>
>>>> 2 SNES Function norm 3.643963189108e-11<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 4.905069034291e+00<br>
>>>> 1 SNES Function norm 8.961489700096e-05<br>
>>>> 2 SNES Function norm 1.008198626576e-09<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> Stepped to time t = 0.350 in 212.7 seconds<br>
>>>> Generating output file: "test2_0007"<br>
>>>> -----------------------------------------------------------------<br>
>>>> 0 SNES Function norm 9.981300326125e+00<br>
>>>> 1 SNES Function norm 1.120717450009e-05<br>
>>>> 2 SNES Function norm 2.503394538611e-11<br>
>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>> 0 SNES Function norm 4.428377952757e+00<br>
>>>> ^CCtrl-C caught... cleaning up processes<br>
>>>><br>
>>>><br>
>>>><br>
>>>><br>
>>>>><br>
>>>>> Barry<br>
>>>>><br>
>>>>> On Oct 9, 2011, at 12:19 PM, Richard Katz wrote:<br>
>>>>><br>
>>>>>> Hi Barry,<br>
>>>>>><br>
>>>>>> I have now run the code with<br>
>>>>>><br>
>>>>>> configure_options = [<br>
>>>>>> '--download-mpich=1',<br>
>>>>>> '--with-cxx=0',<br>
>>>>>> '--with-fortran=0',<br>
>>>>>> 'PETSC_ARCH=debug_mpich',<br>
>>>>>> ]<br>
>>>>>><br>
>>>>>> and<br>
>>>>>><br>
>>>>>> $PETSC_DIR/bin/petscmpiexec -valgrind -n 2 ../porbandsim -nd 3 -alpha 2 -theta 45 -output_file test2<br>
>>>>>><br>
>>>>>> The memory error goes away.<br>
>>>>>><br>
>>>>>> The error in the calculated solution does not go away.<br>
>>>>>><br>
>>>>>> Thanks,<br>
>>>>>> Rich<br>
>>>>>><br>
>>>>>><br>
>>>>>><br>
>>>>>><br>
>>>>>> On 9 Oct 2011, at 17:01, Barry Smith wrote:<br>
>>>>>><br>
>>>>>>><br>
>>>>>>> Rich,<br>
>>>>>>><br>
>>>>>>> Could you do an identical run (to the valgrind run you made below) but using another PETSC_ARCH that uses --download-mpich<br>
>>>>>>><br>
>>>>>>> I am suspecting openmpi because it is a bad read of something INSIDE OpenMPI that was just outside the range of memory allocated BY OpenMPI, it is hard to see how PETSc could induce this particular memory problem.<br>
>>>>>>><br>
>>>>>>> Of course, it is also possible there are two bugs, one in OpenMPI and another in PETSc related to the DMDA.<br>
>>>>>>><br>
>>>>>>> Barry<br>
>>>>>>><br>
>>>>>>> On Oct 9, 2011, at 10:46 AM, Richard Katz wrote:<br>
>>>>>>><br>
>>>>>>>> Hi all,<br>
>>>>>>>><br>
>>>>>>>> I'm beginning to wonder if there is a bug in the filling of ghost-points for 3-dimensional DMDAs under petsc-3.2. I am working on a finite-volume code that uses a DMDACreate3d() to create several DMDAs, all with a stencil width of two.<br>
>>>>>>>><br>
>>>>>>>> * The 2-dimensional mode for the same code runs without problems in parallel.<br>
>>>>>>>><br>
>>>>>>>> * The initial SNES solve converges well, with no seg-faults, but when there are processor bounds in the k-direction, the solution is clearly incorrect at/near the boundaries.<br>
>>>>>>>><br>
>>>>>>>> * When I run the same code in serial, or in parallel but with no processor boundaries in the k-direction, it works without problems.<br>
>>>>>>>><br>
>>>>>>>> * When I run the code with processor boundaries in the k-direction, but zeroing-out the bits of the stencil that are at stencil-radius of 2, the processor-boundary errors don't appear.<br>
>>>>>>>><br>
>>>>>>>> * Extending the stencil width to 3 has no effect on the problem.<br>
>>>>>>>><br>
>>>>>>>> * There does seem to be a memory-corruption problem, but I'm not sure where it comes from, because it appears very high in the stack. The following is valgrind output for a run with 2 processors. The processor boundary is in the k-direction. This memory error ONLY appears on the first SNES solve, but not on any of the following ones. Also, the memory error seems to be associated with OpenMPI; it doesn't occur when running with MPICH2<br>
>>>>>>>><br>
>>>>>>>> [magmox:test]$ $PETSC_DIR/bin/petscmpiexec -valgrind -n 2 ../porbandsim -malloc off -nd 3 -alpha 2 -theta 45 -output_file test2<br>
>>>>>>>><br>
>>>>>>>> [...]<br>
>>>>>>>><br>
>>>>>>>> Continuation in lambda for initial condition<br>
>>>>>>>> 0 SNES Function norm 2.512586736188e+01<br>
>>>>>>>> 1 SNES Function norm 2.149396609776e-03<br>
>>>>>>>> ==4576== Invalid read of size 8<br>
>>>>>>>> ==4576== at 0x7FFFFFE00BAC: ???<br>
>>>>>>>> ==4576== by 0x101A34213: ompi_convertor_pack (in /usr/lib/libmpi.0.dylib)<br>
>>>>>>>> ==4576== by 0x10858909E: mca_pml_ob1_send_request_start_copy (in /usr/lib/openmpi/mca_pml_ob1.so)<br>
>>>>>>>> ==4576== by 0x1085840D2: mca_pml_ob1_isend (in /usr/lib/openmpi/mca_pml_ob1.so)<br>
>>>>>>>> ==4576== by 0x1085D11C3: ompi_coll_tuned_allreduce_intra_recursivedoubling (in /usr/lib/openmpi/mca_coll_tuned.so)<br>
>>>>>>>> ==4576== by 0x101A4419B: MPI_Allreduce (in /usr/lib/libmpi.0.dylib)<br>
>>>>>>>> ==4576== by 0x1001D3F3E: VecMDot_MPI (pvec2.c:20)<br>
>>>>>>>> ==4576== by 0x100204669: VecMDot (rvector.c:1144)<br>
>>>>>>>> ==4576== by 0x10097E28B: KSPGMRESClassicalGramSchmidtOrthogonalization (borthog2.c:66)<br>
>>>>>>>> ==4576== by 0x100975B33: GMREScycle (gmres.c:159)<br>
>>>>>>>> ==4576== by 0x10097656E: KSPSolve_GMRES (gmres.c:231)<br>
>>>>>>>> ==4576== by 0x1009A8606: KSPSolve (itfunc.c:423)<br>
>>>>>>>> ==4576== by 0x100A0271A: SNES_KSPSolve (snes.c:3396)<br>
>>>>>>>> ==4576== by 0x100A17F79: SNESSolve_LS (ls.c:190)<br>
>>>>>>>> ==4576== by 0x1009FB635: SNESSolve (snes.c:2676)<br>
>>>>>>>> ==4576== by 0x100008104: InitSolveElliptic (initialize.c:216)<br>
>>>>>>>> ==4576== by 0x10000538A: FormInitialCondition (initialize.c:31)<br>
>>>>>>>> ==4576== by 0x100001410: main (porbandsim.c:28)<br>
>>>>>>>> ==4576== Address 0x10289f778 is 0 bytes after a block of size 136 alloc'd<br>
>>>>>>>> ==4576== at 0x1010B05BF: malloc (vg_replace_malloc.c:266)<br>
>>>>>>>> ==4576== by 0x1085D0F98: ompi_coll_tuned_allreduce_intra_recursivedoubling (in /usr/lib/openmpi/mca_coll_tuned.so)<br>
>>>>>>>> ==4576== by 0x101A4419B: MPI_Allreduce (in /usr/lib/libmpi.0.dylib)<br>
>>>>>>>> ==4576== by 0x1001D3F3E: VecMDot_MPI (pvec2.c:20)<br>
>>>>>>>> ==4576== by 0x100204669: VecMDot (rvector.c:1144)<br>
>>>>>>>> ==4576== by 0x10097E28B: KSPGMRESClassicalGramSchmidtOrthogonalization (borthog2.c:66)<br>
>>>>>>>> ==4576== by 0x100975B33: GMREScycle (gmres.c:159)<br>
>>>>>>>> ==4576== by 0x10097656E: KSPSolve_GMRES (gmres.c:231)<br>
>>>>>>>> ==4576== by 0x1009A8606: KSPSolve (itfunc.c:423)<br>
>>>>>>>> ==4576== by 0x100A0271A: SNES_KSPSolve (snes.c:3396)<br>
>>>>>>>> ==4576== by 0x100A17F79: SNESSolve_LS (ls.c:190)<br>
>>>>>>>> ==4576== by 0x1009FB635: SNESSolve (snes.c:2676)<br>
>>>>>>>> ==4576== by 0x100008104: InitSolveElliptic (initialize.c:216)<br>
>>>>>>>> ==4576== by 0x10000538A: FormInitialCondition (initialize.c:31)<br>
>>>>>>>> ==4576== by 0x100001410: main (porbandsim.c:28)<br>
>>>>>>>> ==4576==<br>
>>>>>>>> ==4577== Invalid read of size 8<br>
>>>>>>>> ==4577== at 0x7FFFFFE00BAC: ???<br>
>>>>>>>> ==4577== by 0x101A34213: ompi_convertor_pack (in /usr/lib/libmpi.0.dylib)<br>
>>>>>>>> ==4577== by 0x10858909E: mca_pml_ob1_send_request_start_copy (in /usr/lib/openmpi/mca_pml_ob1.so)<br>
>>>>>>>> ==4577== by 0x1085840D2: mca_pml_ob1_isend (in /usr/lib/openmpi/mca_pml_ob1.so)<br>
>>>>>>>> ==4577== by 0x1085D11C3: ompi_coll_tuned_allreduce_intra_recursivedoubling (in /usr/lib/openmpi/mca_coll_tuned.so)<br>
>>>>>>>> ==4577== by 0x101A4419B: MPI_Allreduce (in /usr/lib/libmpi.0.dylib)<br>
>>>>>>>> ==4577== by 0x1001D3F3E: VecMDot_MPI (pvec2.c:20)<br>
>>>>>>>> ==4577== by 0x100204669: VecMDot (rvector.c:1144)<br>
>>>>>>>> ==4577== by 0x10097E28B: KSPGMRESClassicalGramSchmidtOrthogonalization (borthog2.c:66)<br>
>>>>>>>> ==4577== by 0x100975B33: GMREScycle (gmres.c:159)<br>
>>>>>>>> ==4577== by 0x10097656E: KSPSolve_GMRES (gmres.c:231)<br>
>>>>>>>> ==4577== by 0x1009A8606: KSPSolve (itfunc.c:423)<br>
>>>>>>>> ==4577== by 0x100A0271A: SNES_KSPSolve (snes.c:3396)<br>
>>>>>>>> ==4577== by 0x100A17F79: SNESSolve_LS (ls.c:190)<br>
>>>>>>>> ==4577== by 0x1009FB635: SNESSolve (snes.c:2676)<br>
>>>>>>>> ==4577== by 0x100008104: InitSolveElliptic (initialize.c:216)<br>
>>>>>>>> ==4577== by 0x10000538A: FormInitialCondition (initialize.c:31)<br>
>>>>>>>> ==4577== by 0x100001410: main (porbandsim.c:28)<br>
>>>>>>>> ==4577== Address 0x10289f778 is 0 bytes after a block of size 136 alloc'd<br>
>>>>>>>> ==4577== at 0x1010B05BF: malloc (vg_replace_malloc.c:266)<br>
>>>>>>>> ==4577== by 0x1085D0F98: ompi_coll_tuned_allreduce_intra_recursivedoubling (in /usr/lib/openmpi/mca_coll_tuned.so)<br>
>>>>>>>> ==4577== by 0x101A4419B: MPI_Allreduce (in /usr/lib/libmpi.0.dylib)<br>
>>>>>>>> ==4577== by 0x1001D3F3E: VecMDot_MPI (pvec2.c:20)<br>
>>>>>>>> ==4577== by 0x100204669: VecMDot (rvector.c:1144)<br>
>>>>>>>> ==4577== by 0x10097E28B: KSPGMRESClassicalGramSchmidtOrthogonalization (borthog2.c:66)<br>
>>>>>>>> ==4577== by 0x100975B33: GMREScycle (gmres.c:159)<br>
>>>>>>>> ==4577== by 0x10097656E: KSPSolve_GMRES (gmres.c:231)<br>
>>>>>>>> ==4577== by 0x1009A8606: KSPSolve (itfunc.c:423)<br>
>>>>>>>> ==4577== by 0x100A0271A: SNES_KSPSolve (snes.c:3396)<br>
>>>>>>>> ==4577== by 0x100A17F79: SNESSolve_LS (ls.c:190)<br>
>>>>>>>> ==4577== by 0x1009FB635: SNESSolve (snes.c:2676)<br>
>>>>>>>> ==4577== by 0x100008104: InitSolveElliptic (initialize.c:216)<br>
>>>>>>>> ==4577== by 0x10000538A: FormInitialCondition (initialize.c:31)<br>
>>>>>>>> ==4577== by 0x100001410: main (porbandsim.c:28)<br>
>>>>>>>> ==4577==<br>
>>>>>>>> 2 SNES Function norm 1.863405745330e-08<br>
>>>>>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>>>>>> 0 SNES Function norm 1.874271649438e+00<br>
>>>>>>>> 1 SNES Function norm 8.642606950729e-06<br>
>>>>>>>> 2 SNES Function norm 7.636114324949e-11<br>
>>>>>>>> Nonlinear solve converged due to CONVERGED_FNORM_ABS<br>
>>>>>>>><br>
>>>>>>>> [...]<br>
>>>>>>>><br>
>>>>>>>><br>
>>>>>>>><br>
>>>>>>>> Any suggestions about how to debug this? How well-tested are the parallel, 3D DMDA routines? Should I do some testing on a PETSc example code?<br>
>>>>>>>><br>
>>>>>>>> I am running petsc-3.2 with<br>
>>>>>>>><br>
>>>>>>>> changeset: 20338:fbd00cf1070e<br>
>>>>>>>> tag: tip<br>
>>>>>>>> parent: 20337:589176156279<br>
>>>>>>>> parent: 20336:7abcced43e69<br>
>>>>>>>> user: Barry Smith <<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>><br>
>>>>>>>> date: Tue Oct 04 16:45:00 2011 -0500<br>
>>>>>>>> summary: commit after merge<br>
>>>>>>>><br>
>>>>>>>> and<br>
>>>>>>>><br>
>>>>>>>> configure_options = [<br>
>>>>>>>> '--download-fftw=true',<br>
>>>>>>>> '--with-cxx=0',<br>
>>>>>>>> '--with-debugging=1',<br>
>>>>>>>> '--with-fortran=0',<br>
>>>>>>>> '--with-mpi-dir=/usr',<br>
>>>>>>>> 'PETSC_ARCH=debug',<br>
>>>>>>>> ]<br>
>>>>>>>><br>
>>>>>>>> or<br>
>>>>>>>><br>
>>>>>>>> configure_options = [<br>
>>>>>>>> '--download-blacs',<br>
>>>>>>>> '--download-fftw=true',<br>
>>>>>>>> '--download-mpich=yes',<br>
>>>>>>>> '--download-mumps=yes',<br>
>>>>>>>> '--download-parmetis',<br>
>>>>>>>> '--download-scalapack',<br>
>>>>>>>> '--with-cxx=0',<br>
>>>>>>>> '--with-debugging=0',<br>
>>>>>>>> 'PETSC_ARCH=optimize_mumps',<br>
>>>>>>>> ]<br>
>>>>>>>><br>
>>>>>>>><br>
>>>>>>>> Thanks,<br>
>>>>>>>> Rich<br>
>>>>>>>><br>
>>>>>>>><br>
>>>>>>><br>
>>>>>><br>
>>>>>> ________________________________<br>
>>>>>> Richard Foa Katz<br>
>>>>>> Dept Earth Sciences, Univ Oxford<br>
>>>>>> <a href="http://foalab.earth.ox.ac.uk" target="_blank">http://foalab.earth.ox.ac.uk</a><br>
>>>>>><br>
>>>>><br>
>>>><br>
>>>> ________________________________<br>
>>>> Richard Foa Katz<br>
>>>> Dept Earth Sciences, Univ Oxford<br>
>>>> <a href="http://foalab.earth.ox.ac.uk" target="_blank">http://foalab.earth.ox.ac.uk</a><br>
>>>><br>
>>><br>
>><br>
>> ________________________________<br>
>> Richard Foa Katz<br>
>> Dept Earth Sciences, Univ Oxford<br>
>> <a href="http://foalab.earth.ox.ac.uk" target="_blank">http://foalab.earth.ox.ac.uk</a><br>
>><br>
><br>
<br>
</blockquote></div>