[FLASH-USERS] DustCollapse Weird Artifacts in Gravitational Acceleration

[Joshua Wall]

Hello Ryan,

      I double checked this morning, and Richard Wunsch's BHTree gravity
includes variables for storing the gravitation acceleration at each cell if
you use his solver with the added option bhtreeAcc at setup like this (make
sure you take out any references to other gravity solvers in the simulation
Config file):

./setup DustCollapse +auto +3d
--with-unit=physics/Gravity/GravityMain/Poisson/BHTree bhtreeAcc=1

If you check in the code in
source/physics/Gravity/GravityMain/Poisson/BHTree/Gravity_accelOneRow.F90
you can see how the accelerations are calculated and stored exactly from
both the grid and any sinks in the code. I've been using his potential
solver for some time (in Flash 4.2.2, including getting the accelerations)
and find it to be both robust and fast. Or perhaps you could use his
implementation as an example of how to structure your own modifications to
Gravity_accelOneRow.F90.

Note that if you do use his code it already includes the gravity accel
variables in source/physics/Gravity/GravityMain/Poisson/BHTree/Config,
which I encourage you to check out as it describes the other options for
the code. The paper for this method is in prep still I believe.

Hope this helps!

Josh

On Sat, Feb 25, 2017 at 10:03 PM Ryan Farber <rjfarber at umich.edu> wrote:

> Dear FLASH users,
>
> I decided to look at the components of the gravitational acceleration for
> the supplied
> DustCollapse problem (FLASH 4.3), and have noticed some weird artifacts on
> the block boundaries.
>
> The issue seems somewhat similar to one Thomas Peters' refinement issue
> <http://flash.uchicago.edu/pipermail/flash-users/2014-May/001440.html> (a
> plotting bug)
> but I chatted with ngoldbaum on the yt IRC channel and the conclusion was
> the artifacts
> are on the FLASH side of things.
>
> The only edits I made to the DustCollapse problem supplied with FLASH 4.3
> were:
>
> * Config: added the following lines
> VARIABLE grac
> VARIABLE poix
> VARIABLE poiy
> VARIABLE poiz
>
> * Driver_initFlash.F90: copied to problem directory and modified the
> original lines:
>   if(.not. dr_restart) then
>      call Grid_getListOfBlocks(LEAF,blockList,blockCount)
>      call Gravity_potentialListOfBlocks(blockCount,blockList)
>      call Particles_initForces()
>   end if
>
>   ! If we want to free any arrays created during simulation
>   ! initialization that are no longer needed, do it here.
>   call Simulation_freeUserArrays()
>
>   call IO_outputInitial(  dr_nbegin, dr_initialSimTime)
>   if(dr_globalMe==MASTER_PE)print*,'Initial plotfile written'
>
>   if(dr_globalMe==MASTER_PE)print*,'Driver init all done'
>
>
> **to**
>   if(.not. dr_restart) then
>      call Grid_getListOfBlocks(LEAF,blockList,blockCount)
>      call Gravity_potentialListOfBlocks(blockCount,blockList)
>      call Particles_initForces()
>
>      call Gravity_accelListOfBlocks(blockCount, blockList, IAXIS, GRAC_VAR)
>      call Gravity_accelListOfBlocks(blockCount, blockList, JAXIS, GRAC_VAR)
>      call Gravity_accelListOfBlocks(blockCount, blockList, KAXIS, GRAC_VAR)
>   end if
>
>   ! If we want to free any arrays created during simulation
>   ! initialization that are no longer needed, do it here.
>   call Simulation_freeUserArrays()
>
>   call IO_outputInitial(  dr_nbegin, dr_initialSimTime)
>   if(dr_globalMe==MASTER_PE)print*,'Initial plotfile written'
>
>   STOP(99)
>
>   if(dr_globalMe==MASTER_PE)print*,'Driver init all done'
>
> * Gravity_accelOneRow: modified the original lines:
>   do ii = iimin+1, iimax-1
>      grav(ii) = grav(ii) + delxinv * (gpot(ii-1) - gpot(ii+1))
>   enddo
>
>   grav(iimin) = grav(iimin+1)     ! this is invalid data - must not be used
>   grav(iimax) = grav(iimax-1)
>
> **to**
>   do ii = iimin+1, iimax-1
>      grav(ii) = grav(ii) + delxinv * (gpot(ii-1) - gpot(ii+1))
>
>      if (sweepDir .eq. SWEEP_X) then
>        solnVec(POIX_VAR, ii, pos(1), pos(2)) = grav(ii)
>      end if
>      if (sweepDir .eq. SWEEP_Y) then
>        solnVec(POIY_VAR, pos(1), ii, pos(2)) = grav(ii)
>      end if
>      if (sweepDir .eq. SWEEP_Z) then
>        solnVec(POIZ_VAR, pos(1), pos(2), ii) = grav(ii)
>      end if
>   enddo
>
>   grav(iimin) = grav(iimin+1)     ! this is invalid data - must not be used
>   grav(iimax) = grav(iimax-1)
>   if (sweepDir .eq. SWEEP_X) then
>     solnVec(POIX_VAR, iimin, pos(1), pos(2)) = grav(iimin+1)
>     solnVec(POIX_VAR, iimax, pos(1), pos(2)) = grav(iimax-1)
>   end if
>   if (sweepDir .eq. SWEEP_Y) then
>     solnVec(POIY_VAR, pos(1), iimin, pos(2)) = grav(iimin+1)
>     solnVec(POIY_VAR, pos(1), iimax, pos(2)) = grav(iimax-1)
>   end if
>   if (sweepDir .eq. SWEEP_Z) then
>     solnVec(POIZ_VAR, pos(1), pos(2), iimin) = grav(iimin+1)
>     solnVec(POIZ_VAR, pos(1), pos(2), iimax) = grav(iimax-1)
>   end if
>
> Slices were taken along the x-axis.
>
> Thanks for your help,
> Ryan
>
-- 
Joshua Wall
Doctoral Candidate
Department of Physics
Drexel University
3141 Chestnut Street
Philadelphia, PA 19104
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://flash.rochester.edu/pipermail/flash-users/attachments/20170227/c126d1c4/attachment.htm>