Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Cannot set linear scaling in certain PhasePlots #3431

Closed
chummels opened this issue Jul 12, 2021 · 4 comments · Fixed by #3432
Closed

Cannot set linear scaling in certain PhasePlots #3431

chummels opened this issue Jul 12, 2021 · 4 comments · Fixed by #3432
Assignees
Labels

Comments

@chummels
Copy link
Member

Bug report

Bug summary

There appears to be a bug in PhasePlot in how it calculates extrema on subsequent passes, and it leads to failures when trying to changing PhasePlot parameters like which fields are set to log scaling or linear scaling.

Code for reproduction

This code creates a sphere object, then it attempts to create a PhasePlot of spherical theta vs spherical radius.

import yt
ds = yt.load_sample('FIRE_M12i_ref11')
_, c = ds.find_max(('gas', 'density'))
sp = ds.sphere(c, (150.0, 'kpc'))
p = yt.PhasePlot(sp, ('PartType0', 'particle_position_spherical_theta'), ('PartType0', 'particle_position_spherical_radius'), ('gas', 'mass'), weight_field=None)
p.set_log(('PartType0', 'particle_position_spherical_theta'), False)
p.save()

Actual outcome

The code fails on the line where we try to change the theta field to linear scaling (the second to last line).

Traceback (most recent call last):
  File "test_profile.py", line 7, in <module>
    p.save()
  File "/Users/chummels/src/yt/yt/visualization/plot_container.py", line 102, in newfunc
    args[0]._recreate_profile()
  File "/Users/chummels/src/yt/yt/visualization/profile_plotter.py", line 1611, in _recreate_profile
    self._profile = create_profile(
  File "/Users/chummels/src/yt/yt/data_objects/profiles.py", line 1443, in create_profile
    raise YTIllDefinedBounds(mi, ma)
yt.utilities.exceptions.YTIllDefinedBounds: The bounds -4.941e-324 and 5.233e+23 are ill-defined. Typically this happens when a log binning is specified and zero or negative values are given for the bounds.

Expected outcome

For it to change the theta field to linear scaling and then output a PhasePlot.

Notes

If one removes the second-to-last line: p.set_log(('PartType0', 'particle_position_spherical_theta'), False), the code snippet above works perfectly fine. In digging into this, the issue appears to in the extrema code here. When yt initially generates the ProfilePlot, extrema is equal to None, so it enters the first conditional and correctly calculates that the extrema of the two fields (theta, radius):

extrema = None
ex = unyt_array([0.00500881, 3.1364112 ], '(dimensionless)'), unyt_array([6.65419743e+19, 5.23294092e+23], 'cm')

And everything is fine. But if the last line p.set_log(('PartType0', 'particle_position_spherical_theta'), False) is included in our script, yt goes back through all of the ProfilePlot code and runs it again (I'm not sure why). This time, though, extrema is not None (but it's just composed of all Nones) so it enters the second conditional, which is a bunch of messy code that I cannot really follow. The outcome is that the extrema minima are now negative, and yt fails:

extrema = {('PartType0', 'particle_position_spherical_theta'): (None, None), ('PartType0', 'particle_position_spherical_radius'): (None, None)}
ex = [[array(-5.e-324), array(3.1364112)], [array(-5.e-324), array(5.23294092e+23)]]

I'm submitting a PR that simply checks to see if extrema is composed of all Nones, and if so, it treats it as though it is None, so the same logic occurs. But maybe there is a better way?

Version Information
dev tip:

commit aaaf4d00141b123c089f59978ac81bdfaf8de023 (HEAD -> main, upstream/main)
Merge: 4ecf7bcfb aa733c0f6
Author: Madicken Munk <madicken.munk@gmail.com>
Date:   Fri Jul 9 16:49:36 2021 -0500

    Merge pull request #3427 from yt-project/revert-3425-3269_pplot_zlabel

    Revert "Use a fancier repr for labels in Profile/PhasePlots. Fixes #3269"
@neutrinoceros
Copy link
Member

neutrinoceros commented Jul 13, 2021

Your test script also fails with the current tip of the main branch (7adc079), but with a completely different error.

edit: the only difference on main since you reported this is a pre-commit config update, clearly not related. I don't understand why we're getting different errors.

yt : [INFO     ] 2021-07-13 18:41:19,343 Sample dataset found in '/Users/robcleme/dev/yt-project/test-data-dir/FIRE_M12i_ref11/snapshot_600.hdf5'
yt : [INFO     ] 2021-07-13 18:41:19,394 Calculating time from 1.000e+00 to be 4.355e+17 seconds
yt : [INFO     ] 2021-07-13 18:41:19,395 Assuming length units are in kpc/h (comoving)
yt : [INFO     ] 2021-07-13 18:41:19,468 Parameters: current_time              = 4.3545571088051405e+17 s
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: domain_dimensions         = [1 1 1]
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: domain_left_edge          = [0. 0. 0.]
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: domain_right_edge         = [60000. 60000. 60000.]
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: cosmological_simulation   = 1
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: current_redshift          = 0.0
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: omega_lambda              = 0.728
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: omega_matter              = 0.272
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: omega_radiation           = 0.0
yt : [INFO     ] 2021-07-13 18:41:19,469 Parameters: hubble_constant           = 0.702
yt : [INFO     ] 2021-07-13 18:41:19,483 Allocating for 4.787e+06 particles
Initializing coarse index : 100%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 10/10 [00:00<00:00, 40.88it/s]
╭─────────────────────────────── Traceback (most recent call last) ────────────────────────────────╮
│ /Users/robcleme/dev/yt-project/yt/yt/geometry/particle_geometry_handler.py:181 in                │
│ _initialize_index                                                                                │
│                                                                                                  │
│   178 │   │   try:                                                                               │
│   179 │   │   │   if dont_load:                                                                  │
│   180 │   │   │   │   raise OSError                                                              │
│ ❱ 181 │   │   │   rflag = self.regions.load_bitmasks(fname)                                      │
│   182 │   │   │   rflag = self.regions.check_bitmasks()                                          │
│   183 │   │   │   self._initialize_frontend_specific()                                           │
│   184 │   │   │   if rflag == 0:                                                                 │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/geometry/particle_oct_container.pyx:1028 in                 │
│ yt.geometry.particle_oct_container.ParticleBitmap.load_bitmasks                                  │
│                                                                                                  │
│   1025 │   │   cdef bytes serial_BAC                                                             │
│   1026 │   │   cdef np.uint64_t ifile                                                            │
│   1027 │   │   f = open(fname,'wb')                                                              │
│ ❱ 1028 │   │   # Header                                                                          │
│   1029 │   │   f.write(struct.pack('Q', _bitmask_version))                                       │
│   1030 │   │   f.write(struct.pack('q', self.file_hash))                                         │
│   1031 │   │   f.write(struct.pack('Q', self.nfiles))                                            │
╰──────────────────────────────────────────────────────────────────────────────────────────────────╯
OSError

During handling of the above exception, another exception occurred:

╭─────────────────────────────── Traceback (most recent call last) ────────────────────────────────╮
│ /Users/robcleme/dev/yt-project/yt/t.py:3 in <module>                                             │
│                                                                                                  │
│   1 import yt                                                                                    │
│   2 ds = yt.load_sample('FIRE_M12i_ref11')                                                       │
│ ❱ 3 _, c = ds.find_max(('gas', 'density'))                                                       │
│   4 sp = ds.sphere(c, (150.0, 'kpc'))                                                            │
│   5 p = yt.PhasePlot(sp, ('PartType0', 'particle_position_spherical_theta'), ('PartType0', '     │
│   6 p.set_log(('PartType0', 'particle_position_spherical_theta'), False)                         │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/data_objects/static_output.py:982 in find_max               │
│                                                                                                  │
│    979 │   │   This is a wrapper around _find_extremum                                           │
│    980 │   │   """                                                                               │
│    981 │   │   mylog.debug("Searching for maximum value of %s", field)                           │
│ ❱  982 │   │   return self._find_extremum(field, "max", source=source, to_array=to_array)        │
│    983 │                                                                                         │
│    984 │   def find_min(self, field, source=None, to_array=True):                                │
│    985 │   │   """                                                                               │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/data_objects/static_output.py:947 in _find_extremum         │
│                                                                                                  │
│    944 │   │   """                                                                               │
│    945 │   │   ext = ext.lower()                                                                 │
│    946 │   │   if source is None:                                                                │
│ ❱  947 │   │   │   source = self.all_data()                                                      │
│    948 │   │   method = {                                                                        │
│    949 │   │   │   "min": source.quantities.min_location,                                        │
│    950 │   │   │   "max": source.quantities.max_location,                                        │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/data_objects/static_output.py:1039 in all_data              │
│                                                                                                  │
│   1036 │   │   all_data is a wrapper to the Region object for creating a region                  │
│   1037 │   │   which covers the entire simulation domain.                                        │
│   1038 │   │   """                                                                               │
│ ❱ 1039 │   │   self.index                                                                        │
│   1040 │   │   if find_max:                                                                      │
│   1041 │   │   │   c = self.find_max("density")[1]                                               │
│   1042 │   │   else:                                                                             │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/data_objects/static_output.py:526 in index                  │
│                                                                                                  │
│    523 │   @property                                                                             │
│    524 │   def index(self):                                                                      │
│    525 │   │   if self._instantiated_index is None:                                              │
│ ❱  526 │   │   │   self._instantiated_index = self._index_class(                                 │
│    527 │   │   │   │   self, dataset_type=self.dataset_type                                      │
│    528 │   │   │   )                                                                             │
│    529 │   │   │   # Now we do things that we need an instantiated index for                     │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/geometry/particle_geometry_handler.py:25 in __init__        │
│                                                                                                  │
│    22 │   │   self.dataset = weakref.proxy(ds)                                                   │
│    23 │   │   self.float_type = np.float64                                                       │
│    24 │   │   super().__init__(ds, dataset_type)                                                 │
│ ❱  25 │   │   self._initialize_index()                                                           │
│    26 │                                                                                          │
│    27 │   def _setup_geometry(self):                                                             │
│    28 │   │   self.regions = None                                                                │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/frontends/sph/data_structures.py:91 in _initialize_index    │
│                                                                                                  │
│    88 │   │   if hasattr(self.io, "_generate_smoothing_length"):                                 │
│    89 │   │   │   self.io._generate_smoothing_length(self)                                       │
│    90 │   │                                                                                      │
│ ❱  91 │   │   super()._initialize_index()                                                        │
│    92 │                                                                                          │
│    93 │   def _generate_kdtree(self, fname):                                                     │
│    94 │   │   from yt.utilities.lib.cykdtree import PyKDTree                                     │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/geometry/particle_geometry_handler.py:188 in                │
│ _initialize_index                                                                                │
│                                                                                                  │
│   185 │   │   │   │   raise OSError                                                              │
│   186 │   │   except (OSError, struct.error):                                                    │
│   187 │   │   │   self.regions.reset_bitmasks()                                                  │
│ ❱ 188 │   │   │   self._initialize_coarse_index()                                                │
│   189 │   │   │   self._initialize_refined_index()                                               │
│   190 │   │   │   # We now update fname since index_order2 may have changed                      │
│   191 │   │   │   fname = _current_fname()                                                       │
│                                                                                                  │
│ /Users/robcleme/dev/yt-project/yt/yt/geometry/particle_geometry_handler.py:230 in                │
│ _initialize_coarse_index                                                                         │
│                                                                                                  │
│   227 │   │   │   # domain, it will correspond to a very very high index order, which            │
│   228 │   │   │   # is a large amount of memory!  Having multiple indexes, one for               │
│   229 │   │   │   # each particle type, would fix this.                                          │
│ ❱ 230 │   │   │   new_order2 = self.regions.update_mi2(max_hsml, ds.index_order[1] + 2)          │
│   231 │   │   │   mylog.info(                                                                    │
│   232 │   │   │   │   "Updating index_order2 from %s to %s", ds.index_order[1], new_order2       │
│   233 │   │   │   )                                                                              │
╰──────────────────────────────────────────────────────────────────────────────────────────────────╯
AttributeError: 'yt.geometry.particle_oct_container.ParticleBitmap' object has no attribute 'update_mi2'

@matthewturk
Copy link
Member

matthewturk commented Jul 13, 2021 via email

@chummels
Copy link
Member Author

@neutrinoceros You need to pip install -e . on your local install to rebuild cython methods to make that error go away.

@neutrinoceros
Copy link
Member

Oh that's it, sorry I should have though about it, thank you !

meeseeksmachine pushed a commit to meeseeksmachine/yt that referenced this issue Oct 15, 2021
neutrinoceros added a commit that referenced this issue Oct 16, 2021
…2-on-yt-4.0.x

Backport PR #3432 on branch yt-4.0.x (Checking if extrema are all None in PhasePlot; fixes #3431)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants