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

Lambert's Omega constant #12

Closed
wants to merge 11 commits into from
Closed

Conversation

alyst
Copy link

@alyst alyst commented Nov 19, 2021

Add Lambert's Omega constant as a part of moving LambertW.jl into SpecialFunctions.jl (see JuliaMath/SpecialFunctions.jl/pull/84).

Also add the related inve = exp(-1) constant, since -exp(-1) is the branching point for W_0(t) and W_{-1}(t).

cc @jlapeyre

useful in many situations, also -exp(-1) is the branching point of
Lambert's W_0(t) and W_{-1}(t)
@codecov
Copy link

codecov bot commented Nov 19, 2021

Codecov Report

Merging #12 (1ee42a5) into main (d226b95) will increase coverage by 80.00%.
The diff coverage is 80.00%.

Impacted file tree graph

@@            Coverage Diff            @@
##           main      #12       +/-   ##
=========================================
+ Coverage      0   80.00%   +80.00%     
=========================================
  Files         0        1        +1     
  Lines         0       10       +10     
=========================================
+ Hits          0        8        +8     
- Misses        0        2        +2     
Impacted Files Coverage Δ
src/lambertw_omega.jl 80.00% <80.00%> (ø)

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update d226b95...1ee42a5. Read the comment docs.

@coveralls
Copy link

coveralls commented Nov 19, 2021

Pull Request Test Coverage Report for Build 1481553290

Warning: This coverage report may be inaccurate.

This pull request's base commit is no longer the HEAD commit of its target branch. This means it includes changes from outside the original pull request, including, potentially, unrelated coverage changes.

Details

  • 10 of 10 (100.0%) changed or added relevant lines in 1 file are covered.
  • No unchanged relevant lines lost coverage.
  • Overall coverage remained the same at 100.0%

Totals Coverage Status
Change from base Build 1333249280: 0.0%
Covered Lines: 10
Relevant Lines: 10

💛 - Coveralls

src/IrrationalConstants.jl Outdated Show resolved Hide resolved
log4π, # log(4π)
invℯ, # 1 / ℯ

LambertW_Ω # Ω exp(Ω) = 1
Copy link
Member

Choose a reason for hiding this comment

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

We should stick with the convention in this package and in base, and use lowercase names and in particular no camel case. I.e., e.g. lambertwω or just ω. Again, also we might want to replace it with or at least add an ASCII alias.

Copy link
Author

Choose a reason for hiding this comment

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

