[yt-users] fields in ghost zones

Jean-Claude Passy jcpassy at gmail.com
Thu Jun 3 04:24:06 PDT 2010


Hi Greg,

that was indeed my first idea.
The issue is I can't trust any simulation done with restarting Enzo. 
Look at 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
>   
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20100603/54d4f088/attachment.html>


More information about the yt-users mailing list