JuliaApproximation / ApproxFun.jl

Julia package for function approximation
http://juliaapproximation.github.io/ApproxFun.jl/
Other
539 stars 70 forks source link

Support (f::Fun)(::Domain) for restricting a function to a sub domain #414

Open BoundaryValueProblems opened 7 years ago

BoundaryValueProblems commented 7 years ago

Hello. I want to solve the Helmholtz equation on the square [-1,1] x [-1,1] with the Dirichlet or Neumann boundary conditions where the boundary values are prescribed by the boundary value of a given 2D function. For example, let's use a Gabor function of the form:

julia> d = Interval()^2
julia> g=Fun((x,y)->exp(-(x^2+3*y^2)*cos(5*pi*x+7*pi*y),d)

Now, I want to use the boundary values of this function g to solve the Helmholtz equation. Following the example in README.md, I tried:

julia> Δ = Laplacian(d)                            # Represent the Laplacian
julia> f = g(∂(d))                                 # f should represent the boundary value of g
julia> u = [dirichlet(d);Δ+500I]\[f;0.]            # Solve the PDE

However, it generates the following error at the 2nd line. What should I do? Thanks! PS: By the way, there is a typo in the Helmholtz solver example in README.md. The third line above in README.md says:

julia> u = [Dirichlet(d); Δ+500I]\[f;0.]            # Solve the PDE

But Dirichlet should be replaced by dirichlet; otherwise it causes an error.

Here is the error generated by:

julia> f = g(∂(d))                                 # f should represent the boundary value of g

MethodError: no method matching start(::ApproxFun.PiecewiseInterval{Complex{Float64}}) Closest candidates are: start(!Matched::SimpleVector) at essentials.jl:170 start(!Matched::Base.MethodList) at reflection.jl:258 start(!Matched::IntSet) at intset.jl:184 ... in append_any(::ApproxFun.PiecewiseInterval{Complex{Float64}}, ::Vararg{ApproxFun.PiecewiseInterval{Complex{Float64}},N}) at ./essentials.jl:98 in evaluate(::ApproxFun.Fun{ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2},Float64}, ::ApproxFun.PiecewiseInterval{Complex{Float64}}) at /Users/xxx/.julia/v0.5/ApproxFun/src/Fun/Fun.jl:141 in (::ApproxFun.Fun{ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2},Float64})(::ApproxFun.PiecewiseInterval{Complex{Float64}}, ::Vararg{ApproxFun.PiecewiseInterval{Complex{Float64}},N}) at...

dlfivefifty commented 7 years ago

The README is for the development version, which hasn't been tagged yet.

The syntax f = g(∂(d)) doesn't work (though I guess it should!). Use instead

Fun((x,y)->g(x,y),∂(d))
BoundaryValueProblems commented 7 years ago

Thanks, dlfivefifty! However, I tried the following two possibilities, and still got errors.

julia> f=Fun((x,y)->g(x,y),∂(d))
julia> f=Fun((x,y)->exp(-(x^2+3*y^2))*cos(5*pi*x+7*pi*y),∂(d))

