[yt-users] fields in ghost zones

Jean-Claude Passy jcpassy at gmail.com
Sun Jun 6 09:05:15 PDT 2010


Hi David, hi Greg,

your guess was totally right. I was computing with 64 bits and writing 
with 32 bits only. I fixed it and restarting Enzo worked well. Great !
Now, I can get back to Greg's suggestion which was using the restarting 
capability of Enzo instead of reading an HDF5 in my routines. That would 
make my code much lighter. However, I need to change a few things and 
would like to have your opinion about it. I want to:
    - change some global parameters like StopCycle, CycleDataDump, 
GlobalDir,... Would it work if I only edit the files 
DDnnnn/CommonEnvelopennnn ? Wouldn't it create a conflict somewhere ?
    - keep all BaryonFields the same except the velocity that I need to 
set up to 0. Is there an easy way to do that ?

Thanks for your help,

JC

David Collins a écrit :
>> The issue is I can't trust any simulation done with restarting Enzo. Look at
>>     
>
> I wouldn't be so sure that the solution you propose isn't going to
> suffer from the same problems that you see with restarting.  The
> restarting machinery doesn't do anything fancy, there aren't many
> places for errors to creep in there.  So baring a few things, this is
> potentially a deeper issue.
>
> One thing to check, which you may have done already, is to ensure that
> you're writing with the same precision that you're computing with.
> Ensure that the top 5 things in make show-config are all 64 (or 32, if
> you must)
>    CONFIG_PRECISION:             64
>    CONFIG_PARTICLES:             64
>    CONFIG_INTEGERS:              64
>    CONFIG_INITS:                 64
>    CONFIG_IO:                    64
> (Inits doesn't matter for you, but the last one, CONFIG_IO, is important)
>
> Have you checked that restarting without changing anything at all
> gives different answers?
> What compiler optimizations are you using?  Try with -O1 or -O0,
> sometimes more aggressive compiler options can change the answer.
>
>
> d.
>
>
>   
>> those two plots: http://drop.io/0etjfmc. I don't want to bother you with too
>> many details but they represent the separation between my two particles -
>> the core of the RGB star and the compact companion.
>>
>> The first run was with dumping every 80 cycles - dumping based on CycleSkip.
>> It stopped at dumping #196. I restarted it and got data #197, #198,...Then,
>> I wanted to check that restarting did not mess things up, so I run the same
>> sim again but with a dumping frequency twice smaller. So data #98 and #99
>> correspond to the same cycle as file #196 and #198 respectively. Files #98
>> and #196 are exactly the same (which is expected) but #99 and #198 are not.
>> They are slightly different and lead at the end to a unphysical behavior -
>> the companion bouncing outwards instead of sinking deeper.
>>
>> I don't what is causing that, may be the new structure of particles I am
>> using. Has anyone already had that kind of problem ?
>> I will inspect the data at some point in the future but in the meantime, I
>> cannot restart Enzo.
>>
>> Cheers,
>>
>> JC
>>
>> Greg Bryan a écrit :
>>
>> Or if you just want to read the data directly back in, why not do a
>> restart ("enzo -r")?
>>
>> Greg
>>
>> On Jun 2, 2010, at 3:36 PM, David Collins wrote:
>>
>>
>>
>> Hi, JC--
>>
>>
>> Would it be easier to skip YT all together, and just have your problem
>> generator read the HDF5 file produced by the 1d run directly?  This
>> would
>> streamline your processes, and you wouldn't have the noise incurred
>> by both the
>> ASCII file format and the two units conversions.
>>
>> I would probably approach your task by looping over only the active
>> zones in enzo, and incrementing an extra counter.
>>
>> int counter = 0;
>> for(int field=0;field<NumberOfBaryonFields;field++)
>> for(int k=GridStartIndex[2];k<=GridEndIndex[2]; k++)
>>  for(int j=GridStartIndex[1];j<=GridEndIndex[1]; j++)
>>    for(int i=GridStartIndex[0];i<=GridEndIndex[0]; i++){
>>      index = i + GridDimension[0]*(j + GridDimension[1]*k);
>>      BaryonField[ DensNum ][ i ] = DensityFromFile[ counter ];
>>      counter++
>>    }
>>
>>
>>
>> This code snippit came from the wiki page on Baryon access, for future
>> reference:
>> http://lca.ucsd.edu/projects/enzo/wiki/Tutorials/BaryonFieldAccess
>>
>> Does that help?
>> d.
>>
>>
>> On Wed, Jun 2, 2010 at 12:07 PM, Jean-Claude Passy
>> <jcpassy at gmail.com> wrote:
>>
>>
>> Hi Matt,
>>
>> I will try to explain what my goal is.
>> I have a 1D profile of a star that is physically at the equilibrium.
>> Regarding the problem initialization, I set up the density and the
>> total
>> specific energy fields according to the model, velocities equal 0.
>> However,
>> the star is not numerically stable so I have to relax it.
>> Therefore, I
>> divide the velocity by 2 at each timestep, and I let Enzo run for a
>> few
>> dynamical times. Then, I read the last data dump and create an
>> relaxed_model.dat:
>>
>> ##########################################################
>> region = pf.h.region([0.5, 0.5, 0.5], [0.0, 0.0, 0.0], [1.0, 1.0,
>> 1.0])
>> output = 'relaxed_model.dat' f = open(output,'w') print >>f, '#
>> Grid Index
>> rho pressure' density = region["Density"] pressure =
>> region["Pressure"] size
>> = NP.size(density)
>> for i in range(0,size,1):             dtmp = density[i]
>> ptmp =
>> pressure[i]             print >>f, "0", i, '%06e' % dtmp, '%06e' %
>> ptmp
>> ##########################################################
>>
>> so I end up with a file like this:
>>
>> ########################################################## #
>> Restart file
>> for only one grid
>> # Grid Index rho pressure 0 0 5.099011e-10 2.568020e+05 0 1
>> 5.099011e-10
>> 2.568020e+05 0 2 5.099011e-10 2.568020e+05
>> ...
>> ##########################################################
>> Finally, I read index, density, pressure from that previous file
>> and do:
>>
>> ##########################################################  for (k
>> = 0; k <
>> GridDimension[2]; k++)       for (j = 0; j < GridDimension[1]; j++)
>> for (i = 0; i < GridDimension[0]; i++) {
>>           index = i + GridDimension[0]*(j + GridDimension[1]*k);
>> density  /= DensityUnits;           pressure /= PressureUnits;
>> BaryonField[0][index] = density;           BaryonField[1][index] =
>> pressure
>> / ((Gamma - 1.0) * density);
>> ##########################################################
>> The problem is GridDimensions contain the ghost zones so the
>> variable index
>> does not match with the actual index read in the file
>> relaxed_model.dat.
>> That is why I wanted to have that file values for the ghost zones
>> as well.
>>
>> Does it make sense ? Do you have any suggestion ?
>>
>> Thanks for your help,
>>
>>
>> JC
>>
>> Matthew Turk a écrit :
>>
>> Hi Jean-Claude,
>>
>> There are a couple aspects to this.  The first is that Enzo doesn't
>> output the ghost zones -- so any ghost zones handled inside yt are
>> generated by yt.  Were Enzo to output the ghost zones, we would
>> probably be able to handle this, but it does not.
>>
>> Derived fields can depend on the generation of ghost zones, but keep
>> in mind that these are ghost zones generated by yt.  These ghost
>> zones
>> are constructed in a similar manner to how Enzo generates them, but
>> there may be minor differences.  You can manually inspect ghost zones
>> on a *grid* by calling retrieve_ghost_zones on that grid.
>>
>> If you could tell us a bit more about your goal, maybe we could help
>> out a bit more?
>>
>> -Matt
>>
>> On Wed, Jun 2, 2010 at 9:32 AM, Jean-Claude Passy
>> <jcpassy at gmail.com> wrote:
>>
>>
>> Dear all,
>>
>> in order to set up my simulation, I need to access the values of
>> certain
>> fields in the ghost zones.
>> If I do something like:
>>
>> pf = load(data) region = pf.h.region(...) x = region[field]
>>
>> x contains the values of field for the physical grid only. Is there
>> a way I
>> can get the same 1D-array but with the ghost zones included as well ?
>>
>> Thanks for your help,
>>
>> Jean-Claude
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>>
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>>
>>
>> --
>> Sent from my Stone Tablet and carried by my Pterodactyl.
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>>     
>
>
>
>   
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20100606/7b389aab/attachment.htm>


More information about the yt-users mailing list