Add new store_density_matricies Result option

[Matt-Ord]

Checklist
Thank you for contributing to QuTiP! Please make sure you have finished the following tasks before opening the PR.

• Please read Contributing to QuTiP Development
• Contributions to qutip should follow the pep8 style.
  You can use pycodestyle to check your code automatically
• Please add tests to cover your changes if applicable.
• If the behavior of the code has changed or new feature has been added, please also update the documentation in the doc folder, and the notebook. Feel free to ask if you are not sure.
• Include the changelog in a file named: doc/changes/<PR number>.<type> 'type' can be one of the following: feature, bugfix, doc, removal, misc, or deprecation (see here for more information).

Delete this checklist after you have completed all the tasks. If you have not finished them all, you can also open a Draft Pull Request to let the others know this on-going work and keep this checklist in the PR description.

Description
Add a new store_density_matricies option to Result

Related issues or PRs
Fixes #2299

[Matt-Ord]

@Ericgig does this look like it is along the right lines - It fixes my issue but I want it signed off before I go adding docs. Also unrelated what is your opinion on type annotations. Would it be accepted if I was to make another pr adding annotations to some functions I use?

[Ericgig]

Thank you for looking at this.

store_state should still be the final factor whether states are stored or not.
It's the one available everywhere and for v4. mesolve with store_state will have result.states being the output density matrix. mcsolve with the same store_state and default options should return a similar result.

So it would be better to have the new option be precompute_average_states and have it only control if the average is compute when all trajectories is stored. The if the average_states property it use, it could compute the average then, but never do so if not used.

Please keep result.states and result.final_states as is. They have the same behaviour as in v4 and we want to keep them for backward compatibility.

runs_... and average_... are new to v5 and could be renamed, but it's not clear that density_matricies is the average while final_states are for each trajectories. smesolve and sometime for mcsolve individual trajectories states are already density matrices.
Thus smesolve's result density_matricies could work as well for all states from all trajectories, causing confusion.

You implemented the density_matricies property twice.

Having the options to keep only the final states without keeping all the trajectories is a good idea. Be careful that the trajectory result have the final state stored

[Matt-Ord]

Thank you for looking at this.

store_state should still be the final factor whether states are stored or not. It's the one available everywhere and for v4. mesolve with store_state will have result.states being the output density matrix. mcsolve with the same store_state and default options should return a similar result.

So it would be better to have the new option be precompute_average_states and have it only control if the average is compute when all trajectories is stored. The if the average_states property it use, it could compute the average then, but never do so if not used.

Please keep result.states and result.final_states as is. They have the same behaviour as in v4 and we want to keep them for backward compatibility.

runs_... and average_... are new to v5 and could be renamed, but it's not clear that density_matricies is the average while final_states are for each trajectories. smesolve and sometime for mcsolve individual trajectories states are already density matrices. Thus smesolve's result density_matricies could work as well for all states from all trajectories, causing confusion.

You implemented the density_matricies property twice.

Having the options to keep only the final states without keeping all the trajectories is a good idea. Be careful that the trajectory result have the final state stored

So I'm not entirely sure what you want the behaviour to be for each set of options. I think it is confusing to return a density matrix if you request store_states - the density matrix is not the same as the state of the system, and surely v5 is the perfect time to clear up issues like this. Personally I found this behaviour very surprising when updating to V5. I do agree that density_matricies should maybe be average_density_matricies. I am happy to change mesolve as part of this for consistency. If you want to keep it backwards compatible would it be possible to send a list of all combination of options and exactly what you want states to return in each case, just to prevent any back and forth.

[Ericgig]

Here is how I see it:

store_states	keep_runs_results	precompute_average_states	->	states	average_states	runs_states
False	Any	Any	->	None	None	None
True	False	Any	->	Dms	Dms	None
True	True	False	->	[runs.states]	Computed in property	[runs.states]
True	True	True	->	[runs.states]	Computed at add	[runs.states]

mcsolve can call mesolve which can call sesolve. The options must match.

I asked around and states being the average did not confuse anyone, a kind of professional deformation.
It would break a lot of code and annoy a lot of users to change it from the only option to not existing.

density_matricies could be added as a new property.
Ideally it should be everywhere, for sesolve it would return the projection of the ket, etc.

[Matt-Ord]

Thanks for the detailed response. One issue - what if a user wants to store only the final state, not the final density matrix for example? For now I've made what I think should be the minimum fix - I've made it so that the density matrix is only pre-computed if keep_all_runs is not specified

[Ericgig]

The approach of only computing average per default after the evolution is good.
But with >300 lines of change, I don't see it as the minimum fix...

For keeping trajectories' last state without saving the full trajectories, I think we need a new option.

[Matt-Ord]

@Ericgig any updates on this? Also can you enable CI?

[Matt-Ord]

@Ericgig could you enable CI again

[coveralls]

coverage: 85.42% (-0.06%) from 85.479%
when pulling d550e7c on Matt-Ord:Split-density-matrix-and-state-results
into c3a1916 on qutip:master.

[Matt-Ord]

@Ericgig Should now be passing all tests, and essentially be ready

[Ericgig]

The code look nice.

Could you revert the names to what they were:
_reduce_density_matrix -> _reduce_states,
_sum_final_density_matrix -> _sum_final_states
etc.
We will be the one having to work with these in years to come, please to not impose your preferences on us.

Ideally I would prefer changes from black, flake8 etc. in their own PR as to not confuse with actual code changes. Also I do not blindly accept black formatting, it must read well for humans as well.

[Matt-Ord]

Thanks for your patience reviewing. I should have fixed all of your comments. Agree with the black formatting - my usual workflow is to format with black and add brackets until both black and I am happy with the result. Seems a few slipped through this time unfortunately.

[Ericgig]

Great work. Thank you.