The error messages in both cases are: v MethodError: no method matching (::##3#4)(::Complex{Float64}) Closest candidates are:

3(::Any, !Matched::Any) at console:1

in zerocfsFun at ApproxFun/src/Fun/constructors.jl:134
in #Fun#45 at ApproxFun/src/Fun/constructors.jl:206
in  at base/<missing>
in #Fun#45 at ApproxFun/src/Fun/constructors.jl:202
in #Fun#48 at ApproxFun/src/Fun/constructors.jl:213
in ApproxFun.Fun{S,T} at ApproxFun/src/Fun/constructors.jl:213

Thanks for your help!

dlfivefifty commented 7 years ago

What version of Julia are you on?

You can you also try on development. The tagged version of PDE solving is less robust. (There was a major revamp.)

Pkg.checkout(“ApproxFun”,”development”)
d=Interval()^2
f=Fun((x,y)->exp(-(x^2+3*y^2))*cos(5*pi*x+7*pi*y),∂(d))
u = linsolve([Dirichlet(d); Δ+500I],[f;0.];tolerance=1E-7)

Note that it is necessary to specify the tolerance as otherwise the solution time will be too long (the solver is still being optimized.)

dlfivefifty commented 7 years ago

If you require multiple RHS, you can do the following the make subsequent solves faster:

d=Interval()^2
f=Fun((x,y)->exp(-(x^2+3*y^2))*cos(5*pi*x+7*pi*y),∂(d))
QR=qrfact([Dirichlet(d); Δ+500I])
@time u = linsolve(QR,[f;0.];tolerance=1E-7)   # 2.392173s
@time u = linsolve(QR,[f;0.];tolerance=1E-7)   # 0.009747s
BoundaryValueProblems commented 7 years ago

The julia version I have is: 0.5.0. As you suggested, I just did:

Pkg.checkout("ApproxFun", "development")
d=Interval()^2
f=Fun((x,y)->exp(-(x^2+3*y^2))*cos(5*pi*x+7*pi*y),∂(d))

Then, I got the same error as I reported earlier. Do I need to use the development version of julia, i.e., 0.6.x, in order for these to work? If possible, I want to stick with 0.5.0 if possible. Thanks!

dlfivefifty commented 7 years ago

No this is on 0.6. 0.5 has some outstanding bugs with conflicting packages which I suspect you are running into. You would have to delete your .julia folder to fix it.

On 4 Nov. 2016, at 8:32 am, BoundaryValueProblems notifications@github.com wrote:

The julia version I have is: 0.5.0. As you suggested, I just did:

Pkg.checkout("ApproxFun", "development") d=Interval()^2 f=Fun((x,y)->exp(-(x^2+3_y^2))_cos(5_pi_x+7_pi_y),∂(d)) Then, I got the same error as I reported earlier. Do I need to use the development version of julia, i.e., 0.6.x, in order for these to work? If possible, I want to stick with 0.5.0 if possible. Thanks!

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/ApproxFun/ApproxFun.jl/issues/414#issuecomment-258280603, or mute the thread https://github.com/notifications/unsubscribe-auth/ABLDqfgiNHaw_UIkQFvJeBnU99drsleMks5q6lL8gaJpZM4Ko5JY.

BoundaryValueProblems commented 7 years ago

Well, I just installed julia version 0.6.0-dev.1178, and did exactly what you suggested including installing the developmental version of ApproxFun. Then, the following line generated tons of errors:

f=Fun((x,y)->exp(-(x^2+3*y^2))*cos(5*pi*x+7*pi*y),∂(d))

ERROR: MethodError: Cannot convert an object of type Int64 to an object of type ApproxFun.Infinity{Bool} This may have arisen from a call to the constructor ApproxFun.Infinity{Bool}(...), since type constructors fall back to convert methods. in mapfoldl_impl(::ApproxFun.#dimension, ::Base.#*, ::ApproxFun.Infinity{Bool}, ::Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}}, ::Int64) at ./reduce.jl:43 in mapfoldl(::ApproxFun.#dimension, ::Function, ::Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}}) at ./reduce.jl:73 in #Fun#81(::String, ::Type{T}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:207 in (::Core.#kw#Type)(::Array{Any,1}, ::Type{ApproxFun.Fun}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at ./:0 in #Fun#81(::String, ::Type{T}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:206 in #Fun#84(::Array{Any,1}, ::Type{T}, ::Function, ::ApproxFun.ProductDomain{Tuple{ApproxFun.Interval{Float64},ApproxFun.Interval{Float64}},Float64,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:217 in ApproxFun.Fun{S,T}(::Function, ::ApproxFun.ProductDomain{Tuple{ApproxFun.Interval{Float64},ApproxFun.Interval{Float64}},Float64,2}) at /Users/xx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:217

dlfivefifty commented 7 years ago

0.6 is not supported

On 4 Nov. 2016, at 9:13 am, BoundaryValueProblems notifications@github.com wrote:

Well, I just installed julia version 0.6.0-dev.1178, and did exactly what you suggested including installing the developmental version of ApproxFun. Then, the following line generated tons of errors:

f=Fun((x,y)->exp(-(x^2+3_y^2))_cos(5_pi_x+7_pi_y),∂(d)) ERROR: MethodError: Cannot convert an object of type Int64 to an object of type ApproxFun.Infinity{Bool} This may have arisen from a call to the constructor ApproxFun.Infinity{Bool}(...), since type constructors fall back to convert methods. in mapfoldl_impl(::ApproxFun.#dimension, ::Base.#*, ::ApproxFun.Infinity{Bool}, ::Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}}, ::Int64) at ./reduce.jl:43 in mapfoldl(::ApproxFun.#dimension, ::Function, ::Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}}) at ./reduce.jl:73 in #Fun#81(::String, ::Type{T}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:207 in (::Core.#kw#Type)(::Array{Any,1}, ::Type{ApproxFun.Fun}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at ./:0 in #Fun#81(::String, ::Type{T}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:206 in #Fun#84(::Array{Any,1}, ::Type{T}, ::Function, ::ApproxFun.ProductDomain{Tuple{ApproxFun.Interval{Float64},ApproxFun.Interval{Float64}},Float64,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:217 in ApproxFun.Fun{S,T}(::Function, ::ApproxFun.ProductDomain{Tuple{ApproxFun.Interval{Float64},ApproxFun.Interval{Float64}},Float64,2}) at /Users/xx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:217

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/ApproxFun/ApproxFun.jl/issues/414#issuecomment-258290188, or mute the thread https://github.com/notifications/unsubscribe-auth/ABLDqaTIuLjRDUMZkRRTKfcxNuQyRZxfks5q6lyWgaJpZM4Ko5JY.

dlfivefifty commented 7 years ago

Sorry, I made a typo. I meant to say “No, this is on 0.5"

On 4 Nov. 2016, at 9:14 am, Sheehan Olver solver@mac.com wrote:

0.6 is not supported

On 4 Nov. 2016, at 9:13 am, BoundaryValueProblems <notifications@github.com mailto:notifications@github.com> wrote:

Well, I just installed julia version 0.6.0-dev.1178, and did exactly what you suggested including installing the developmental version of ApproxFun. Then, the following line generated tons of errors:

f=Fun((x,y)->exp(-(x^2+3_y^2))_cos(5_pi_x+7_pi_y),∂(d)) ERROR: MethodError: Cannot convert an object of type Int64 to an object of type ApproxFun.Infinity{Bool} This may have arisen from a call to the constructor ApproxFun.Infinity{Bool}(...), since type constructors fall back to convert methods. in mapfoldl_impl(::ApproxFun.#dimension, ::Base.#*, ::ApproxFun.Infinity{Bool}, ::Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}}, ::Int64) at ./reduce.jl:43 in mapfoldl(::ApproxFun.#dimension, ::Function, ::Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}}) at ./reduce.jl:73 in #Fun#81(::String, ::Type{T}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:207 in (::Core.#kw#Type)(::Array{Any,1}, ::Type{ApproxFun.Fun}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at ./:0 in #Fun#81(::String, ::Type{T}, ::Function, ::ApproxFun.TensorSpace{Tuple{ApproxFun.Chebyshev{ApproxFun.Interval{Float64}},ApproxFun.Chebyshev{ApproxFun.Interval{Float64}}},ApproxFun.RealBasis,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:206 in #Fun#84(::Array{Any,1}, ::Type{T}, ::Function, ::ApproxFun.ProductDomain{Tuple{ApproxFun.Interval{Float64},ApproxFun.Interval{Float64}},Float64,2}) at /Users/xxx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:217 in ApproxFun.Fun{S,T}(::Function, ::ApproxFun.ProductDomain{Tuple{ApproxFun.Interval{Float64},ApproxFun.Interval{Float64}},Float64,2}) at /Users/xx/.julia/v0.6/ApproxFun/src/Fun/constructors.jl:217

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/ApproxFun/ApproxFun.jl/issues/414#issuecomment-258290188, or mute the thread https://github.com/notifications/unsubscribe-auth/ABLDqaTIuLjRDUMZkRRTKfcxNuQyRZxfks5q6lyWgaJpZM4Ko5JY.

dlfivefifty commented 7 years ago

OK, I remembered in the last tagged version it uses complex numbers instead of x, y. So try this:

f=Fun(z->(x=real(z);y=imag(z);exp(-(x^2+3*y^2))*cos(5*pi*x+7*pi*y)),∂(d))
Δ = Laplacian(d) 
u = pdesolve([dirichlet(d);Δ+500I],[f;0.],200)
dlfivefifty commented 7 years ago

Or, on development in Julia 0.5, you can try moving your .julia folder aside.

dlfivefifty commented 7 years ago

Sorry for the typo....really annoying one.

BoundaryValueProblems commented 7 years ago

Well, before I installed 0.6, I tried the developmental version of ApproxFun as you suggested, and as I reported earlier, I got the error. Please clarify what setting I should use: (1) Julia 0.5.0; tagged version of ApproxFun.jl; with complex numbers as you just suggested or (2) Julia 0.5.0; development version of ApproxFun.jl; without using complex numbers as you suggested. In either case, should I rename my ~/.julia and rebuild a new ~/.julia? I would greatly appreciate if you could clarify these once and for all. Thanks!

dlfivefifty commented 7 years ago

Yes. In case (1) you might be able to get away with not moving away .julia.

There's been several people who've reported weird bugs where moving away .julia fixes it. I believe it's due to the open issue https://github.com/JuliaLang/julia/issues/18465 but that's just a guess.

dlfivefifty commented 7 years ago

Note that (1) and (2) have very different solvers. The solver in (1) is faster if you are changing k. The solver in (2) is faster if yo are changing the RHS.

The solver in (1) has since been moved to https://github.com/ApproxFun/PDESchurFactorization.jl but has not been fixed up for the current development version.

dlfivefifty commented 7 years ago

Ah, for (2) you might have to also do Pkg.checkout("BandedMatrices"). I'll tag that now to remove one step.

BoundaryValueProblems commented 7 years ago

Thanks a lot! I tried (1) on OS X and (2) on Windows, and both worked! Please let me know when you plan to tag the current development version of ApproxFun so that I do not need to use two different ways to specify the boundary values from a given 2D function. Thank you very much for your effort of creating extremely useful package!!

dlfivefifty commented 7 years ago

I'll try to tag it by next week. I was hoping to speed up the solution to Fourier PDEs first, but I think that will have to be postponed.

BoundaryValueProblems commented 7 years ago

Thanks! When you commit the tagged version next week, what I should do with the development version I just installed on OS X via Pkg.checkout("ApproxFun", "development")? Can I just do Pkg.update() at that point and don't need to do anything else? Or do I have to remove the development version of ApproxFun first, then upgrade it via Pkg.update()?

dlfivefifty commented 7 years ago
help?> Pkg.free
  free(pkg)

  Free the package pkg to be managed by the package manager again. It calls
  Pkg.resolve() to determine optimal package versions after. This is an
  inverse for both Pkg.checkout and Pkg.pin.

  You can also supply an iterable collection of package names, e.g.,
  Pkg.free(("Pkg1", "Pkg2")) to free multiple packages at once.