Skip to content

Conversation

@AnatoleStorck
Copy link

PR Summary

Currently, only scalar quantites of particles can be deposited to the grid (such as mass), with final shape (Ncells,) I want to be able to deposit a vector quantity (stellar spectra of each star particle) onto the grid with shape (Ncells, spectra_length), with the main goal to build spectra for each cell (nebular continuum + emission lines, and stellar spectra.) I have not been able to make the deposition of the vector fields work yet, so any help is greatly appreciated.

PR Checklist

  • New features are documented, with docstrings and narrative docs
  • Adds a test for any bugs fixed. Adds tests for new features.

@AnatoleStorck
Copy link
Author

Just to add, I am currently only considering method="sum" for deposition of vector fields

@chrishavlin
Copy link
Contributor

@cphyc some octree enhancements/changes here that you might be interested in checking out when you're back from holiday :)

@jzuhone
Copy link
Contributor

jzuhone commented Aug 21, 2025

@AnatoleStorck would be happy to look more at this after the test failures are addressed.

@chummels
Copy link
Member

I'd love to see this functionality in yt. Let me know if there are ways I can assist in testing this or debugging.

return self.oct_handler.domain_ind(self.selector)

def deposit(self, positions, fields=None, method=None, kernel_name="cubic"):
def deposit(self, positions, fields=None, method=None, kernel_name="cubic", vector_field=False):
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You're missing the corresponding item in the docstring to explain what vector_field means.


source = source.reshape(np.shape(source) + (dims,))
# If source is a vector field, we need to handle it differently.
is_vector_field = len(np.shape(source)) > 5
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe use source.ndim == 5 instead?

source = source.reshape((source.shape[0], source.shape[1],
source.shape[2], source.shape[3], dims))
if is_vector_field:
dest = np.zeros((len(dest), np.shape(source)[-2], dims),
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am not sure this is correct - if dest is not None, we should modify it in place instead of creating a new array. It is possible that the upstream code relies on dest being modified (if passed).

In addition, np.shape(source)[-2] is cryptic and requires knowing about https://github.com/yt-project/yt/pull/5227/files#diff-3cef3e330a489f3d2975d39b0f17b55a03f3134b27492604af30bc23f1bf7fa7R179. What about something like:

*_, ncells, ndim = source.shape
# assert ndim == dims?
dest = np.zeros((len(dest), ncells, dims), dtype=source.dtype, order='C')

Comment on lines 84 to 85
#print("IN COPY2DARRAYF64")
#print(self.dest.shape, self.source.shape)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
#print("IN COPY2DARRAYF64")
#print(self.dest.shape, self.source.shape)

# even though in other instances it varies the fastest. To see this in
# action, try looking at the results of an n_ref=256 particle CIC plot,
# which shows it the most clearly.
return ((i*dims[1])+j)*dims[2]+k
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Shouldn't this be (i*dims[0] + j) * dims[1] + k?

arr2.shape = arr2.shape + (1,) * (naxes - arr2.ndim)
return arr2

cdef class ParticleDepositOperation:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could it superclass the one defined in particle_deposit.pyx to minimize code duplication?

@cphyc
Copy link
Member

cphyc commented Nov 10, 2025

Note for @AnatoleStorck and @cphyc, the particle_deposit_vector.pyx and .pxd should be removed, and merged into particle_deposit.pyx. The latter should be modified to accept a vector field in general, with the specific case of a scalar field being a vector field of dimension 1.

AnatoleStorck and others added 8 commits November 11, 2025 11:58
This works, but is still using the 'particle_deposit_vector.pyx' code.
The latter is a duplication from 'particle_deposit.pyx', which
only works for scalar fields but could be replaced with the more
generic vector version.
Assume that we're depositing vector, and treat scalar fields
as a special case of vectors
@cphyc cphyc force-pushed the add_vector_field_support_to_deposit branch from 34fba04 to 3526408 Compare November 11, 2025 11:03
@cphyc
Copy link
Member

cphyc commented Nov 11, 2025

@chummels this should support SPH datasets, but I have never used an SPH dataset, and I'm not sure what the particle deposition would mean for SPH? Onto what are you deposing?

Given my ignorance, I haven't really tested this, but it'd be amazing if you could give it a go

@cphyc
Copy link
Member

cphyc commented Nov 11, 2025

We still have failures for AMR (non-octree) datasets. It seems that vector fields end up being stored vector-first, i.e., with shape (vector_size, nx, ny, nz, num_zones) whereas oct-based dataset have shape (nx, ny, nz, num_zones, vector_size). I'm sure there must be a transpose somewhere, but couldn't find it yet. Any help is welcome :)

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 this pull request may close these issues.

5 participants