ENH: avoid unnecessary array copies through np.ndarray.astype where p…

…ossible

yt/data_objects/construction_data_containers.py

@@ -153,7 +153,7 @@ def _get_cut_mask(self, grid):
         for mi, (i, pos) in enumerate(zip(pids, self.positions[points_in_grid])):
             if not points_in_grid[i]:
                 continue
-            ci = ((pos - grid.LeftEdge) / grid.dds).astype("int")
+            ci = ((pos - grid.LeftEdge) / grid.dds).astype("int64")
             if grid.child_mask[ci[0], ci[1], ci[2]] == 0:
                 continue
             for j in range(3):
@@ -1581,7 +1581,7 @@ def _minimal_box(self, dds):
             # How many root cells do we occupy?
             end_index = np.rint(cell_end).astype("int64")
             dims = end_index - start_index + 1
-        return start_index, end_index.astype("int64"), dims.astype("int32")
+        return start_index, end_index, dims.astype("int32")
 
     def _update_level_state(self, level_state):
         ls = level_state

yt/data_objects/index_subobjects/octree_subset.py

@@ -234,7 +234,7 @@ def mesh_sampling_particle_field(self, positions, mesh_field, lvlmax=None):
             positions.shape[0],
             positions.shape[0] ** 0.3333333,
         )
-        pos = positions.to("code_length").value.astype("float64")
+        pos = positions.to_value("code_length").astype("float64", copy=False)
 
         op.process_octree(
             self.oct_handler,

yt/data_objects/level_sets/contour_finder.py

@@ -40,7 +40,7 @@ def identify_contours(data_source, field, min_val, max_val, cached_fields=None):
             g.id, [contour_ids.view("float64")], mask, LE, RE, dims.astype("int64")
         )
         contours[nid] = (g.Level, node.node_ind, pg, sl)
-    node_ids = np.array(node_ids).astype("int64")
+    node_ids = np.array(node_ids, dtype="int64")
     if node_ids.size == 0:
         return 0, {}
     trunk = data_source.tiles.tree.trunk

yt/data_objects/level_sets/tests/test_clump_finding.py

@@ -128,8 +128,8 @@ def test_clump_tree_save():
     t2 = list(ds2.tree)
     mt1 = ds.arr([c.info["cell_mass"][1] for c in t1])
     mt2 = ds2.arr([c["clump", "cell_mass"] for c in t2])
-    it1 = np.array(np.argsort(mt1).astype(int))
-    it2 = np.array(np.argsort(mt2).astype(int))
+    it1 = np.argsort(mt1).astype("int64")
+    it2 = np.argsort(mt2).astype("int64")
     assert_array_equal(mt1[it1], mt2[it2])
 
     for i1, i2 in zip(it1, it2):
@@ -143,8 +143,8 @@ def test_clump_tree_save():
     c2 = list(ds2.leaves)
     mc1 = ds.arr([c.info["cell_mass"][1] for c in c1])
     mc2 = ds2.arr([c["clump", "cell_mass"] for c in c2])
-    ic1 = np.array(np.argsort(mc1).astype(int))
-    ic2 = np.array(np.argsort(mc2).astype(int))
+    ic1 = np.argsort(mc1).astype("int64")
+    ic2 = np.argsort(mc2).astype("int64")
     assert_array_equal(mc1[ic1], mc2[ic2])
 
     os.chdir(curdir)

yt/data_objects/particle_trajectories.py

@@ -46,7 +46,7 @@ class ParticleTrajectories:
     ... ]
     >>> ds = load(my_fns[0])
     >>> init_sphere = ds.sphere(ds.domain_center, (0.5, "unitary"))
-    >>> indices = init_sphere[("all", "particle_index")].astype("int")
+    >>> indices = init_sphere[("all", "particle_index")].astype("int64")
     >>> ts = DatasetSeries(my_fns)
     >>> trajs = ts.particle_trajectories(indices, fields=fields)
     >>> for t in trajs:
@@ -284,8 +284,8 @@ def _get_data(self, fields):
                             pfield[field],
                             self.num_indices,
                             cube[fds[field]],
