[yt-users] cloud in cell mesh construction for particle data

Brendan Griffen brendan.f.griffen at gmail.com
Thu Jun 5 18:34:23 PDT 2014


Thanks very much Nathan. I tried to load in my own data but I think there
are too many particles or I have to specifically set the domain size.

In this area:

data = {'particle_position_x': pos[:,0],
        'particle_position_y': pos[:,1],
        'particle_position_z': pos[:,2],
        'particle_mass': np.array([mpart]*npart)}

ds = yt.load_particles(data, length_unit=3.08e24, mass_unit=1.9891e36)
ad = ds.all_data()
print ad['deposit', 'all_cic']

In [3]: run ytcic.py
yt : [INFO     ] 2014-06-05 21:29:06,183 Parameters: current_time
   = 0.0
yt : [INFO     ] 2014-06-05 21:29:06,183 Parameters: domain_dimensions
    = [2 2 2]
yt : [INFO     ] 2014-06-05 21:29:06,184 Parameters: domain_left_edge
   = [ 0.  0.  0.]
yt : [INFO     ] 2014-06-05 21:29:06,185 Parameters: domain_right_edge
    = [ 1.  1.  1.]
yt : [INFO     ] 2014-06-05 21:29:06,185 Parameters:
cosmological_simulation   = 0.0
yt : [INFO     ] 2014-06-05 21:29:06,188 Allocating for 1.342e+08 particles
---------------------------------------------------------------------------
YTDomainOverflow                          Traceback (most recent call last)
/nfs/blank/h4231/bgriffen/data/lib/yt-x86_64/lib/python2.7/site-packages/IPython/utils/py3compat.pyc
in execfile(fname, *where)
    202             else:
    203                 filename = fname
--> 204             __builtin__.execfile(filename, *where)

/nfs/blank/h4231/bgriffen/work/projects/caterpillar/c2ray/cic/ytcic.py in
<module>()
     52
     53 ad = ds.all_data()
---> 54 print ad['deposit', 'all_cic']
     55 slc = yt.SlicePlot(ds, 2, ('deposit', 'all_cic'))
     56 slc.set_figure_size(4)

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/data_containers.pyc
in __getitem__(self, key)
    205         Returns a single field.  Will add if necessary.
    206         """
--> 207         f = self._determine_fields([key])[0]
    208         if f not in self.field_data and key not in self.field_data:
    209             if f in self._container_fields:

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/data_containers.pyc
in _determine_fields(self, fields)
    453                     raise YTFieldNotParseable(field)
    454                 ftype, fname = field
--> 455                 finfo = self.pf._get_field_info(ftype, fname)
    456             else:
    457                 fname = field

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/static_output.pyc
in _get_field_info(self, ftype, fname)
    445     _last_finfo = None
    446     def _get_field_info(self, ftype, fname = None):
--> 447         self.index
    448         if fname is None:
    449             ftype, fname = "unknown", ftype

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/data_objects/static_output.pyc
in index(self)
    277                 raise RuntimeError("You should not instantiate
Dataset.")
    278             self._instantiated_index = self._index_class(
--> 279                 self, dataset_type=self.dataset_type)
    280             # Now we do things that we need an instantiated index
for
    281             # ...first off, we create our field_info now.

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/frontends/stream/data_structures.pyc
in __init__(self, pf, dataset_type)
    942     def __init__(self, pf, dataset_type = None):
    943         self.stream_handler = pf.stream_handler
--> 944         super(StreamParticleIndex, self).__init__(pf, dataset_type)
    945
    946     def _setup_data_io(self):

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
in __init__(self, pf, dataset_type)
     48         self.directory = os.path.dirname(self.index_filename)
     49         self.float_type = np.float64
---> 50         super(ParticleIndex, self).__init__(pf, dataset_type)
     51
     52     def _setup_geometry(self):

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/geometry_handler.pyc
in __init__(self, pf, dataset_type)
     54
     55         mylog.debug("Setting up domain geometry.")
---> 56         self._setup_geometry()
     57
     58         mylog.debug("Initializing data grid data IO")

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
in _setup_geometry(self)
     52     def _setup_geometry(self):
     53         mylog.debug("Initializing Particle Geometry Handler.")
---> 54         self._initialize_particle_handler()
     55
     56

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
in _initialize_particle_handler(self)
     87                 pf.domain_left_edge, pf.domain_right_edge,
     88                 [N, N, N], len(self.data_files))
---> 89         self._initialize_indices()
     90         self.oct_handler.finalize()
     91         self.max_level = self.oct_handler.max_level

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/geometry/particle_geometry_handler.pyc
in _initialize_indices(self)
    109             npart = sum(data_file.total_particles.values())
    110             morton[ind:ind + npart] = \
--> 111                 self.io._initialize_index(data_file, self.regions)
    112             ind += npart
    113         morton.sort()

/bigbang/data/bgriffen/lib/yt-x86_64/src/yt-hg/yt/frontends/stream/io.pyc
in _initialize_index(self, data_file, regions)
    144                 raise YTDomainOverflow(pos.min(axis=0),
pos.max(axis=0),
    145
 data_file.pf.domain_left_edge,
--> 146
 data_file.pf.domain_right_edge)
    147             regions.add_data_file(pos, data_file.file_id)
    148             morton.append(compute_morton(

YTDomainOverflow: Particle bounds [ 0.  0.  0.] and [ 99.99999237
 99.99999237  99.99999237] exceed domain bounds [ 0.  0.  0.] code_length
and [ 1.  1.  1.] code_length



On Thu, Jun 5, 2014 at 8:22 PM, Nathan Goldbaum <nathan12343 at gmail.com>
wrote:

> Here's a worked out example that does what you're looking for using a fake
> 1 million particle dataset:
>
> http://nbviewer.ipython.org/gist/ngoldbaum/546d37869aafe71cfe38
>
> In this notebook I make use of two key yt features: `load_particles`, and
> `covering_grid`.
>
> load_particles creates a "stream" dataset based on in-memory data fed in
> as a numpy array.  This dataset acts just like an on-disk simulation
> dataset, but doesn't come with the baggage of needing to write a custom
> frontend to read a specific data format off disk.
>
> covering_grid is a way to generate uniform resolution data from an AMR
> dataset. It acts like a python dictionary where the keys are field names
> and returns 3D numpy arrays of whatever uniform resolution you specify when
> you create the covering_grid.
>
> Note that if you're using load_particles all of your data needs to live in
> memory.  If your data is too big for that you'll need to write a frontend
> for your data format or use a memmap to an on-disk file somehow.  I'm not
> an expert on that but others on the list should be able to help out.
>
> Hope that gets you well on your way :)
>
> -Nathan
>
>
> On Thu, Jun 5, 2014 at 5:04 PM, Desika Narayanan <dnarayan at haverford.edu>
> wrote:
>
>> Hey Brendan,
>>
>> A couple of extra tools you might find helpful in conjunction with
>> Nathan's example of depositing the particles onto an octree are at:
>>
>> http://paste.yt-project.org/show/4737/
>>
>> Where I load a gadget snapshot, and then recover the coordinates and
>> width of each cell.
>>
>> In response to your last question - the particles are deposited into an
>> octree grid (so, you'll see that the cell sizes aren't all the same size).
>>   I don't know if depositing onto a regular NxNxN mesh is possible, though
>> would be interested to hear if so.
>>
>> -d
>>
>>
>> On Thu, Jun 5, 2014 at 7:58 PM, Brendan Griffen <
>> brendan.f.griffen at gmail.com> wrote:
>>
>>> Thanks. I'll get the "bleeding edge" version first then try your
>>> suggestions. Though I want to return the NxNxN array and be able to write
>>> this mesh to a file. It is *only* using the cic part of yt and it should
>>> return the mesh to be written? Just wanted to clarify?
>>>
>>> Thanks.
>>> Brendan
>>>
>>>
>>> On Thu, Jun 5, 2014 at 6:49 PM, Nathan Goldbaum <nathan12343 at gmail.com>
>>> wrote:
>>>
>>>>
>>>>
>>>>
>>>> On Thu, Jun 5, 2014 at 3:36 PM, John ZuHone <jzuhone at gmail.com> wrote:
>>>>
>>>>> Hi Brendan,
>>>>>
>>>>> Which version of yt are you using?
>>>>>
>>>>> If you're using 3.0, this is actually fairly easy. If you look in
>>>>> yt.fields.particle_fields.py, around line 85, you can see how this is
>>>>> done for the "particle_density" and "particle_mass" fields. Basically you
>>>>> can call a "deposit" method which takes the particle field quantity you
>>>>> want deposited and deposits it into cells. The underlying calculation is
>>>>> done using Cython, so it's fast.
>>>>>
>>>>
>>>> And you shouldn't ever actually need to call these "deposit" functions,
>>>> since "deposit" is exposed as a field type for all datasets that contain
>>>> particles.
>>>>
>>>> Here is a notebook that does this for Enzo AMR data:
>>>>
>>>> http://nbviewer.ipython.org/gist/ngoldbaum/5e19e4e6cc2bf330149c
>>>>
>>>> This dataset contains about a million particles and generates a CIC
>>>> deposition for the whole domain in about 6 seconds from a cold start.
>>>>
>>>>
>>>>>
>>>>> If you're using 2.x, then you can do the same thing, but it's not as
>>>>> straightforward. You can see how this works in
>>>>> yt.data_objects.universal_fields.py, around line 986, where the
>>>>> "particle_density" field is defined. Basically, it calls CICDeposit_3,
>>>>> which is in yt.utilities.lib.CICDeposit.pyx.
>>>>>
>>>>> Let me know if you need any more clarification.
>>>>>
>>>>> Best,
>>>>>
>>>>> John Z
>>>>>
>>>>> On Jun 5, 2014, at 6:07 PM, Brendan Griffen <
>>>>> brendan.f.griffen at gmail.com> wrote:
>>>>>
>>>>> > Hi,
>>>>> >
>>>>> > I was wondering if there were any Cython routines within yt which
>>>>> takes particle data and converts it into a cloud-in-cell based mesh which
>>>>> can be written to a file of my choosing.
>>>>>
>>>>
>>>> What sort of mesh were you looking for?  yt will internally construct
>>>> an octree if it is fed particle data.  I'm not sure whether this octree can
>>>> be saved to disk for later analysis.
>>>>
>>>>  It's also possible to create a uniform resolution covering grid
>>>> containing field data for a deposited quantity, which can be quite easily
>>>> saved to disk in a number of ways.
>>>>
>>>>
>>>>> I heard a while ago there was some such functionality but it could be
>>>>> too far down the yt rabbit hole to be used as a standalone? Is this true? I
>>>>> have my own Python code for doing it but it just isn't fast enough and
>>>>> thought I'd ask the yt community if there were any wrapper tools available
>>>>> to boost the speed.
>>>>> >
>>>>> > Thanks.
>>>>> > Brendan
>>>>> > _______________________________________________
>>>>> > 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
>>>>
>>>>
>>>
>>> _______________________________________________
>>> 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/20140605/2151a804/attachment.htm>


More information about the yt-users mailing list