mnt: changed index ordering of spin-charge quantities

[zerothi]

Before the resulting charge arrays were not strictly the same order. Now we always have indices such as:

states, spin, orbitals

which means that spin is not the last index anymore.

Furthermore for polarized calculations we will return spin dimension as:

total, Sz

to make it more compatible with the SOC cases.
This will perhaps confuse users in the beginning (those who were used to the old ordering), but in the long run this should be easier to handle since users don't need to put in checks for polarized vs non-colinear calculations.

This closes #501

• Closes Spin-polarized quantities return order #501
• Tests added
• Documentation for functionality
• User visible changes (including notable bug fixes) are documented in CHANGELOG

[codecov]

Codecov Report

Merging #516 (dc0afcf) into main (ecd42bd) will decrease coverage by 0.07%.
The diff coverage is 90.62%.

❗ Current head dc0afcf differs from pull request most recent head 77d73d1. Consider uploading reports for the commit 77d73d1 to get more accurate results

@@            Coverage Diff             @@
##             main     #516      +/-   ##
==========================================
- Coverage   87.16%   87.08%   -0.08%     
==========================================
  Files         353      353              
  Lines       45277    45407     +130     
==========================================
+ Hits        39466    39543      +77     
- Misses       5811     5864      +53     

Impacted Files	Coverage Δ
sisl/geom/surfaces.py	96.64% <ø> (ø)
sisl/io/siesta/pdos.py	32.19% <37.50%> (+0.81%)	⬆️
sisl/io/tbtrans/delta.py	92.73% <89.36%> (-1.11%)	⬇️
sisl/io/tbtrans/tests/test_delta.py	99.17% <96.87%> (-0.83%)	⬇️
sisl/io/siesta/out.py	75.07% <100.00%> (ø)
sisl/io/siesta/tests/test_dm.py	94.11% <100.00%> (ø)
sisl/physics/densitymatrix.py	90.87% <100.00%> (+0.10%)	⬆️
sisl/physics/electron.py	77.17% <100.00%> (ø)
sisl/physics/tests/test_density_matrix.py	100.00% <100.00%> (ø)
sisl/physics/tests/test_energy_density_matrix.py	100.00% <100.00%> (ø)
... and 11 more

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

[pfebrer]

Isn't it a bit weird for PDOS? Because if you calculate it in sisl you need to to the calculation separately per spin channel and therefore you'll get spin 0 and spin 1. Don't you think this will create inconsistencies?

You could create a class for storing PDOS that would retain the information of what each data corresponds to, then you could have a method that converts from one to the other. A PDOSData could be useful for other things as well.

Also, with the current ordering, you have to implement checks for non polarized calculations. For example, every time you want to select:

nonpol = PDOS.ndim == 2

if nonpol:
   selected = PDOS[2, 3]
else:
   selected = PDOS[2, :, 3]

If the spin index was at the end or the beggining you could select in all cases with PDOS[..., 3, 2] or PDOS[2, 3, ...]. If the PDOS was returned as a dataarray you wouldn't have this kind of problems though :) And perhaps you wouldn't need to have a PDOSData structure. You would also be able to merge with other dataarrays that have orbitals, (spin) coordinates to create a dataset, which could be nice. But yeah perhaps it can be confusing for users to work with dataarrays instead of numpy arrays.

[zerothi]

I think you are misreading? The PDOS have ordering spin, orbital, Energy

While I agree that the non-polarized check is necessary (should we add a dimension of 1?) the others are as you suggest? Spin first... Or did I miss something?

[zerothi]

As for the polarized case, I have been wondering about some flag for BZ.apply that also loops spin and merges the quantities... But I am a bit at a loss as to how it should be done without complicating it for end-users... hmm...

[pfebrer]

Oh great, I didn't look at the ordering, I just went by your comment here:

Before the resulting charge arrays were not strictly the same order. Now we always have indices such as:

states, spin, orbitals

So I thought in general you were doing something like (*dimensions, spin, orbitals)

[zerothi]