-                            np.array(grid.LeftEdge).astype(np.float64),
-                            np.array(grid.ActiveDimensions).astype(np.int32),
+                            np.array(grid.LeftEdge, dtype="float64"),
+                            np.array(grid.ActiveDimensions, dtype="int32"),
                             grid.dds[0],
                         )
             sto.result_id = ds.parameter_filename

yt/data_objects/time_series.py

@@ -417,7 +417,7 @@ def particle_trajectories(
         ... ]
         >>> ds = load(my_fns[0])
         >>> init_sphere = ds.sphere(ds.domain_center, (0.5, "unitary"))
-        >>> indices = init_sphere[("all", "particle_index")].astype("int")
+        >>> indices = init_sphere[("all", "particle_index")].astype("int64")
         >>> ts = DatasetSeries(my_fns)
         >>> trajs = ts.particle_trajectories(indices, fields=fields)
         >>> for t in trajs:

yt/fields/fluid_vector_fields.py

@@ -36,21 +36,21 @@ def setup_fluid_vector_fields(
     sl_center = slice(1, -1, None)
 
     def _baroclinic_vorticity_x(field, data):
-        rho2 = data[ftype, "density"].astype(np.float64) ** 2
+        rho2 = data[ftype, "density"].astype("float64", copy=False) ** 2
         return (
             data[ftype, "pressure_gradient_y"] * data[ftype, "density_gradient_z"]
             - data[ftype, "pressure_gradient_z"] * data[ftype, "density_gradient_z"]
         ) / rho2
 
     def _baroclinic_vorticity_y(field, data):
-        rho2 = data[ftype, "density"].astype(np.float64) ** 2
+        rho2 = data[ftype, "density"].astype("float64", copy=False) ** 2
         return (
             data[ftype, "pressure_gradient_z"] * data[ftype, "density_gradient_x"]
             - data[ftype, "pressure_gradient_x"] * data[ftype, "density_gradient_z"]
         ) / rho2
 
     def _baroclinic_vorticity_z(field, data):
-        rho2 = data[ftype, "density"].astype(np.float64) ** 2
+        rho2 = data[ftype, "density"].astype("float64", copy=False) ** 2
         return (
             data[ftype, "pressure_gradient_x"] * data[ftype, "density_gradient_y"]
             - data[ftype, "pressure_gradient_y"] * data[ftype, "density_gradient_x"]
@@ -260,23 +260,23 @@ def _vorticity_growth_timescale(field, data):
     ########################################################################
 
     def _vorticity_radiation_pressure_x(field, data):
-        rho = data[ftype, "density"].astype(np.float64)
+        rho = data[ftype, "density"].astype("float64", copy=False)
         return (
             data[ftype, "radiation_acceleration_y"] * data[ftype, "density_gradient_z"]
             - data[ftype, "radiation_acceleration_z"]
             * data[ftype, "density_gradient_y"]
         ) / rho
 
     def _vorticity_radiation_pressure_y(field, data):
-        rho = data[ftype, "density"].astype(np.float64)
+        rho = data[ftype, "density"].astype("float64", copy=False)
         return (
             data[ftype, "radiation_acceleration_z"] * data[ftype, "density_gradient_x"]
             - data[ftype, "radiation_acceleration_x"]
             * data[ftype, "density_gradient_z"]
         ) / rho
 
     def _vorticity_radiation_pressure_z(field, data):
-        rho = data[ftype, "density"].astype(np.float64)
+        rho = data[ftype, "density"].astype("float64", copy=False)
         return (
             data[ftype, "radiation_acceleration_x"] * data[ftype, "density_gradient_y"]
             - data[ftype, "radiation_acceleration_y"]

yt/frontends/adaptahop/data_structures.py

@@ -309,7 +309,7 @@ def ihalo(self):
             return ihalo
         else:
             halo_id = self.particle_identifier
-            halo_ids = self.halo_ds.r["halos", "particle_identifier"].astype(int)
+            halo_ids = self.halo_ds.r["halos", "particle_identifier"].astype("int64")
             ihalo = np.searchsorted(halo_ids, halo_id)
 
             assert halo_ids[ihalo] == halo_id
@@ -339,7 +339,7 @@ def _set_halo_member_data(self):
             f *= 1.1
             sph = parent_ds.sphere(center, f * radius)
 
-            part_ids = sph[ptype, "particle_identity"].astype(int)
+            part_ids = sph[ptype, "particle_identity"].astype("int64")
 
             ok = len(np.lib.arraysetops.setdiff1d(members, part_ids)) == 0
 
@@ -348,7 +348,7 @@ def _set_halo_member_data(self):
 
         # Build subregion that only contains halo particles
         reg = sph.cut_region(
-            ['np.in1d(obj[("io", "particle_identity")].astype(int), members)'],
+            ['np.in1d(obj[("io", "particle_identity")].astype("int64"), members)'],
             locals={"members": members, "np": np},
         )
 

yt/frontends/adaptahop/tests/test_outputs.py

@@ -68,10 +68,10 @@ def test_get_halo(brick):
     members = np.sort(halo.member_ids)
 
     # Check sphere contains all the members
-    id_sphere = halo.sphere["io", "particle_identity"].astype(int)
+    id_sphere = halo.sphere["io", "particle_identity"].astype("int64")
     assert len(np.lib.arraysetops.setdiff1d(members, id_sphere)) == 0
 
     # Check region contains *only* halo particles
-    id_reg = np.sort(halo["io", "particle_identity"].astype(int))
+    id_reg = np.sort(halo["io", "particle_identity"].astype("int64"))
 
     np.testing.assert_equal(members, id_reg)

yt/frontends/arepo/io.py

@@ -16,7 +16,7 @@ def _get_smoothing_length(self, data_file, position_dtype, position_shape):
         ind = int(ptype[-1])
         si, ei = data_file.start, data_file.end
         with h5py.File(data_file.filename, mode="r") as f:
-            pcount = f["/Header"].attrs["NumPart_ThisFile"][ind].astype("int")
+            pcount = f["/Header"].attrs["NumPart_ThisFile"][ind].astype("int64")
             pcount = np.clip(pcount - si, 0, ei - si)
             # Arepo cells do not have "smoothing lengths" by definition, so
             # we compute one here by finding the radius of the sphere

yt/frontends/athena/data_structures.py

@@ -99,11 +99,11 @@ def parse_line(line, grid):
         grid["level"] = int(str23(splitup[time_index + 3]).rstrip(","))
         grid["domain"] = int(str23(splitup[time_index + 5]).rstrip(","))
     elif chk23("DIMENSIONS") in splitup:
-        grid["dimensions"] = np.array(str23(splitup[-3:])).astype("int")
+        grid["dimensions"] = np.array(str23(splitup[-3:]), dtype="int")
     elif chk23("ORIGIN") in splitup:
-        grid["left_edge"] = np.array(str23(splitup[-3:])).astype("float64")
+        grid["left_edge"] = np.array(str23(splitup[-3:]), dtype="float64")
     elif chk23("SPACING") in splitup:
-        grid["dds"] = np.array(str23(splitup[-3:])).astype("float64")
+        grid["dds"] = np.array(str23(splitup[-3:]), dtype="float64")
     elif chk23("CELL_DATA") in splitup or chk23("POINT_DATA") in splitup:
         grid["ncells"] = int(str23(splitup[-1]))
     elif chk23("SCALARS") in splitup:
@@ -147,7 +147,7 @@ def _detect_output_fields(self):
             chkp = chk23("POINT_DATA")
             if chkd in splitup:
                 field = str23(splitup[-3:])
-                grid_dims = np.array(field).astype("int")
+                grid_dims = np.array(field, dtype="int64")
                 line = check_readline(f)
             elif chkc in splitup or chkp in splitup:
                 grid_ncells = int(str23(splitup[-1]))
@@ -317,7 +317,7 @@ def _parse_index(self):
         # know the extent of all the grids.
         glis = np.round(
             (glis - self.dataset.domain_left_edge.ndarray_view()) / gdds
-        ).astype("int")
+        ).astype("int64")
         new_dre = np.max(gres, axis=0)
         dre_units = self.dataset.domain_right_edge.uq
         self.dataset.domain_right_edge = np.round(new_dre, decimals=12) * dre_units
@@ -329,7 +329,7 @@ def _parse_index(self):
         )
         self.dataset.domain_dimensions = np.round(
             self.dataset.domain_width / gdds[0]
-        ).astype("int")
+        ).astype("int64")
 
         if self.dataset.dimensionality <= 2:
             self.dataset.domain_dimensions[2] = 1
