[FLASH-USERS] Gravity, GPOL_VAR, and GPOT_VAR

[Cole Miller]

Hi all,

I've been trying to modify the Cellular test problem to include
gravity (ultimately with success, I think, but see below).  Along
the way I ran into a couple of things that might be useful to have
on record.  I'm using a version of FLASH 3.0 beta that I downloaded
on September 11, 2007, so perhaps these have been changed already.

1. There is a typo in Simulation_initBlock.F90 in
   source/Simulation/SimulationMain/Cellular.  In this particular
   test problem, one sets up an ambient medium with a perturbation
   in it, and watches the burning progress.  In the 
   Simulation_initBlock.F90 code we have the lines
           else if (NDIM .EQ. 2) then
              if (sim_usePseudo1d) then
                 dist = xx - sim_xCenterPerturb
              else
                 dist = sqrt((xx - sim_xCenterPerturb)**2 + &
                      &                 (xx - sim_yCenterPerturb)**2)
              endif

   Clearly, the second to last line should have
   (yy - sim_yCenterPerturb)**2
   instead.

2. On a much more complicated note, for a long time I was unsuccessful
   in adding constant gravity to this problem.  I did all the obvious 
   things, e.g., adding REQUIRES physics/Gravity/ to the Config file,
   along with gconst and gdirec, making sure that there was a line
   "useGravity = .true." in flash.par, running it in three dimensions
   instead  of two, and so on.  No luck.  Eventually, with the help of
   Sean O'Neill here at Maryland, I tracked down the problem.

   In source/physics/Hydro/HydroMain/split/PPM/PPMKernel/
   there is a program called hydro_1d.F90.  This program reads the
   gravity and applies the force, if any, to the fluid.  In it, we have the
   following lines:
   #ifdef GPOL_VAR
   	  call Gravity_accelOneRow (pos, xyzswp, blockID, numIntCells8,ograv,GPOL_VAR)
   #ifndef GPOT_VAR
   	  call Driver_abortFlash("Shouldn't have gpol defined without gpot")
   #endif
   #endif
   #ifdef GPOT_VAR
   	  call Gravity_accelOneRow (pos, xyzswp, blockID, numIntCells8,grav,GPOT_VAR)
   #endif
   	  dtfac = dt/dt_old

   	  do i = 1,numIntCells8
   	     dg       = dtfac*(grav(i) - ograv(i))
   	     hgrav(i) = grav(i) + 0.5e0*dg
   	     ngrav(i) = grav(i) +	dg
   	  enddo

   and then hgrav is passed to the Riemann solver as a force.

   The problem is that nowhere that I saw is GPOL_VAR or GPOT_VAR defined.
   Not in Flash.h, not in any of the files after gmake is run.  I tried
   several variants of the setup script, with the same result.  Therefore,
   since prior to the lines quoted above the arrays grav() and ograv() are
   both set to zero, hgrav is as well, so there is no gravity.

   However, on further examination, in the program Gravity_accelOneRow
   the last variable (GPOL_VAR or GPOT_VAR here) is optional, and is in
   fact not used at all for constant gravity.  Those two variables aren't
   used anywhere else, either.  Therefore, just before the dtfac line I added:

   #ifndef GPOT_VAR
   	  call Gravity_accelOneRow (pos, xyzswp, blockID, numIntCells8,grav,-1)
   #endif

   This appears to work.  

   Is there a more elegant solution?  Is there, for example, a specific
   option in ./setup that will define GPOL_VAR and/or GPOT_VAR?  If not,
   perhaps it would be useful to look into this.

Cole