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

Adds Velocity, ShearStress, StressDiagonal, observables #285

Open
wants to merge 37 commits into
base: main
Choose a base branch
from

Conversation

harveydevereux
Copy link
Collaborator

@harveydevereux harveydevereux commented Aug 21, 2024

This adds a Velocity observable capable of averaging over supplied components x, y, z and atom lists. The behaviour is abstracted so also for stress we can obtain shear stress correlations $\langle \sigma_{xy}\sigma_{xy}\rangle+\langle \sigma_{yz}\sigma_{yz}\rangle+\langle \sigma_{zx}\sigma_{zx}\rangle$ with the built in ShearStress === Stress(['xy', 'yz', 'zx'])

  • ci tests now test against post-process vaf, dependent on changes to it post_process uses full correlation, fix filter_atoms #284 .

  • I used the ShearStress observable to get experimental viscosity for LiF which could be an example perhaps along with the VAF.

  • still need to fix the docs/ update the developer guide.

  • perhaps we could have a helper function to create Velocity-Velocity correlations from a structure (say, creating correlation_kwargs for different atom types automatically). As it stands Observable does not attach to a particular atoms object and so does not have this information at construction. Except implicitly via atoms: list[int].

  • Velocity can take a SliceLike or list[int]. So open-ended ranges can be specified i.e. slice(0, None), but also non-monotonic, decreasing or whatever sequences can be specified that cannot be represented as a slice.

bitmap

Copy link
Member

@ElliottKasoar ElliottKasoar left a comment

Choose a reason for hiding this comment

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

Looks like a good start, thanks @harveydevereux!

I suspect I knew a lot of the answers to my questions when we discussed this when you started, but maybe it's helpful coming to it with fresh eyes, so we can make things as clear as possible.

