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

Distribution function #141

Merged
merged 3 commits into from
Oct 6, 2017

Conversation

antoinetavant
Copy link
Contributor

First contribution for me !
To start with #101, I added a 1D Maxwellian distribution function calculation.
Please tell me if you like the interface, numbre of arguments and so on.

Needs to be done :

  • Do the same with Kappa
  • implement a 2D and/or 3D version of it.
  • add particle generator
  • maybe : add symbolic formulation

@pep8speaks
Copy link

pep8speaks commented Oct 4, 2017

Hello @antoinelpp! Thanks for updating the PR.

Cheers ! There are no PEP8 issues in this Pull Request. 🍻

Comment last updated on October 06, 2017 at 15:48 Hours UTC

Copy link
Member

@StanczakDominik StanczakDominik left a comment

Choose a reason for hiding this comment

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

Hey, thanks for this one! I like it. There are some minor improvements I could see being useful.

I can't find a good place to put comments on the Jupyter notebook, so I'll put them here:

assert np.isclose(v_vect[Maxwellian_1D(v_vect,T=T_e, part='e',V_drift = V_drift ).argmax()].value,V_drift.value)

#integral of the distribution over v_vect
integral = (Maxwellian_1D(v_vect,T= 30000*u.K, part='e')).sum()*dv
Copy link
Member

Choose a reason for hiding this comment

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

If we used scipy.integrate here (I don't think adding scipy to requirements is an issue given we're working within the scientific Python stack), we'd also get an error estimate on the integral we could use for the next line!

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I'n not so used to Scipy, I didn't knew we could do it ! I'll check it out.

from ..distribution import (Maxwellian_1D)

T_e = 30000*u.K;
V_drift = 10*u.m/u.s;
Copy link
Member

Choose a reason for hiding this comment

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

I haven't looked at the speed of this test, but if it's not too slow (it's an integral, it shouldn't be) this would be a perfect place to test a few parameter combinations via pytest fixtures. I'd be happy to add them!

The temperature in Kelvin

part: string, optional
Representation of the ion species (e.g., 'p' for protons, 'D+'
Copy link
Member

Choose a reason for hiding this comment

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

Could we rename this to particle, and say Representation of the particle species?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Ok for the argument name !


else :
try:
m_s = ion_mass(ion)
Copy link
Member

Choose a reason for hiding this comment

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

This should be possible to clean up a bit once #140 goes through :)

Representation of the ion species (e.g., 'p' for protons, 'D+'
for deuterium, or 'He-4 +1' for singly ionized helium-4),
which defaults to electrons.

Copy link
Member

Choose a reason for hiding this comment

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

We could probably use a latex (mathrm) representation of the formula. I can get that done for you.

Come to think of it, AFAIK SymPy can autogenerate those representations... #136 ? :D

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I'm not so sure to understand this comment. I tried something tell me if this is what you think of.

Copy link
Member

Choose a reason for hiding this comment

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

My apologies, I guess I wasn't clear at all here! I added a commit to show what I mean (the actual formula for the Maxwellian distribution) within the docstring.

from ..utils import _check_quantity, check_relativistic, check_quantity


def Maxwellian_1D(v,T,part ="e", V_drift = 0* units.m/units.s):
Copy link
Member

Choose a reason for hiding this comment

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

Seems like a nice place for @astropy.units.quantity_input :D

Copy link
Contributor Author

Choose a reason for hiding this comment

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

good point !

@StanczakDominik
Copy link
Member

That's a weird test failure:

063     Example
064     -------
065     >>> from plasmapy.physics import Maxwellian_1D
066     >>> from astropy import units as u
067     >>> Maxwellian_1D(v=1*u.m/u.s, T= 30000*u.K, particle='e',V_drift=0*u.m/u.s)
Expected:
    <Quantity 5.916 m/s >
Got:
    <Quantity 5.916329687405701e-07 s / m>

so you expect something on the order of 6m/s and get 6 * 10^-7 s/m ? This seems like some bug in the formula, I'll take a look.

@StanczakDominik
Copy link
Member

Okay, so I just noticed that this returns return f.to(units.s/units.m). Is this correct? It's 1/velocity, so given this is a probability distribution function it probably is correct as you can integrate it over velocities and get a 1 without units. The formula seems all right, though.

I also noticed I forgot to add k_B for each temperature in the formula docstring, sorry :(

@antoinetavant
Copy link
Contributor Author

antoinetavant commented Oct 4, 2017

so you expect something on the order of 6m/s and get 6 * 10^-7 s/m ? This seems like some bug in the formula, I'll take a look.

No that's my mistake !

@antoinetavant
Copy link
Contributor Author

Sorry all for the PR spam, but I couldn't figure out how to run the tests without pushing...
Anyone could tell me the best way to run them from my repo ?

I tried to rebase it on Upstream, but the PR is now a mess... Shouldn't I do it ? Or I maybe I did it wrong...

@StanczakDominik
Copy link
Member

I wouldn't worry about rebasing for now. I think we can live with a lengthy history.

We really need to look into this: run pytest plasmapy from repository root. Could someone start an issue for this? I can't right now.

@antoinetavant
Copy link
Contributor Author

Ho !
Actually pytest plasmapy works fine ! I didn't knew about it.
I'm sorry :s

@StanczakDominik
Copy link
Member

No worries, glad it works for you :)

I'll try to do a final review on this but I won't be able to do that by Friday, probably. If anyone wants to step in, I'd be much obliged!

Copy link
Member

@StanczakDominik StanczakDominik left a comment

Choose a reason for hiding this comment

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

You know, I checked and this did indeed turn out to be a lot of commits. The review screen on GitHub at least turned out weird and it looks like you're kind of duplicating changes from #140. That may be a glitch. However, could you please still try rebasing onto current master and squashing some of those smaller commits after all?

I'll go through the notebook as soon as time permits :)

@@ -26,3 +21,5 @@
from .relativity import Lorentz_factor

from .transport import Coulomb_logarithm

from .distribution import (Maxwellian_1D)
Copy link
Member

Choose a reason for hiding this comment

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

Why the ( ) parentheses here?

Quantity)

