by

Where communities thrive


  • Join over 1.5M+ people
  • Join over 100K+ communities
  • Free without limits
  • Create your own community
People
Repo info
Activity
    Kirill Ignatiev
    @ikirill
    It would make the self-adjoint ODE actually have a symmetric matrix when solving the generalized eigenvalue problem, but I couldn't find this.
    Sheehan Olver
    @dlfivefifty
    It's just an operator that hasn't been added yet. Probably there should be a Normalized(Jacobi(1,1)) space and this done as a conversion operator.
    (Before you ask "why not always use normalized", its because the unnormalized versions have rational formulae which are much faster.)
    Mikael Slevinsky
    @MikaelSlevinsky
    I showed how to get a symmetric-definite + banded GEP for associated Legendre functions here (section 3.1 https://arxiv.org/abs/1711.07866). I have some ApproxFun-like code that I wrote to verify the correctness of the formulas, (but not the FMM speedup)
    Mikael Slevinsky
    @MikaelSlevinsky
    Sorry, that should say "one way" to get banded symmetry
    Ashton Bradley
    @AshtonSBradley
    if I have 2D cartesian grid data at a certain resolution, is there any way to make a high precision Chebychev interpolant using ApproxFun?
    Ashton Bradley
    @AshtonSBradley
    ah, I see a comment above that answers this!
    Soham Mukherjee
    @soham1112
    Hi, I was just trying out an example from https://github.com/JuliaApproximation/ApproxFunExamples/blob/master/PDEs/Rectangle%20PDEs.ipynb. When trying the example for solving the convection equation (see In[10]), I end up getting an error if I end up changing the interval for t, say, for example, from 0 .. 2 to 2 .. 3. ERROR: LoadError: LoadError: Cannot convert coefficients from Chebyshev(-1.0..1.0) to Chebyshev(-1.0..1.0)⊗ConstantSpace(Point(2.0))
    I don't understand why this should be giving an error, since I should be able to choose an arbitrary time-interval. Am I missing something?
    Sheehan Olver
    @dlfivefifty
    It’s using an embedding of a 1D interval in the 2D plane where the interval is assumed to be at t = 0. You can use a fully 2D version: u0 = Fun((x,_) -> ..., Chebyshev(-1.0..1.0)⊗ConstantSpace(Point(2.0)))
    Soham Mukherjee
    @soham1112
    Ah! I see. Thanks--that worked like a charm! However, I'm curious if there's a reason for this default behaviour?
    Mitko Georgiev
    @mitkoge
    Hi,
    i would like to ask a question about joined spaces
    I tried tzhe following
    using ApproxFun
    sp1= Chebyshev(0.0..50.0) ∪ Chebyshev(50.0..60.0)
    x1=Fun(sp1); g1 = Fun(exp(-x1^2),sp1)
    s1= sum(g1) #NaN
    sp2= Chebyshev(0.0..5.0) ∪ Chebyshev(5.0..60.0)
    x2=Fun(sp2); g2 = Fun(exp(-x2^2),sp2)
    s2= sum(g2) #0.8862269254527582
    i expected s1~= s2. but it seems it is not.
    Mikael Slevinsky
    @MikaelSlevinsky
    I pushed a fix for this on this branch JuliaApproximation/ApproxFun.jl@07f50d0. If the tests pass, then it could be merged.
    Mitko Georgiev
    @mitkoge
    Wow i am impressed by your reactivity! Fixed, pushed and merged in a day! Thank you!
    Mitko Georgiev
    @mitkoge

    May i ask another question?
    I try to reproduce an example from this chat given by @cortner

    using ApproxFun
    dom = Interval(0.001, 1) * PeriodicInterval(-pi, pi)
    space = Space(dom)
    Dr = Derivative(space, [1,0])
    Dθ = Derivative(space, [0,1])
    Mr = Multiplication(Fun( (r, θ) -> r, space ), space)
    rDr = Mr * Dr
    L = rDr * rDr + Dθ * Dθ

    but PeriodicInterval(-pi, pi) is no more recognized.
    What alternative can be used? Circle?

    Sheehan Olver
    @dlfivefifty
    It’s been renamed PeriodicSegment as it’s really an oriented line segment
    Mitko Georgiev
    @mitkoge

    Thank you,
    i tried

    julia> dom = Interval(0.001, 1) * PeriodicSegment(-pi, pi)
    ERROR: MethodError: no method matching *(::Interval{:closed,:closed,Float64}, ::PeriodicSegment{Floa
    t64})

    and

    julia> dom = Segment(0.001, 1) * PeriodicSegment(-pi, pi)
    ERROR: MethodError: no method matching *(::Segment{Float64}, ::PeriodicSegment{Float64})
    Sheehan Olver
    @dlfivefifty
    Ah and * has changed to imes
    imes
    \times
    Mitko Georgiev
    @mitkoge

    thank you.

    using ApproxFun,LinearAlgebra
    dm = Segment(0.001,1) × PeriodicSegment(-pi,pi)
    sp = Space(dm)
    Dr = Derivative(sp, [1,0]); Dθ = Derivative(sp, [0,1])
    Mr = Multiplication(Fun( (r, θ) -> r, sp ), sp)
    rDr = Mr * Dr
    Lr = rDr * rDr; L = Lr + Dθ * Dθ

    i modified it and now can play with.

    Simon Etter
    @ettersi

    Not sure whether this is the right place to ask since technically my question concerns SingularIntegralEquations.jl, but that package seems to be related enough to ApproxFun.jl that I'll go ahead anyway.

    In the code of SingularIntegralEquations it says

    # stieltjesintegral is an indefinite integral of stieltjes
    # normalized so that there is no constant term

    Can someone elaborate on that normalisation condition?

    Sheehan Olver
    @dlfivefifty
    Sure, fine to ask here. This normalisation condition comes from the fact that Jacobi polynomials have nice integration formulae coming from the weighted differentiation d/dx[(1-x)^a *(1+x)^b P_n^(a,b)] = C (1-x)^{a-1} (1+x)^{b-1} P_{n+1}^(a-1,b-1). Roughly speaking, for functions that vanish at ±1, derivatives and Stieltjes transforms commute, so provided a,b > 0 we can use this to write the integral of a Stieltjes transforms in terms of the integral of the integrand.
    Note that the n = 0 case is special, and dictates the normalisation constant (as the construction for other n will automatically decay something like O(z^(-n)) )
    For the cases actually implemented, I think its normalized to go to zero at ∞.
    Sorry, that's not right.
    You always have logarithmic growth (unless the integrand integrates to zero).
    Sheehan Olver
    @dlfivefifty
    On -1..1 it just vanishes apart from the log term:
    julia> x = Fun();
    
    julia> f = exp(x)/sqrt(1-x^2);
    
    julia> z = 10_000; stieltjesintegral(f,z) - log(z) * sum(f)
    -0.00017756097918919522
    julia> f = exp(x) * sqrt(1-x^2);
    
    julia> z = 10_000; stieltjesintegral(f,z) - log(z) * sum(f)
    -4.264886882765495e-5
    Oh, wait, that's true for other intervals too:
    
    julia> x = Fun(1..2);
    
    julia> f = exp(x) * sqrt((x-1)*(2-x));
    
    julia> z = 100_000; stieltjesintegral(f,z) - log(z) * sum(f)
    -2.8356239955229512e-5
    Sorry for making it seem more complicated than it actually it is. It's really simple: it's normalized so that stieltjes(f,z) = log(z) *sum(f) + o(1).
    Which happens to be consistent with logkernel (whose normalization is imposed by the actual definition):
    julia> logkernel(f, z) -log(abs(z))/π * sum(f)
    -9.026071512430178e-6
    Simon Etter
    @ettersi
    If I understand correctly, I can summarise what you said as: stieltjes(f,z) = integral log(z - w) f(w) dw + C with Csuch that stieltjes(f,z) = log(z) *sum(f) + o(1) for z to Inf. Correct?
    Sheehan Olver
    @dlfivefifty
    Yes
    Simon Etter
    @ettersi
    Great, thanks a lot!
    Mitko Georgiev
    @mitkoge

    Hi i was playing with modified example from Quantum states.jl

    Using  ApproxFun
    fc1 = Fun(1, 0..1); fc2 = Fun(2, 1..2); fc3 = Fun(3, 2..3)
    
    fc= fc1+fc2+fc3
    D = Derivative(); sp = space(fc); x= Fun(sp);
    B = [Dirichlet(sp); continuity(sp,0:1)]
    λ, v = ApproxFun.eigs(B, -D^2 + fc, 500,tolerance=1E-10)  
    #returns result
    
    λ, v = ApproxFun.eigs(B, -D^2 - 1/x*D + fc, 500,tolerance=1E-10)
    # says Not implemented  but if fc=fc1+fc2 it manages to return result

    Is 1/x*D kind of a problematic composition?

    Sheehan Olver
    @dlfivefifty
    Yes. It would be easy to work around by multiplying by x and updating the generalised eigenvalue call appropriately in the implementation of eigs.
    Mitko Georgiev
    @mitkoge

    Thanks,
    you mean calling like this

    λ, v = ApproxFun.eigs_overx(B, -x*D^2 - D + x*fc, 500,tolerance=1E-10)

    and somehow divide C1 to x before calling to

        λ,V = eigen(A1,C1)

    inside the copy of eigs?

    Sheehan Olver
    @dlfivefifty
    No multiply by x: instead of solving A*v = λ*C*v we solve x * A * v = λ * x * C* v. That is, in line 34 change to C = Conversion(ds,rangespace(A)); C = x * C should be enough.
    Note eigs is quite stale and needs to be rewritten, but that's low priority for me. I'll accept any PRs.
    Mitko Georgiev
    @mitkoge
    Thanks for the great support,
    Mitko Georgiev
    @mitkoge

    i followed you and copy paste defined

        function eigs(Bcs_in::Operator,A_in::Operator,fmul::Fun,n::Integer;tolerance::Float64=100eps())

    with the correction

            C = Conversion(ds,rangespace(A)); C = fmul * C

    now have only to check if it works correct. :-)

    λ, v = ApproxFun.eigs(B, L, x, 500,tolerance=1E-10)
    Mitko Georgiev
    @mitkoge
    :point_up: 22. November 2018 12:35
    I intendet to help with replacing it in the docs but as i see it is already fixed.
    Unfortunately it might be not at http://juliaapproximation.github.io/ApproxFun.jl/latest
    Ashton Bradley
    @AshtonSBradley
    is it possible to solve a nonlinear PDE in approxfun? for example I would like to add a simple nonlinear term to this problem
    https://github.com/JuliaApproximation/ApproxFun.jl#solving-partial-differential-equations
    , presumably switching to newton solver, but not sure how to inform it about boundary conditions etc in 2D
    Sheehan Olver
    @dlfivefifty
    Not yet: I want to make the PDE solver faster before trying to implement that sort of thing. For now you can try doing Newton iteration by hand.