[yt-users] fields in ghost zones
Jean-Claude Passy
jcpassy at gmail.com
Wed Jun 2 14:58:06 PDT 2010
Hi David,
thanks for you help. I think this is a good strategy !
I will let you know how it goes.
Thanks again,
JC
David Collins a écrit :
> Hi, Jean-Claude--
>
> (I copied this to the list, so all may enjoy my wisdom.)
>
> HDF5 brings some extra work over straight ascii files, but it's 1000%
> worth it in the long run. It's your call if you think it's right for
> this problem. I'm the kind of person who will go far out of my way to
> be lazy later, so this is what I'd do. (Even if it's not right for
> this problem, HDF5 is worth knowing. I hate all datasets that aren't
> in HDF5, now that I've gotten comfortable with it.)
>
> The HDF5 website is here:
>
> http://www.hdfgroup.org/HDF5/
>
> I've been told that the HDF5 Lite interface is pretty good:
>
> http://www.hdfgroup.org/HDF5/Tutor/h5lite.html
>
> though I haven't used it myself, the syntax is supposed to be simpler.
> It's worth going through a few examples on the website to get used to
> it.
>
> I would do it using a similar syntax that's in Grid_Group_ReadGrid.C,
> though this file has a bunch of enzo baggage in it, so it's kind of
> difficult to read. One thing to note is that HDF5 is a C library, but
> Enzo arrays are stored in Fortran ordering. So be careful that the
> array will have the X and Z indices swapped if you go to 3d. Another
> thing is to watch out for data types.
>
> The short version: HDF5 files have a directory structure like a unix
> directory tree. Enzo output files have a format like
>
> File/Grid Group/ Dataset
>
> Each processor does it's own output, so there should be one cpu????
> file for each MPI task. So you need to do the following in your
> code:
>
> 1.) Open the file (data0000.cpu0001)
> 2.) Open the "Group", in this case the grid ID, which will be
> something like Grid00001. h5ls will show you the contents of the file
> 3.) Open the dataset ('Density', 'x-velocity', etc)
> 4.) Read the dataset
> 5.) Close the dataset, group, file
>
> Syntax and details can be found on the hdf5 websites or in the Enzo
> Grid_Group_ReadGrid.C files.
>
> Hope that helps!
> d.
>
> On Wed, Jun 2, 2010 at 2:06 PM, Jean-Claude Passy <jcpassy at gmail.com> wrote:
>
>> Hi David,
>>
>> it looks like a good idea, thanks.
>> Could you just be more specific and explain briefly how to read the HDF5
>> file, especially in a C++ program ? I am sorry but I am also not familiar
>> with the h5ls command.
>>
>> Thanks a lot,
>>
>> Jean-Claude
>>
>> David Collins a écrit :
>>
>> 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
>>
>>
>>
>>
>>
>>
>
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20100602/fa94d5aa/attachment.htm>
More information about the yt-users
mailing list