Oh great, I didn't look at the ordering, I just went by your comment here:

Before the resulting charge arrays were not strictly the same order. Now we always have indices such as:
states, spin, orbitals

So I thought in general you were doing something like (*dimensions, spin, orbitals)

Yes, this is actually the only case where I am not sure exactly what to do. Take a look at spin_moment for instance. I would also be talked into moving the spin dimension first there? @ahkole, @juijan, you have probably used these routines, comments?

[zerothi]

ok, moved it to the start ;)

[pfebrer]

By the way, this will break the pdos plot and spin textured bands in sisl.viz. Let me know before you merge it so that I can make the appropiate changes.

[zerothi]

By the way, this will break the pdos plot and spin textured bands in sisl.viz. Let me know before you merge it so that I can make the appropiate changes.

Ok, thanks for the heads up, I think once @ahkole and @juijan confirms, this should just go in.

[ahkole]

Oh great, I didn't look at the ordering, I just went by your comment here:

Before the resulting charge arrays were not strictly the same order. Now we always have indices such as:
states, spin, orbitals

So I thought in general you were doing something like (*dimensions, spin, orbitals)

Yes, this is actually the only case where I am not sure exactly what to do. Take a look at spin_moment for instance. I would also be talked into moving the spin dimension first there? @ahkole, @juijan, you have probably used these routines, comments?

This one feels more weird to me because in the result value of spin moment the spin dimension refers to the Cartesian component of the spin and I am used to the Cartesian dimension being the last as it is usually the case in Python for quantities such as forces, positions, cell vectors, k points, etc. Although you usually write vectors in linear algebra as columns where the Cartesian dimension would be the first index.

[zerothi]

So would you rather have spin last? Always?

[ahkole]

Now I'm starting to doubt, because when spin refers to spin up or down instead of Cartesian component I like having it as the first index. What do you think @juijan ?

[nils-wittemeier]

I don't have a strong preference for either choice.

One way to look at it is that values returned for [total, x, y, z] are separate quantities, in which case having the spin as the first dimension seems like a natural choice, which would also support q, s_x, s_y, s_z = ...function(). I think it would be consistent to apply this logic to the spin_moment function, the return values are three expectation values of different operators.

If spin is not the last index, I would lean toward adding a dimension of size 1 in the unpolarized case, to always have the same number of dimensions.

At the end of the day, there probably is no perfect solution, and there is also a case to be made for keeping the current order and not breaking existing scripts.

A more consistent solution would always be to return arrays of the shape used in spin-orbit case and using zero as a fill value for spin (un)polarized case. In this way, a user wouldn't have to think about what the spin type of the calculation is. However, this seems like a waste of memory.

[zerothi]

Thanks for your comments!

I can also see your points about different operators and thus always having it the first. This would make things much simpler, and easier to explain in documentation.

So what would you say to that solution?

spin, *rest

[nils-wittemeier]

I'm in favor of spin, *rest.

[zerothi]

I am beginning to think the spin, *rest in the appropriate choice. Any other comments before I commence these changes?

[zerothi]

This will also force the velocity calculations to return cartesian, *rest since they are individual operators. Might be a bit weird, but I think this is the more stringent usage.

[zerothi]

The latest commit should have fixed all remaining routines. I'll give this until next week, then I'll merge. @pfebrer you should get a ping once it is merged, then you can fix the viz stuff?

[zerothi]

It might be a bit weird in the beginning, but I think it might be better to have consistency across the code base. So while it might be a bit unnatural that the velocity has dimensions (cart, states), I think it will be easier in the future.

[zerothi]

@pfebrer I have tried to see if I could understand how you did it in viz, I think I grasp some of it, but also lots I don't ;)

1. I believe that now the PDOS.xml file and PDOS calculation from H both return the same dimensionality of things
2. the TBT.nc handling is however not fixed, in principle the TBT_UP.nc and TBT_DN.nc would be created, but the code does not accommodate reading both files, only 1, could this be enabled to have both in one plot as normally?

[pfebrer]

