<div dir="ltr">Hi Klaus, and others listening. Thanks for your help so far.<div><br></div><div>   I've investigated a bit further. Some confusing things are happening. </div><div><br></div><div>  I made a new build that should amount to the same initial conditions, with density set everywhere to 4g/cc, tele=trad=100K, box on [1.5e7,-3.e8,-1.5e8]:[6.15e8,3.e8,1.5e8], with lrefine_min=lrefine_max=4, and 16=NXB=NYB=2NZB. To make sure it wasn't due to bad compiling with too high an optimization, I had -O0 throughout for this new build (but -O3 for previous communication). </div><div><br></div><div>  Looking at the physical quantities in the 0000 checkpoint file, everything was the same with this new build compared with the old run except for gpot. In the Original build posted before (which is based on my existing code which then overwrote every cell at init with rho,T=4,100) I was getting gpot= -3.06777770730863e+27 in the first cell and it doesn't vary much (dgpot/gpot ~ 1.e-10). </div><div><br></div><div>  In my new build I was getting gpot=2760347382.52056 and varied highly from cell to cell. </div><div><br></div><div>  A quick calculation of replacing each cell with a pointmass makes the latter seems close to what I'd expect. </div><div><br></div><div>In both of these cases the gravitational boundary condition was set to isolated. However I noticed that in my new build I had my x and y hydro boundaries as periodic while in my original run I had all boundaries as diode. So I reran the new build with all diode and now I get: +5.97937076239123e26. This seemed really odd; why would changing the fluid boundaries, but not the gravitational boundaries nor the density or cell spacing, change the gravitational potential? Also, once the potential was at this magnitude, then HYPRE was suffering from the non-convergence case again.</div><div><br></div><div>I then reran my original run with -O0 (instead of -O3 which I had originally) and now it also gets +5.97937076239123e26.</div><div><br></div><div>I then reran my new build with -O3 (instead of -O0) and now it gets -3.06777770730863e+27.</div><div><br></div><div>So: Once everything was the same in the new build, then they agree, whether with -O0 or -O3. However :</div><div><br></div><div>  1) Why would changing only the hydro boundary conditions impact gpot (and by ~20 orders of magnitude). </div><div>  2) Has anyone else seen changing the optimization level impact the gpot solution.  </div><div>  3) Clearly these values are much larger than expected: as I said a simple N-body calculation gets me around -1e11 for the potential. Why is this? And does it explain the previous issues I was seeing with HYPRE large RHS?</div><div><br></div><div>Cheers,</div><div>   -Mark</div><div><br></div><div>PS:  Other relevant information: </div><br><p class="">Number of MPI tasks:                 32<br> MPI version:                          3<br> MPI subversion:                       0<br> Dimensionality:                       3<br> Max Number of Blocks/Proc:          400<br> Number x zones:                      16<br> Number y zones:                      16<br> Number z zones:                       8<br></p><div class="gmail_extra">System info:     <br>Linux <a href="http://login3.stampede.tacc.utexas.edu">login3.stampede.tacc.utexas.edu</a> 2.6.32-431.17.1.el6.x86_64 #1 SMP Wed May<br> Version:</div><div class="gmail_extra"> FLASH 4.3_release                                                              </div><div class="gmail_extra"><br></div><div class="gmail_extra">FLASH Units used:<br>   Driver/DriverMain/Unsplit<br>   Grid/GridAllocatableScratches<br>   Grid/GridBoundaryConditions<br>   Grid/GridMain/paramesh/interpolation/Paramesh4/prolong<br>   Grid/GridMain/paramesh/interpolation/prolong<br>   Grid/GridMain/paramesh/paramesh4/Paramesh4dev/PM4_package/headers<br>   Grid/GridMain/paramesh/paramesh4/Paramesh4dev/PM4_package/mpi_source<br>   Grid/GridMain/paramesh/paramesh4/Paramesh4dev/PM4_package/source<br>   Grid/GridMain/paramesh/paramesh4/Paramesh4dev/PM4_package/utilities/multigrid<br>   Grid/GridSolvers/HYPRE/paramesh<br>   Grid/GridSolvers/IsoBndMultipole<br>   Grid/GridSolvers/Multigrid/fft<br>   IO/IOMain/hdf5/serial/PM<br>   Multispecies/MultispeciesMain<br>   PhysicalConstants/PhysicalConstantsMain<br>   RuntimeParameters/RuntimeParametersMain<br>   Simulation/SimulationMain/ChondruleMappingMSMT_NoTill_ConstOpac<br>   flashUtilities/contiguousConversion<br>   flashUtilities/general<br>   flashUtilities/interpolation/oneDim<br>   flashUtilities/nameValueLL<br>   flashUtilities/sorting/quicksort<br>   flashUtilities/system/memoryUsage/legacy<br>   monitors/Logfile/LogfileMain<br>   monitors/Timers/TimersMain/MPINative<br>   physics/Diffuse/DiffuseMain/Unsplit<br>   physics/Eos/EosMain/multiTemp/MatRad3/Multigamma<br>   physics/Eos/EosMain/multiTemp/Multigamma<br>   physics/Gravity/GravityMain/Poisson/Multigrid<br>   physics/Hydro/HydroMain/unsplit_rad/Hydro_Unsplit<br>   physics/Hydro/HydroMain/unsplit_rad/multiTemp<br>   physics/RadTrans/RadTransMain/MGD<br>   physics/materialProperties/Conductivity/ConductivityMain/Constant<br>   physics/materialProperties/Opacity/OpacityMain/Constcm2g<div class="gmail_quote"><br></div>Setup options: -auto -3d -nxb=16 -nyb=16 -nzb=8 -maxblocks=400 species=rock,watr +uhd3tr mgd_meshgroups=1<br><div class="gmail_quote"><br></div><div class="gmail_quote"><br></div><div class="gmail_quote">On 19 February 2016 at 00:15, Klaus Weide <span dir="ltr"><<a href="mailto:klaus@flash.uchicago.edu" target="_blank">klaus@flash.uchicago.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Hello Mark,<br>
<br>
Some further comments on your last message blow.<br>
<br>
<br>
On Mon, 15 Feb 2016, Mark Richardson wrote:<br>
<br>
 .....<br>