from ..constants import (m_p, m_e, c, mu0, k_B, e, eps0, pi, e)
from ..atomic import (ion_mass, charge_state)
Copy link
Member

Choose a reason for hiding this comment

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

Once again, ( ) ?

Copy link
Member

Choose a reason for hiding this comment

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

I've made a habit of using parentheses when doing long import statements, as per PEP 328. At least to me, this makes it easier to read. It is at its most useful when import statements span multiple lines.

Copy link
Member

Choose a reason for hiding this comment

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

Oh right. I'm sorry, I forgot about this and I was checking this in a slight rush. Multilines, sure. Here, not particularly certain they're quite necessary!


# Get mass

if particle == "e":
Copy link
Member

Choose a reason for hiding this comment

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

Another way to do this is to import _is_electron from atomic and use that here:

    if _is_electron(particle):
        m_s = m_e

This will catch when e- or electron is used to represent electrons as well.

Thank you for doing this!

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Good point !
I'm not fluant with the atomic package... I should take a better look into it.

@StanczakDominik
Copy link
Member

@antoinelpp I took a closer look at the branch history and moanad, because something must have gone horribly wrong there. That tree was duplicated in some weird ways. I tried everything, but in the end the only thing that made sense was to fix some of that tangled mess via rebase - you can see the result of that here, but at the same time I may have messed up some of the chronological ordering of the changes.

If you think we can let it go and let that previous messed up version rest under the sea, I made another version where it's all squashed into a single commit. Doing it that way means no worries about chronology for this PR.

@antoinetavant
Copy link
Contributor Author

Hi @StanczakDominik !

I took a closer look at the branch history and moanad, because something must have gone horribly wrong there.

My bad :( I sorry about it.

we can let it go and let that previous messed up version rest under the sea,

No no no, I'd like not to mess with the branch !

I made another version where it's all squashed into a single commit.

Ok that looks way better ! I can try to squash it myself (on my repo) in order to fix the other issues. Or I can remove this PR, clean my branch and do it again... What do you think ?

Again sorry

@StanczakDominik
Copy link
Member

StanczakDominik commented Oct 6, 2017 via email

@StanczakDominik
Copy link
Member

StanczakDominik commented Oct 6, 2017 via email

created distribution.py
created test_distribution.py
change `part` argument to `particle` and minor docstrings
Added is_electron


Update formula docstring
@StanczakDominik
Copy link
Member

Awesome, I'll just do a short stylistic pass on the example (TBD within two hours) and we're good to go :)

@StanczakDominik
Copy link
Member

All right, take a look at this and once it's done I'd say we can merge this one :)

@antoinetavant
Copy link
Contributor Author

Thanks @StanczakDominik !
I just changed the V_drift parameter in the last exemple. That's look good to me. I think we can merge it as well !

@StanczakDominik StanczakDominik merged commit 53da3c5 into PlasmaPy:master Oct 6, 2017
@StanczakDominik
Copy link
Member

And I merged it, thanks a bunch for this one!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
plasmapy.formulary Related to the plasmapy.formulary subpackage
Projects
Development

Successfully merging this pull request may close these issues.

None yet

4 participants