I wouldn't bother understanding the current framework, as it will work in a very different (hopefully simpler) way :)

Here I would just do the changes so that current functionality is not broken, and any new functionality should be added to the new design. I have not changed anything in the viz module since I created the PR because it is hard to keep track of changes to both designs at the same time.

Otherwise if you think the PR can be merged before the next release the changes can be done there.

Apart from the changes that you did, I think that there should be changes to the handling of PDOS requests, because currently the spin key accepts 0 or 1 in the spin polarized case, which up until now meant up and down. I'll try to look at it tomorrow, but I'm not sure if I will have time.

[zerothi]

I wouldn't bother understanding the current framework, as it will work in a very different (hopefully simpler) way :)

Ok. :)

Here I would just do the changes so that current functionality is not broken, and any new functionality should be added to the new design. I have not changed anything in the viz module since I created the PR because it is hard to keep track of changes to both designs at the same time.

Agreed. Lets wait with tbtrans then, but a check for the dimensions of the changes I did would be great.

Otherwise if you think the PR can be merged before the next release the changes can be done there.

I plan on getting this PR in before releasing, basically all I need is this PR 🤞

Apart from the changes that you did, I think that there should be changes to the handling of PDOS requests, because currently the spin key accepts 0 or 1 in the spin polarized case, which up until now meant up and down. I'll try to look at it tomorrow, but I'm not sure if I will have time.

Yes, this was where i was not sure how to do it, for instance the density plots still does up/down, whereas PDOS would now be total/z. And viz/input_fields/spin does not seem to distinguish this. What would you suggest?

[pfebrer]

The PDOS file read returns now total/z, right? So to be consistent one would have to calculate total/z from the Hamiltonian and the WFSX file reads, no?

[pfebrer]

The PDOS file read returns now total/z, right? So to be consistent one would have to calculate total/z from the Hamiltonian and the WFSX file reads, no?

Could there be a utility function to convert from one to the other (I know it's a simple operation but I think the code will look easier to understand if I call a function with an explicit name)

[zerothi]

The PDOS file read returns now total/z, right? So to be consistent one would have to calculate total/z from the Hamiltonian and the WFSX file reads, no?

yes, do PDOS for both spin, and convert to total/z

[pfebrer]

I can not push to this branch, if you give me rights I can try to do the changes

[zerothi]

you can make a PR against this branch, that should be easier for now. I cannot make a single branch open, AFAIK, since this is my private repo it complicates things, I should have begun with an organisational repo...

[pfebrer]

The PDOS file's read_data returns (no, nE) if the information is spin unpolarized. It should return (1, no, nE) for consistency, no?

[zerothi]

The PDOS file's read_data returns (no, nE) if the information is spin unpolarized. It should return (1, no, nE) for consistency, no?

thanks, fixed!

[zerothi]

@pfebrer now I get this error when running:

Here are the errors for each source:
	- bands file: ValueError.A path was provided, therefore we can not use the .bands file even if there is one
	- aiida bands: AttributeError.'NoneType' object has no attribute '_get_bandplot_data'
	- wfsx file: TypeError.expected str, bytes or os.PathLike object, not NoneType
	- band structure: ValueError.conflicting sizes for dimension 'axis': length 4 on 'spin_moments' and length 3 on {'k': 'k', 'band': 'band', 'spin': 'E', 'axis': 'axis'}. 

I tried to figure out where the 2nd array came from, but I can't seem to dig it out...

[pfebrer]

I think my changes were gone when you rebased and force pushed (?)

[zerothi]

that is weird, I thought rebase should always do the merges as well.. Ok, I'll have a look

[pfebrer]

But did you pull my changes locally?

[zerothi]

But did you pull my changes locally?

1. I merged your PR remotely
2. I updated my local branches
3. I then rebased on main...

At least this is what I think I did! :( I will merge your 501-viz-fixes branch locally, and push it, you don't need to do anything. Thanks!

[zerothi]

Just don't delete your branch ;)

[pfebrer]

Ok!

[zerothi]

I fixed this, cherry-picked your commit... Thanks!