JuliaApproximation / SemiclassicalOrthogonalPolynomials.jl

A Julia repository for semiclassical orthogonal polynomials
MIT License
7 stars 3 forks source link

Lowering connection cfs from raising connection cfs #112

Closed TSGut closed 2 months ago

TSGut commented 2 months ago

Not sure whether we will end up merging this or your PR @DanielVandH but the point of this is to show how a lowering implementation would look. It actually works already but it's not efficient.

julia> @testset "basic lowering" begin
           # lowering
           t = 2.2
           P = SemiclassicalJacobi(t, 2, 2, 2)
           R_0 = SemiclassicalJacobi(t, 1, 2, 2) \ P
           R_1 = SemiclassicalJacobi(t, 2, 1, 2) \ P
           R_t = SemiclassicalJacobi(t, 2, 2, 1) \ P

           R_0inv = P \ SemiclassicalJacobi(t, 1, 2, 2)
           R_1inv = P \ SemiclassicalJacobi(t, 2, 1, 2)
           R_tinv = P \ SemiclassicalJacobi(t, 2, 2, 1)

           @test ApplyArray(inv,R_0inv)[1:20,1:20] ≈ R_0[1:20,1:20]
           @test ApplyArray(inv,R_1inv)[1:20,1:20] ≈ R_1[1:20,1:20]
           @test ApplyArray(inv,R_tinv)[1:20,1:20] ≈ R_t[1:20,1:20]
       end
Test Summary:  | Pass  Total  Time
basic lowering |    3      3  0.0s

The reason for the inefficiency lies entirely with three lines where I 'cheat'. The iterative connection coefficients above are computed correctly, that's not the problem. The problem is in the computation we need Q.X which the user could either supply always (this would correspond to @dlfivefifty 's suggestion of changing all conventions to be an semiclassical_ldiv(Q,X) format or alternatively we can implement the reverse Cholesky Jacobi matrix which would go in the big marked by TODO section, i.e. this bit

 elseif isone(-Δa) && iszero(Δb) && iszero(Δc)  # raising by 1
        # TODO: This is re-constructing. It should instead use reverse Cholesky (or an alternative)!
        semiclassical_jacobimatrix(Q.t,a,b,c)
    elseif iszero(Δa) && isone(-Δb) && iszero(Δc)
        # TODO: This is re-constructing. It should instead use reverse Cholesky (or an alternative)!
        semiclassical_jacobimatrix(Q.t,a,b,c)
    elseif iszero(Δa) && iszero(Δb) && isone(-Δc)
        # TODO: This is re-constructing. It should instead use reverse Cholesky (or an alternative)!
        semiclassical_jacobimatrix(Q.t,a,b,c)

If the above was reverse Cholesky jacobimatrix then this PR would do everything you want I think. Without that, we have a design decision to make.

TSGut commented 2 months ago

Guess something else broke as a result, I will fix that later today.

DanielVandH commented 2 months ago

Thanks a lot @TSGut, this seems to be a better approach to it.

The problem is in the computation we need Q.X which the user could either supply always (this would correspond to @dlfivefifty 's suggestion of changing all conventions to be an semiclassical_ldiv(Q,X) format or alternatively we can implement the reverse Cholesky Jacobi matrix

I think it would be easier to work with just being able to keep Q around. Do you have a preference? Maybe makes this easier in the future too if ever we need to come back to these functions and again have a need for Q.

TSGut commented 2 months ago

I think either is fine, there are some cost savings here and there depending on the use case with either approach but probably all negligible.

TSGut commented 2 months ago

Alright, tests passing (though I may add more).

In this package I am usually happy to just merge things myself but since you guys are actively developing with this, I would prefer not to mess unless this is wanted, so please take a look @dlfivefifty and @DanielVandH.


Changes and new features:

TSGut commented 2 months ago

Any idea why codecov went and died on this repo? It would be nice to add more tests to this but it isn't being displayed and manually checking https://app.codecov.io/gh/JuliaApproximation/SemiclassicalOrthogonalPolynomials.jl/pull/112 shows some error in the config that started a few PRs ago

dlfivefifty commented 2 months ago

At some point I had to add keys for codexov, probably I haven’t done this yet. I’ll check when I’m home

dlfivefifty commented 2 months ago

codecov shouldb e fixded

dlfivefifty commented 2 months ago

can you remove the WIP when this is ready?

codecov[bot] commented 2 months ago

Codecov Report

Attention: Patch coverage is 98.76543% with 1 line in your changes missing coverage. Please review.

Project coverage is 92.52%. Comparing base (5d48cc5) to head (f4e9683). Report is 3 commits behind head on master.

Files Patch % Lines
src/SemiclassicalOrthogonalPolynomials.jl 98.76% 1 Missing :warning:
Additional details and impacted files ```diff @@ Coverage Diff @@ ## master #112 +/- ## ========================================== + Coverage 92.40% 92.52% +0.11% ========================================== Files 3 3 Lines 527 575 +48 ========================================== + Hits 487 532 +45 - Misses 40 43 +3 ```

:umbrella: View full report in Codecov by Sentry.
:loudspeaker: Have feedback on the report? Share it here.

dlfivefifty commented 2 months ago

also get the codecov up

TSGut commented 2 months ago

Ok I think this is ready if you're happy with it @dlfivefifty. Codecov is good and WIP tag removed. Version was already previously bumped.

Once reverse Cholesky jacobimatrix is implemented, which if nobody beats me to it I will get to this summer, I will also update this package to remove the fallback option.

I am excited to see the negative parameter stuff @DanielVandH! Let me know if I can help any more regarding this package.

DanielVandH commented 2 months ago

Great! Thanks so much @TSGut for the help here.