@@ -382,7 +382,7 @@ def _parse_index(self):
             gdds = (self.grid_right_edge - self.grid_left_edge) / self.grid_dimensions
             glis = np.round(
                 (self.grid_left_edge - self.ds.domain_left_edge) / gdds
-            ).astype("int")
+            ).astype("int64")
             for i in range(self.num_grids):
                 self.grids[i] = self.grid(
                     i,

yt/frontends/athena_pp/data_structures.py

@@ -74,7 +74,7 @@ def _initialize_mesh(self):
         nlevel = self._handle.attrs["MaxLevel"] + 1
 
         nb = np.array([nbx, nby, nbz], dtype="int64")
-        self.mesh_factors = np.ones(3, dtype="int64") * ((nb > 1).astype("int") + 1)
+        self.mesh_factors = np.ones(3, dtype="int64") * ((nb > 1) + 1)
 
         block_grid = -np.ones((nbx, nby, nbz, nlevel), dtype="int64")
         block_grid[log_loc[:, 0], log_loc[:, 1], log_loc[:, 2], levels[:]] = np.arange(

yt/frontends/fits/data_structures.py

@@ -485,7 +485,7 @@ def _determine_nprocs(self):
         if self.specified_parameters["nprocs"] is None:
             nprocs = np.around(
                 np.prod(self.domain_dimensions) / 32**self.dimensionality
-            ).astype("int")
+            ).astype("int64")
             self.parameters["nprocs"] = max(min(nprocs, 512), 1)
         else:
             self.parameters["nprocs"] = self.specified_parameters["nprocs"]
@@ -724,7 +724,7 @@ def _domain_decomp(self):
         dz = self.ds.quan(1.0, "code_length") * self.ds.spectral_factor
         self.grid_dimensions[:, 2] = np.around(
             float(self.ds.domain_dimensions[2]) / self.num_grids
-        ).astype("int")
+        ).astype("int64")
         self.grid_dimensions[-1, 2] += self.ds.domain_dimensions[2] % self.num_grids
         self.grid_left_edge[0, 2] = self.ds.domain_left_edge[2]
         self.grid_left_edge[1:, 2] = (
@@ -811,7 +811,7 @@ def _parse_parameter_file(self):
     def _determine_nprocs(self):
         # If nprocs is None, do some automatic decomposition of the domain
         if self.specified_parameters["nprocs"] is None:
-            nprocs = np.around(self.domain_dimensions[2] / 8).astype("int")
+            nprocs = np.around(self.domain_dimensions[2] / 8).astype("int64")
             self.parameters["nprocs"] = max(min(nprocs, 512), 1)
         else:
             self.parameters["nprocs"] = self.specified_parameters["nprocs"]

yt/frontends/fits/misc.py

@@ -186,8 +186,8 @@ def ds9_region(ds, reg, obj=None, field_parameters=None):
         prefix = ""
 
     def _reg_field(field, data):
-        i = data[prefix + "xyz"[ds.lon_axis]].d.astype("int") - 1
-        j = data[prefix + "xyz"[ds.lat_axis]].d.astype("int") - 1
+        i = data[prefix + "xyz"[ds.lon_axis]].d.astype("int64") - 1
+        j = data[prefix + "xyz"[ds.lat_axis]].d.astype("int64") - 1
         new_mask = mask[i, j]
         ret = np.zeros(data[prefix + "x"].shape)
         ret[new_mask] = 1.0

yt/frontends/gadget/io.py

@@ -70,7 +70,7 @@ def _read_particle_coords(self, chunks, ptf):
     def _yield_coordinates(self, data_file, needed_ptype=None):
         si, ei = data_file.start, data_file.end
         f = h5py.File(data_file.filename, mode="r")
-        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int")
+        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int64")
         np.clip(pcount - si, 0, ei - si, out=pcount)
         pcount = pcount.sum()
         for key in f.keys():
@@ -236,7 +236,7 @@ def _read_particle_data_file(self, data_file, ptf, selector=None):
     def _count_particles(self, data_file):
         si, ei = data_file.start, data_file.end
         f = h5py.File(data_file.filename, mode="r")
-        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int")
+        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int64")
         f.close()
         if None not in (si, ei):
             np.clip(pcount - si, 0, ei - si, out=pcount)

yt/frontends/ramses/fields.py

@@ -320,7 +320,7 @@ def _species_mass(field, data, species: str):
         def gen_pdens(igroup):
             def _photon_density(field, data):
                 # The photon density depends on the possibly level-dependent conversion factor.
-                ilvl = data["index", "grid_level"].astype(int)
+                ilvl = data["index", "grid_level"].astype("int64")
                 dc = dens_conv[ilvl]
                 rv = data["ramses-rt", f"Photon_density_{igroup + 1}"] * dc
                 return rv

yt/frontends/swift/io.py

@@ -40,7 +40,7 @@ def _read_particle_coords(self, chunks, ptf):
     def _yield_coordinates(self, sub_file, needed_ptype=None):
         si, ei = sub_file.start, sub_file.end
         f = h5py.File(sub_file.filename, mode="r")
-        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int")
+        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int64")
         np.clip(pcount - si, 0, ei - si, out=pcount)
         pcount = pcount.sum()
         for key in f.keys():
@@ -63,7 +63,7 @@ def _get_smoothing_length(self, sub_file, pdtype=None, pshape=None):
         ind = int(ptype[-1])
         si, ei = sub_file.start, sub_file.end
         with h5py.File(sub_file.filename, mode="r") as f:
-            pcount = f["/Header"].attrs["NumPart_ThisFile"][ind].astype("int")
+            pcount = f["/Header"].attrs["NumPart_ThisFile"][ind].astype("int64")
             pcount = np.clip(pcount - si, 0, ei - si)
             # we upscale to float64
             hsml = f[ptype]["SmoothingLength"][si:ei, ...]
@@ -116,7 +116,7 @@ def _read_particle_data_file(self, sub_file, ptf, selector=None):
     def _count_particles(self, data_file):
         si, ei = data_file.start, data_file.end
         f = h5py.File(data_file.filename, mode="r")
-        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int")
+        pcount = f["/Header"].attrs["NumPart_ThisFile"][:].astype("int64")
         f.close()
         # if this data_file was a sub_file, then we just extract the region
         # defined by the subfile

yt/frontends/ytdata/data_structures.py

@@ -514,7 +514,7 @@ def _parse_parameter_file(self):
                 )
                 self.domain_dimensions = (
                     (self.domain_right_edge - self.domain_left_edge) / dx
-                ).astype(int)
+                ).astype("int64")
             else:
                 self.domain_right_edge = self.parameters["right_edge"]
                 self.domain_dimensions = self.parameters["ActiveDimensions"]