Is it "ω" or "Ω" (I've taken the latter from https://en.wikipedia.org/wiki/Omega_constant)? I think upper/lowercase for constants has a difference (π vs Π).
I also prefer Julian convention to avoid underscore in the names, but in this case mixing latin and greek doesn't contribute to readability. What about lambertw_Ω?
I also think ASCII alias (lambertw_Omega) would be very useful, what's the conventional way to define it?

Copy link
Member

Choose a reason for hiding this comment

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

Since it is called "Omega constant" I think we should just use Ω (or ω if we want to stick with the convention in Base and this package of only using lowercase names). And correspondingly, the ASCII alias would just be Omega or omega.

Copy link
Author

Choose a reason for hiding this comment

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

In the revised version it's called lambertw_Ω and lambertw_Omega alias is added. I think plain Ω would be confusing, since it's less well known and/or may interfere with variable name in the user code.
lambertwΩ IMO is quite unreadable.

Copy link

@jlapeyre jlapeyre Dec 31, 2021

Choose a reason for hiding this comment

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

All of the constants in MathConstants and IrrationalConstants follow a naming convention. The omega constant should follow the same convention. None of your proposals follow this convention. But, in any case, it should stay in LambertW.jl.

Comment on lines 2 to 19
const LambertW_Omega_BigFloat256 = Ref{BigFloat}()

# compute BigFloat Omega constant at arbitrary precision
function compute_LambertW_Omega()
# initialize Omega_BigFloat256
isassigned(LambertW_Omega_BigFloat256) ||
(LambertW_Omega_BigFloat256[] = BigFloat("0.5671432904097838729999686622103555497538157871865125081351310792230457930866845666932194"))
o = LambertW_Omega_BigFloat256[] # initial value
precision(BigFloat) <= 256 && return o
# iteratively improve the precision of the constant
myeps = eps(BigFloat)
for _ in 1:100
o_ = (1 + o) / (1 + exp(o))
abs(o - o_) <= myeps && break
o = o_
end
return o
end
Copy link
Member

Choose a reason for hiding this comment

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

I'm not sure if it's worth caching the value in LambertW_Omega_BigFloat256. I don't think this is done for other constants? Also it seems it would return a value of a different precision than requested in the case precision(BigFloat) < 256.

Copy link
Author

Choose a reason for hiding this comment

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

The other constants either have corresponding mpfr_const_xxx() functions or simple formulas.
But I agree it would be more straightforward to have just

    o = BigFloat("0.5671432904097838729999686622103555497538157871865125081351310792230457930866845666932194") # initial value with 256-bit precision
    precision(BigFloat) <= 256 && return o

if BigFloat("...") is fast enough.

Copy link
Author

Choose a reason for hiding this comment

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

Now the returned BigFloat constant should respect the current precision better, and the appropriate tests were added.
I still have kept the lambertw_Omega_BigFloat256 as I assume it's faster to reuse the initialized BigFloat than parse it each time.

Choose a reason for hiding this comment

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

This is much slower than the code in LambertW.jl. The two should be benchmarked explicitly before a PR would be accepted.

But, I agree with the others that it does not make sense to move the omega constant out of LambertW.jl. It should stay in LambertW.jl and the latter should be moved into SpecialFuncions

Comment on lines 23 to 30
"""
Lambert's Omega (Ω) constant, such that Ω exp(Ω) = 1.

*W(Ω) = 1*, where *W(t) = t exp(t)* is the *Lambert's W function*.

# See
https://en.wikipedia.org/wiki/Omega_constant
"""
Copy link
Member

Choose a reason for hiding this comment

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

The docstring should follow the conventions and recommendations in the Julia documentation: https://docs.julialang.org/en/v1/manual/documentation/

Copy link
Author

Choose a reason for hiding this comment

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

Could you please suggest your variant?

Copy link
Author

Choose a reason for hiding this comment

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

I have improved the docstring. But if you think it does not conform to the guidelines, please let me know of the specific changes you want to implement or just fix it yourself.

https://en.wikipedia.org/wiki/Omega_constant
"""
LambertW_Ω

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

Copy link
Member

@devmotion devmotion left a comment

Choose a reason for hiding this comment

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

Can you remind me why this constant should be added to IrrationalConstants? Isn't it sufficient to add it to SpecialFunctions where it is needed?

precision(o) <= 256 && return o
# iteratively improve the precision of the constant
myeps = eps(BigFloat)
for _ in 1:100
Copy link
Member

Choose a reason for hiding this comment

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

Is this enough for higher precision? It seems a bit dangerous to upper bound the number of iterations.

Suggested change
for _ in 1:100
while true

or alternatively move the stopping criterion here.

Copy link
Author

Choose a reason for hiding this comment

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

Wiki says (I've just added the reference) that this is the quadratic method, the number of correct digits is doubled at each iteration. So 100 should be safe for any reasonable precision.
while true sounds a bit dangerous.
I can increase it to 1000 and add a warning that the precision was not reached.

Copy link
Member

Choose a reason for hiding this comment

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

Maybe just use the stopping criterion, i.e., while abs(next_o - o) > eps?

Copy link
Author

Choose a reason for hiding this comment

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

I want to avoid dead loop if there's some subtle bug in BigFloat implementation or something like that.

Copy link
Member

Choose a reason for hiding this comment

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

What bug could cause a dead loop? Generally, we should assume that BigFloat does not contain any bugs - and if we notice any, they should be fixed upstream.

BTW just came across the following while true loop, I still think it would be a simple and straightforward implementation here as well: https://github.com/JuliaLang/julia/blob/3d11f7db65a3461320542aee3b0f26619c4e65e3/base/irrationals.jl#L54

@alyst
Copy link
Author

alyst commented Dec 30, 2021

Can you remind me why this constant should be added to IrrationalConstants? Isn't it sufficient to add it to SpecialFunctions where it is needed?

I guess it's optional to add it to IrrationalConstants, but isn't it the exact purpose of this package to have all well-accepted constants here? (given that the code to calculate them is simple enough and doesn't introduce additional dependencies)
Omega doesn't introduce any dependencies to IrrationalConstants. Neither putting it here would require introducing SpecialFunctions -> IrrationalConstants dependency: it is already there.

I don't have any real use-cases in mind (besides some educational), but MathWorld page shows that the constant possesses the properties outside of the scope of Lambert's W (golden ratio of exponents, power tower).

@devmotion
Copy link
Member

isn't it the exact purpose of this package to have all well-accepted constants here?

The purpose is to avoid duplicate definitions in multiple packages by defining commonly used constants in a lightweight package. If a constant is only used in a single package there is no need to add it to IrrationalConstants - there is no benefit over defining them in the package where it is used but it might increase compilation and loading times for all other packages that depend on IrrationalConstants.

I think we need at least 2, better 3, packages that would make use of it.

@alyst
Copy link
Author

alyst commented Dec 30, 2021

The purpose is to avoid duplicate definitions in multiple packages by defining commonly used constants in a lightweight package.

I think I'm repeating myself, but

  • with this PR the package stays lightweight, and there is no duplication of code for calculation of lambertw() and Omega constant
  • with this PR no additional dependencies are introduced anywhere
  • there are potential uses of this constant outside of lambertw()

Do you think being so conservative about this package is the right thing for the Julia ecosystem?
Not to mention that you could have said that you are hesitant to introduce the constant when I have just submitted the PR, so that I don't waste my time improving it.
I'm working on the updated PR for lambertw() in SpecialFunctions.jl, so it would be nice to resolve this issue rather soon.

cc @lindahua @JeffreySarnoff

@devmotion
Copy link
Member

Do you think being so conservative about this package is the right thing for the Julia ecosystem?

Yes. I'm fine with adding new constants (proposed and added some myself) if they are used in multiple downstream packages. The package is used by large parts of the Julia ecosystem, both directly and indirectly, so I don't think increased compilation and load times are acceptable if it is only used in a single package.

there are potential uses of this constant outside of lambertw()

Are there any concrete packages in addition to SpecialFunctions where you plan to use it?

Not to mention that you could have said that you are hesitant to introduce the constant when I have just submitted the PR, so that I don't waste my time improving it.

You would want to use these improvements in any case, regardless of whether the code ends up in SpecialFunctions or IrrationalConstants. I don't think it was a waste of time to improve it.

@jlapeyre
Copy link

I agree with @devmotion . This constant belongs in LambertW.jl.

If LambertW.jl is moved to SpecialFunctions, I'm not sure how to organize the name. Currently, I'd recommend LambertW.omega. But, SpecialFunctions.omega is a step down in clarity and conservative name spacing.

@alyst
Copy link
Author

alyst commented Dec 31, 2021

Currently, I'd recommend LambertW.omega. But, SpecialFunctions.omega is a step down in clarity and conservative name spacing.

LambertW.Omega is a good idea. I guess it should be possible to do it here or in SpecialFunctions.jl with

....
export LambertW
module LambertW
    @irrational Omega ....
    ....
end
....

@devmotion
Copy link
Member

@irrational Omega ... will define Base.Irrational{:Omega} regardless of the module where it is defined in, scoping does not matter. (Therefore actually every use of @irrational outside of Julia base is type piracy but this is something we accept in this package so far.)

@alyst
Copy link
Author

alyst commented Dec 31, 2021

@irrational Omega ... will define Base.Irrational{:Omega}

ah, true. My bad

@alyst
Copy link
Author

alyst commented Dec 31, 2021

then it could be something like

....
@irrational lambertw_Omega .... # private declaration
export LambertW
module LambertW
    const Omega = ..lambertw_Omega # public alias
    ....
end
....

@devmotion
Copy link
Member

@irrational lambertw_Omega .... # private declaration

@irrational is never a private declaration, it can be accessed with Irrational{:lambertw_Omega}() regardless of exports and imports whenever the defining package is loaded (possibly even only indirectly by a dependency).

@alyst
Copy link
Author

alyst commented Dec 31, 2021

@irrational is never a private declaration, it can be accessed with Irrational{:lambertw_Omega}()

Sure. My idea is that with the suggested scheme the identifier used by Irrational{} is more specific than :Omega (avoiding potential clashes with the other packages).
It is "private" in the sense that it is not advertised to be used directly.

@alyst
Copy link
Author

alyst commented Jan 2, 2022

Since the discussion about Lambert W function has been reactivated (JuliaMath/SpecialFunctions.jl/issues/84), and I have prepared an updated version of that PR (JuliaMath/SpecialFunctions.jl/issues/371) that expects Omega in IrrationalConstants.jl, I think we need to make a decision rather soon.

I am fine with moving Omega to SpecialFunctions.jl or LambertW.jl, I just find that the admittance bar (several packages depending on the constant) is quite high, so only the constants that are transformations of Pi or Euler would make it here.
At the same time for "discoverability" it would be nice to have Julian "repository of constants" as a part of its ecosystem infrastructure.

cc @ViralBShah

@devmotion
Copy link
Member

My stance is clear: I'm fine with adding the constant if it is used by other packages apart from LambertW/SpecialFunctions (which obviously don't depend on these packages). Otherwise there is no need to add them to IrrationalConstants since the purpose of this package is not to define all irrational constants used in Julia packages but to avoid redundant definitions in independent packages.

(BTW I also don't think there's a need to rush this PR or the ones in SpecialFunctions since there is no missing functionality or bugs addressed by them - one can use the Lambert W function with the LambertW package without any issues as far as I know.)

@ViralBShah
Copy link
Member

I don't have a strong view one way or another. I was mainly going through old issues and PRs and seeing what could be closed, salvaged, etc.

It seems like the preference is to have a separate package.

@ViralBShah
Copy link
Member

I prefer packages to be in orgs, so that more people can help with ongoing maintenance tasks. Naturally, the package does need the original author's energy and contributions for the most part.

@alyst alyst mentioned this pull request Jan 3, 2022
@alyst
Copy link
Author

alyst commented Jan 3, 2022

Ok, I'm closing this PR and opening the one for adding just inve = 1/e.

@alyst alyst closed this Jan 3, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Jan 3, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Jan 3, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Jan 3, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request May 23, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Jun 24, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Aug 2, 2022
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Mar 7, 2023
alyst added a commit to alyst/SpecialFunctions.jl that referenced this pull request Mar 24, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

5 participants