<span class="">><br>
> When I was filling the guard cells, this was accomplished in<br>
> Simulation_initBlock, using the exact same method as filling the rest of<br>
</span>> the cells. Thus, using *gcell = .true., *the coordinates of all cells was<br>
<div><div class="h5">> found with:<br>
><br>
> ! get the integer index information for the current block<br>
>   call Grid_getBlkIndexLimits(blockId,blkLimits,blkLimitsGC)<br>
><br>
>   ! Determine size of interior of block<br>
>   sizeX = blkLimitsGC(HIGH,IAXIS) - blkLimitsGC(LOW,IAXIS) + 1<br>
>   sizeY = blkLimitsGC(HIGH,JAXIS) - blkLimitsGC(LOW,JAXIS) + 1<br>
>   sizeZ = blkLimitsGC(HIGH,KAXIS) - blkLimitsGC(LOW,KAXIS) + 1<br>
><br>
>   ! allocate memory for iCen variables<br>
>   allocate(xCen(sizeX),stat=istat)<br>
>   allocate(yCen(sizeY),stat=istat)<br>
>   allocate(zCen(sizeZ),stat=istat)<br>
><br>
>   ! Initialize<br>
>   xCen = 0.0d0<br>
>   yCen = 0.0d0<br>
>   zCen = 0.0d0<br>
><br>
>   ! Get interior cell coordinates<br>
>   if (NDIM == 3) &<br>
>     call Grid_getCellCoords(KAXIS, blockId, CENTER ,gcell, zCen, sizeZ)<br>
>   if (NDIM >= 2) &<br>
>     call Grid_getCellCoords(JAXIS, blockId, CENTER,gcell, yCen, sizeY)<br>
>   call Grid_getCellCoords(IAXIS, blockId, CENTER,gcell, xCen, sizeX)<br>
><br>
>   ! Get lower and upper coordinates<br>
>   blkLow(1)   = xCen(1)<br>
>   blkHigh(1) = xCen(sizeX)<br>
><br>
>   blkLow(2)   = yCen(1)<br>
>   blkHigh(2) = yCen(sizeY)<br>
><br>
>   blkLow(3)   = zCen(1)<br>
>   blkHigh(3) = zCen(sizeZ)<br>
><br>
> ...<br>
><br>
> Then later, I fill each cell in the loop led with:<br>
><br>
>   do k = blkLimitsGC(LOW,KAXIS),blkLimitsGC(HIGH,KAXIS)<br>
>      do j = blkLimitsGC(LOW,JAXIS),blkLimitsGC(HIGH,JAXIS)<br>
>         do i = blkLimitsGC(LOW,IAXIS),blkLimitsGC(HIGH,IAXIS)<br>
<br>
</div></div>Yes, that looks correct (if gcell==.true.). However, ...<br>
<br>
> Originally *gcell* was *.false.* and the *size{X,Y,Z}* variables were set<br>
> by the non-GC versions, and then the *k,j,i* loops were over the non GC<br>
<span class="">> versions. By switching to the GC version, the small CFL went away.<br>
<br>
</span>When you replace blkLimitsGC with blkLimits everywhere (after the<br>
Grid_getBlkIndexLimits call), and set gcell=.false., to restore your<br>
"original" code, this may become wrong.  It depends on whether you use<br>
the xCen / yCen / zCen arrays anywhere in your loop, indexed by i / j / k.<br>
That's because<br>
<br>
 * with blkLimitsGC : xCen(1) corresponds to blkData(1,j,k), etc.<br>
 * with blkLimitsGC : xCen(1) corresponds to blkData(5,j,k), etc.<br>
   (assuming NGUARD=4).<br>