janus_core/helpers/correlator.py Outdated Show resolved Hide resolved
janus_core/helpers/correlator.py Outdated Show resolved Hide resolved
janus_core/helpers/janus_types.py Outdated Show resolved Hide resolved
janus_core/helpers/janus_types.py Outdated Show resolved Hide resolved
janus_core/helpers/correlator.py Outdated Show resolved Hide resolved
@ElliottKasoar ElliottKasoar added the enhancement New/improved feature or request label Oct 17, 2024
@ElliottKasoar ElliottKasoar linked an issue Oct 17, 2024 that may be closed by this pull request
janus_core/helpers/correlator.py Outdated Show resolved Hide resolved
janus_core/helpers/correlator.py Outdated Show resolved Hide resolved
janus_core/helpers/correlator.py Outdated Show resolved Hide resolved
janus_core/helpers/observables.py Outdated Show resolved Hide resolved
janus_core/helpers/observables.py Outdated Show resolved Hide resolved
janus_core/processing/observables.py Outdated Show resolved Hide resolved
list[int]
The indices for each self._components.
"""
return [self._components[c] for c in self.components]
Copy link
Collaborator

Choose a reason for hiding this comment

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

Is this just list(self._components.values())?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

this may lead to duplicate indices. e.g. stress supports both "xy" and "yx" as self._components. self.components are the particular (say stress) components an observable "tracks" hence the need to access those particular values from the dict.

I've changed this so there is no allowed_components property, just a _allowed_components member to track this. Hopefully makes this clearer.

janus_core/processing/observables.py Outdated Show resolved Hide resolved


StressDiagonal = Stress(components=["xx", "yy", "zz"])
ShearStress = Stress(components=["xy", "yz", "zx"])
Copy link
Collaborator

Choose a reason for hiding this comment

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

Don't like that this isn't StressShear to go with others. May also want StressHydrostatic rather than Diagonal?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

those are both good ideas,StressHydrostatic sounds good, unless someone can think of better I'll go with it

janus_core/processing/observables.py Outdated Show resolved Hide resolved
Comment on lines 141 to 145
for component in self.allowed_components.keys() - components.keys():
if component not in self.allowed_components:
component_names = list(self._components.keys())
raise ValueError(
f"'{component}' invalid, must be '{', '.join(component_names)}'"
)
Copy link
Collaborator

@oerc0122 oerc0122 Nov 4, 2024

Choose a reason for hiding this comment

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

Blindly applying this isn't right.

You no longer need to check the if, as that's done by the operation. component_names is now just undefined and need to refer to the components which are not defined.

if any(self.allowed_components.keys() - components.keys()):
    raise ValueError(
                    f"'{','.join(self.allowed_components.keys() - components.keys())}' invalid, must be in '{', '.join(self.allowed_components)}'"
    )

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

yes, unfortunately I did check the changes locally but then added them all together anyway... the code was odd due to the line length but this can be solved

Copy link
Collaborator

Choose a reason for hiding this comment

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

It's my problem really, I click "suggestion" because I'm too lazy to match the formatting and then it renders with "add this" even when it's invalid code. I should change the type to Python

janus_core/helpers/janus_types.py Outdated Show resolved Hide resolved
janus_core/helpers/utils.py Show resolved Hide resolved
if any(components - self._allowed_components.keys()):
raise ValueError(
f"'{components-self.allowed_components.keys()}'"
" invalid, must be '{', '.join(self._components)}'"
Copy link
Collaborator

Choose a reason for hiding this comment

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

Suggested change
" invalid, must be '{', '.join(self._components)}'"
f" invalid, must be '{', '.join(self._components)}'"

Comment on lines 212 to 217
return (
atoms.get_stress(include_ideal_gas=self.include_ideal_gas, voigt=True)[
self._index
]
sliced_atoms.get_stress(
include_ideal_gas=self.include_ideal_gas, voigt=True
)
/ units.GPa
)
)[self._indices]
Copy link
Collaborator

Choose a reason for hiding this comment

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

I think this is a bit messy and confusing, personally would prefer something more like:

        stresses = sliced_atoms.get_stress(
                include_ideal_gas=self.include_ideal_gas, voigt=True
            ) / units.GPa
        return stresses[self._indices]

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I've split the index onto the return line as above, it is clearer.

Comment on lines 290 to 292
if isinstance(self.atoms_slice, list):
return len(self.atoms_slice) * self.dimension
return slicelike_len_for(self.atoms_slice, self.n_atoms) * self.dimension
Copy link
Collaborator

@oerc0122 oerc0122 Nov 5, 2024

Choose a reason for hiding this comment

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

It might be best to move the list case into slicelike_len_for so we can have a clean interface here. Maybe selector_len or something in that case?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

done, but this also made me realise slicelike_len_for is not defined due to the import guarding (ultimately a circular import issue from before). @ElliottKasoar I have moved CorrelationKwargs into correlator.py so there is no longer an issue with Observable and SliceLike circular imports.

docs/source/conf.py Outdated Show resolved Hide resolved
tests/test_utils.py Outdated Show resolved Hide resolved
if isinstance(index, (slice, range)):
return (index.start, index.stop, index.step)

return index
Copy link
Member

Choose a reason for hiding this comment

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

Is it worth checking if index is a StartStopStep (or any other valid input)?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

do you mean throwing an error if it is not actually a SliceLike?

i.e. these cases are no good but "run"

>>> slicelike_to_startstopstep([1,2])
[1, 2]
>>> slicelike_to_startstopstep((None, None, None))
(None, None, None)

its probably a good move

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

We could do

def validate_slicelike(maybe_slicelike: SliceLike):
    # ...
    if isinstance(maybe_slicelike, (slice, range, int)):
        return
    if isinstance(maybe_slicelike, tuple) and len(maybe_slicelike) == 3:
        start, stop, step = maybe_slicelike
        if (
            isinstance(start, Optional[int])
            and isinstance(stop, Optional[int])
            and isinstance(step, int)
        ):
            return

    raise ValueError(f"{maybe_slicelike} is not a valid SliceLike")

?

Copy link
Member

Choose a reason for hiding this comment

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

Yes, that sort of thing would be great

docs/source/developer_guide/tutorial.rst Outdated Show resolved Hide resolved
docs/source/developer_guide/tutorial.rst Outdated Show resolved Hide resolved
docs/source/developer_guide/tutorial.rst Outdated Show resolved Hide resolved
@ElliottKasoar
Copy link
Member

ElliottKasoar commented Nov 5, 2024

Not really work for this PR as such, but how easy would it be to integrate/use these in post processing, rather than on-the-fly calculations?

@harveydevereux
Copy link
Collaborator Author

Not really work for this PR as such, but how easy would it be to integrate/use these in post processing, rather than on-the-fly calculations?

All it needs is to be called on some sequence of atoms objects. But if you are post-processing why not just use post-process? The only advantage of these is they extract the data from Atoms in this case

@ElliottKasoar
Copy link
Member

All it needs is to be called on some sequence of atoms objects. But if you are post-processing why not just use post-process? The only advantage of these is they extract the data from Atoms in this case

Most significantly, I'm thinking about how easy it is to extend postprocess to observables we've defined, but aren't supported yet.

At the moment we can only do RDFs and VAFs, so if we didn't calculate any other observable on the fly, we might still want to be able to use an Observable.

It could also fit into unifying the interfaces e.g. for the VAF, depending on how their efficiencies compare, which would avoid questions about their consistency any time we make a change.

@harveydevereux harveydevereux force-pushed the add_velocity_observable branch 2 times, most recently from 5e25e15 to 47cf3bc Compare November 8, 2024 11:25
Copy link
Collaborator

@oerc0122 oerc0122 left a comment

Choose a reason for hiding this comment

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

Just so we know it's being run for side-effects.

janus_core/helpers/utils.py Outdated Show resolved Hide resolved
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement New/improved feature or request
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Add velocity observables
3 participants