[FLASH-USERS] energy inaccuracy in Sedov simulation

Joshua Wall joshua.e.wall at gmail.com
Thu Jul 28 07:56:00 EDT 2016


Hello everyone,

    I got several responses about this bit of code, so we have decided to
share it with the community at large (for science!) so everyone can use it.
It is a simple "uniform sampling" Monte Carlo integration of the overlap
between a sphere or a right cylinder, as was first implemented in ZEUS-MP
by David Clarke. Note that you could likely make this faster and/or more
accurate by making the sampling random. There are also a few known tricks
for random sampling in the "right" region, check your local Numerical
Methods text under MC integration. But especially when you just want to
initialize once or just want an approximation, this method is very nice. I
hope everyone finds it as useful as I have.

Cordially,

Joshua Wall


On Wed, Jul 27, 2016 at 9:06 AM Slavin, Jonathan <jslavin at cfa.harvard.edu>
wrote:

> Hi Josh,
>
> That's an interesting approach and one that's nicely containable within
> Simulation_initBlock.  Would you be willing to share that code?  I was
> wrong, by the way, about the on/off switch regarding pressure.  In fact
> what is done is a linear interpolation (see around line 250 in
> Simulation_initiBlock.F90 in the Sedov directory).  However, the inadequacy
> of that approximation is shown by the substantial discrepancy between the
> actual energy and sim_expEnergy.
>
> An alternative approach that I thought of that could be applied to the
> thermal energy only case is to, after all blocks have been initialized, 1)
> total up the thermal energy, 2) calculate the ratio  sim_expEnergy/E_total
> and 3) adjust the pressure in all the zones in which the pressure is above
> the ambient pressure by multiplying by that factor.
>
> While I think I could write the code to calculate the ratio and apply it,
> I don't fully understand the structure of FLASH, so I'm not sure where it
> would need to be called from.  One advantage to this approach is that one
> is pretty much guaranteed that the explosion energy will be what is desired.
>
> Any suggestions on this approach would be appreciated.
>
> Thanks,
> Jon
>
> On Tue, Jul 26, 2016 at 10:04 PM, Joshua Wall <joshua.e.wall at gmail.com>
> wrote:
>
>> Hello Jon,
>>
>>     After I saw your statement:
>>
>> On Tue, Jul 26, 2016 at 11:48 AM Slavin, Jonathan <
>> jslavin at cfa.harvard.edu> wrote:
>>
>>>
>>>
>>> One thing that occurred to me is that Simulation_initBlock uses a simple
>>> on/off switch to determine if the pressure in a zone should be set to the
>>> high explosion pressure or the ambient pressure.  If the grid does not
>>> match up well with the sim_rInit, that is, if the edge of a zone does not
>>> align with sim_rInit, then part of a zone that should contain a pressure
>>> above ambient could be left out of the explosion.  I'm looking over the
>>> code now to see if that could be the case.
>>>
>>>
>>
>>     I thought I would mention that one way I've used to inject energy
>> (thermal and kinetic) into the Cartesian grid in spherical regions at the
>> edge of the sphere was to do a simple Monte Carlo integration of the
>> overlapping volume by just sampling the cell over a uniform mesh to see if
>> its inside or outside the sphere. Its one way to quickly get the overlap
>> fraction, which I've used to scale thermal energy or momentum I've put down
>> on the grid. Maybe you try something similar to put the pressure or thermal
>> energy on the grid to "smooth" it out at the edges?
>>
>> Cordially,
>>
>> Josh
>>
>>
>>
>>
>>
>>>
>>> On Tue, Jul 26, 2016 at 11:00 AM, Robert Fisher <rfisher1 at umassd.edu>
>>> wrote:
>>>
>>>> Hi Jonathan :
>>>>
>>>>   I've recently been looking into the Sedov problem as a test case for
>>>> some remnant work I am doing. The Sedov problem is easy to set up as a
>>>> point blast and get "sort of" correct, but because of its self-similar
>>>> nature, it is singular at [image: t = 0], and consequently it is
>>>> extremely challenging to get precisely correct in all its details.
>>>>
>>>>   By default, FLASH initializes a point-like energy deposition in the
>>>> central zones, in an approximation to a delta function. Firstly, if you
>>>> carefully read through the code, you will find that the explosion energy
>>>> parameter is not actually used in this case at all -- the code is simply
>>>> initializing density and pressure. Additionally, adding more resolution is
>>>> a somewhat better approximation to the delta function, but is still not
>>>> quite what you may need if you are looking to do a very accurate
>>>> initialization. A better procedure is to set the initial time (tinitial
>>>> in flash.par) to a small but non-zero value. This will kick in a
>>>> second piece of the initialization code in Simulation_initBlock.F90, which
>>>> initializes the exact solution at the specified time, and actually *does
>>>> use the explosion energy parameter sim_expEnergy*, though this code is
>>>> currently hard-coded to work only for [image: \gamma = 7/5] and for
>>>> 3D. This initialization procedure has the advantage of avoiding the initial
>>>> singularity of the solution. In this case I believe you will find better
>>>> convergence behavior with added resolution, but there are limitations in
>>>> the hardcoded solution.
>>>>
>>>>   Lastly, depending on what sort of accuracy you may be demanding,
>>>> there are also some errors in the well-known texts on this problem which
>>>> will reveal themselves as you test against the exact solutions. There is
>>>> much more on these issues on Frank Timmes website :
>>>> http://cococubed.asu.edu/research_pages/sedov.shtml, where you can
>>>> also find a highly-accurate initialization method which can be used as an
>>>> initialization.
>>>>
>>>>   Best wishes,
>>>>
>>>>   Bob
>>>>
>>>> On Tue, Jul 26, 2016 at 10:48 AM, Marissa Adams <
>>>> madams at pas.rochester.edu> wrote:
>>>>
>>> Hi Jonathan,
>>>>>
>>>>> I have also noticed this problem.
>>>>>
>>>>> The energy one sees in the .dat file is typically less than the input
>>>>> energy for sim_expEnergy. Also noticed that this is geometry dependent.
>>>>>
>>>>> I input on the order of 1e10 ergs, and in cartesian I get out 9.2e9;
>>>>> for cylindrical I get even less! 3.2e8... which is huge. Note I am going
>>>>> based off of memory here. I cannot get into my HPC system now.
>>>>>
>>>>> If I am remembering things correctly: seems that in the source code's
>>>>> volume function in Simulation_init.f90 doesn't always choose the volume
>>>>> you're choosing to distribute the energy over.
>>>>>
>>>>> I think it is the section where you'll find a comment that notes the
>>>>> code will calculate the initial volume and interior pressure. It is
>>>>> assuming pi r^2 for to get the energy density... where as you may want
>>>>> something different for other chosen geometries.
>>>>>
>>>>> If I remember correctly the output energy also changes depending where
>>>>> you put the hotspot.
>>>>>
>>>>> Sorry for the haphazard e-mail, but I am experiencing similar issues!
>>>>>
>>>>> Best,
>>>>> Marissa
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On Tue, Jul 26, 2016 at 10:30 AM, Slavin, Jonathan <
>>>>> jslavin at cfa.harvard.edu> wrote:
>>>>>
>>>>>> Hi all,
>>>>>>
>>>>>> I've been running supernova remnant evolution simulations.  The
>>>>>> initialization is done with code taken from the supplied Sedov simulation
>>>>>> with certain small changes.  One issue that I've run into is that the total
>>>>>> energy is not very accurate.  That is, though the energy is conserved, the
>>>>>> total energy, initially input as thermal energy, is below that which is
>>>>>> given in the sim_expEnergy parameter by ~5%.  Has anyone else run into
>>>>>> this?  I don't think that it has to do with any of the changes that I've
>>>>>> made to Simulation_initBlock or Simulation_init.  I'm wondering if
>>>>>> increasing the sim_nsubzones parameter might help.  Any advice would be
>>>>>> appreciated.
>>>>>>
>>>>>> Thanks,
>>>>>> Jon
>>>>>>
>>>>>> --
>>>>>> ________________________________________________________
>>>>>> Jonathan D. Slavin                 Harvard-Smithsonian CfA
>>>>>> jslavin at cfa.harvard.edu       60 Garden Street, MS 83
>>>>>> phone: (617) 496-7981       Cambridge, MA 02138-1516
>>>>>> cell: (781) 363-0035             USA
>>>>>> ________________________________________________________
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> Marissa Adams                         E-mail: madams at pas.rochester.edu
>>>>> Graduate Student                      Ph: (585) 402-5779
>>>>> Department of Physics & Astronomy     Website:
>>>>> http://www.pas.rochester.edu/~madams/
>>>>> University of Rochester
>>>>> 478 Bausch & Lomb Hall
>>>>> P.O. Box 270171
>>>>>
>>>>> --
>>>>> Dr. Robert Fisher
>>>>> Associate Professor / Graduate Program Director
>>>>> University of Massachusetts/Dartmouth
>>>>> Department of Physics
>>>>> 285 Old Westport Road
>>>>> North Dartmouth, Massachusetts 02747
>>>>> robert.fisher at umassd.edu
>>>>> http://www.novastella.org
>>>>>
>>>> --
>> Joshua Wall
>> Doctoral Candidate
>> Department of Physics
>> Drexel University
>> 3141 Chestnut Street
>> Philadelphia, PA 19104
>>
>
>
>
> --
> ________________________________________________________
> Jonathan D. Slavin                 Harvard-Smithsonian CfA
> jslavin at cfa.harvard.edu       60 Garden Street, MS 83
> phone: (617) 496-7981       Cambridge, MA 02138-1516
> cell: (781) 363-0035             USA
> ________________________________________________________
>
> --
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/20160728/fb956fa8/attachment.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: overlap.F90
Type: text/x-fortran
Size: 2547 bytes
Desc: not available
URL: <http://flash.rochester.edu/pipermail/flash-users/attachments/20160728/fb956fa8/attachment.bin>


More information about the flash-users mailing list