yt/geometry/grid_geometry_handler.py

@@ -417,9 +417,9 @@ def _chunk_io(
             chunk_ngrids = len(gobjs)
             if chunk_ngrids > 0:
                 nproc = int(ytcfg.get("yt", "internals", "global_parallel_size"))
-                chunking_factor = np.ceil(
-                    self._grid_chunksize * nproc / chunk_ngrids
-                ).astype("int")
+                chunking_factor = np.int64(
+                    np.ceil(self._grid_chunksize * nproc / chunk_ngrids)
+                )
                 size = max(self._grid_chunksize // chunking_factor, 1)
             else:
                 size = self._grid_chunksize

yt/geometry/oct_geometry_handler.py

@@ -35,8 +35,8 @@ def _cell_index(field, data):
             # Get the position of the particles
             pos = data[ptype, "particle_position"]
             Npart = pos.shape[0]
-            ret = np.zeros(Npart)
-            tmp = np.zeros(Npart)
+            ret = np.zeros(Npart, dtype="float64")
+            tmp = np.zeros(Npart, dtype="float64")
 
             if isinstance(data, FieldDetector):
                 return ret
@@ -71,7 +71,7 @@ def _cell_index(field, data):
                 remaining[remaining] = np.isnan(tmp[:Nremaining])
                 Nremaining = remaining.sum()
 
-            return data.ds.arr(ret.astype(np.float64), units="1")
+            return data.ds.arr(ret, units="1")
 
         def _mesh_sampling_particle_field(field, data):
             """