<br>
So you would have to replace any xCen(i) in the loop body with<br>
xCen(i-NGUARD), etc.<br>
<br>
(alternatively, change the bounds of the allocated coordiante arrays<br>
so that the numbering agrees with that used by the look indices i,j,k -<br>
I.e., replace<br>
<br>
>   allocate(xCen(sizeX),stat=istat)<br>
<br>
with<br>
<br>
>   allocate(xCen(blkLimits(LOW,IAXIS):blkLimits(HIGH,IAXIS)))<br>
<br>
etc. (and appropriate changes where you set your blkLow and blkHigh).<br>
<br>
Maybe this is obvious to you, but I thought it might be related to why<br>
your original version was misbehaving.<br>
<span class=""><br>
<br>
<br>
> I'm having one other problem, which is maybe related. Originally I was<br>
> running this with a self-made Eos (Tillotson for planetary bodies, that<br>
> includes radiative transfer). I was getting the following errors<br>
> immediately after initialisation:<br>
><br>
>  [ 02-15-2016  09:57:48.810 ] [gr_hypreInit]: HYPRE_RELEASE_NAME=hypre,<br>
> HYPRE_RELEASE_VERSION=2.10.0b, HYPRE_RELEASE_DATE=2015/01/22<br>
<br>
</span>.....<br>
<span class=""><br>
>  [ 02-15-2016  09:57:58.820 ] [gr_hypreSolve]: Nonconv./failure: ierr=<br>
>       0<br>
>  [ 02-15-2016  09:57:58.820 ] [gr_hypreSolve]: Nonconv./failure: component=<br>
>           0<br>
>  [ 02-15-2016  09:57:58.820 ] [gr_hypreSolve]: Nonconv./failure: converged=<br>
> F<br>
>  [ 02-15-2016  09:57:58.820 ] [gr_hypreSolve]: Nonconv./failure: |initial<br>
> RHS|^2=  1.2520E+35<br>
>  [ 02-15-2016  09:57:58.820 ] [gr_hypreSolve]: Nonconv./failure:<br>
> num_iterations=         500<br>
>  [ 02-15-2016  09:57:58.820 ] [gr_hypreSolve]: Nonconv./failure:<br>
> final_res_norm=   88143026.2676523<br>
<br>
</span>Not sure what is going on here.<br>
<br>
Hopefully it's not the HYPRE version (I think we still haven't tested with 2.10.0b).<br>
<br>
The high |initial RHS^2| value may indicate that you are feeding HYPRE<br>
an unphysically large r.h.s.  (The B in A X = B); whether 1.2520E+35<br>
should be considered unphysical depenss on your problem (initial conditions,<br>
number of cells...)<br>
<span class=""><br>
<br>
> Where hypre isn't converging. I'm not sure why this would be, as I<br>
> eventually tested it with a uniform density and temperature (with temp_ele<br>
> = temp_rad) for the MatRad3 case of no ions. Also, I'm surprised Flash<br>
> doesn't crash when non-convergence occurs in the hypre grid solver.<br>
<br>
</span>Maybe FLASH should abort.  But sometimes we have failures in HYPRE to<br>
be transitory, and to disappear a few time steps.  In case of severe<br>
and persisting problems - they usually cause a serious blowup in a<br>
variable which will normally not go undetected for long.<br>
<span class=""><br>
> I've now tested it with the MatRad3 Multigamma EOS, and get the same error,<br>
> so I'm unsure how this would be related to my EOS. Has anyone had similar<br>
> issues with the hypre grid solver?<br>
><br>
> These may be the most relevant parts of my flash.par file:<br>
> # Flux control<br>
> hy_fPresInMomFlux = 1.0      # Percentage of pressure gradient for the<br>
<br>
</span>You should set the above to 0.0 if you are using the unsplit_rad Hydro code.<br>
(If not, leave at at 1.0.)<br>
<br>
If you are using MatRad3 EOS but with the "unsplit" (not _rad) Hydro code,<br>
you are testing a combination that is MEANT to work, bu we haven't tested it much.<br>
<span class=""><br>
<br>
> ppmEintFluxConstructionMeth     = 0<br>
> ppmEintCompFluxConstructionMeth = 0 #(or 1) #Advect like any old mass scalar<br>
> ppmEnerFluxConstructionMeth     = 0<br>
<br>
</span>Those three dont't matter.<br>
<span class=""><font color="#888888"><br>
Klaus<br>
</font></span></blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><br>Mark Richardson<br><div>Beecroft Postdoctoral Fellow, Department of Physics<div>Denys Wilkinson Building - 555F</div><div>University of Oxford<div><a href="mailto:Mark.Richardson@physics.ox.ac.uk" target="_blank">Mark.Richardson@physics.ox.ac.uk</a><br><div>+44 1865 280793</div></div></div></div></div></div></div></div></div></div></div></div></div></div></div></div>
</div></div>