Library

Library

Constructing a Fun

GaussWeight(Hermite(L), L) is a space spanned by exp(-Lx²) * H_k(sqrt(L) * x) where H_k(x)'s are Hermite polynomials.

GaussWeight() is equivalent to GaussWeight(Hermite(), 1.0) by default.

Fun(s::Space,coefficients::AbstractVector)

returns a Fun with the specified coefficients in the space s

source
Fun(f,s::Space)

return a Fun representing the function, number, or vector f in the space s. If f is vector-valued, it returns a vector-valued analogue of s.

source
Fun(f,d::Domain)

returns Fun(f,Space(d)), that is, it uses the default space for the specified domain.

source
Fun(s::Space)

returns Fun(identity,s)

source
Fun(f)

returns Fun(f,Chebyshev())

source
Fun()

returns Fun(identity,Chebyshev()).

source
Base.onesMethod.
ones(d::Space)

Return the Fun that represents the function one on the specified space.

source
Base.zerosMethod.
zeros(d::Space)

Return the Fun that represents the function one on the specified space.

source

Domains

Arc(c,r,(θ₁,θ₂))

represents the arc centred at c with radius r from angle θ₁ to θ₂.

Circle(c,r,o)

represents the circle centred at c with radius r which is positively (o=true) or negatively (o=false) oriented.

ApproxFun.CurveConstant.

Curve Represents a domain defined by the image of a Fun. Example usage would be

x=Fun(1..2)
Curve(exp(im*x))  # represents an arc
source
Disk(c,r)

represents the disk centred at c with radius r.

Segment(a,b)

represents a line segment from a to b. In the case where a and b are real and a < b, then this is is equivalent to an Interval(a,b).

An Interval{L,R}(left, right) where L,R are :open or :closed is an interval set containg x such that

  1. left ≤ x ≤ right if L == R == :closed
  2. left < x ≤ right if L == :open and R == :closed
  3. left ≤ x < right if L == :closed and R == :open, or
  4. left < x < right if L == R == :open
Line{a}(c)

represents the line at angle a in the complex plane, centred at c.

PeriodicSegment(a,b)

represents a periodic interval from a to b, that is, the point b is identified with a.

Point(x)

represents a single point at x.

A ProductDomain represents the cartesian product of other domains.

A product domain has two eltypes, an internal type S and an external type T. The internal type S is a tuple containing the eltypes of the elements of the product domain. The external eltype T is a type whose associated space is isomorphic to that of S, but which has been simplified. (See also simplify_product_eltype).

For example, if S is Tuple{Float64,Float64}, then T is SVector{2,Float64}.

Ray{a}(c,o)

represents a ray at angle a starting at c, with orientation out to infinity (o = true) or back from infinity (o = false).

A UnionDomain represents the union of a set of domains.

Missing docstring.

Missing docstring for . Check Documenter's build log for details.

Accessing information about a spaces

canonicalspace(s::Space)

returns a space that is used as a default to implement missing functionality, e.g., evaluation. Implement a Conversion operator or override coefficients to support this.

source
itransform(s::Space,coefficients::AbstractVector)

Transform coefficients back to values. Defaults to using canonicalspace as in transform.

source
transform(s::Space,vals::Vector)

Transform values on the grid specified by points(s,length(vals)) to coefficients in the space s. Defaults to coefficients(transform(canonicalspace(space),values),canonicalspace(space),space)

source
dimension(s::Space)

returns the dimension of s, which is the maximum number of coefficients.

source

Inbuilt spaces

SequenceSpace is the space of all sequences, i.e., infinite vectors. Also denoted ℓ⁰.

ConstantSpace is the 1-dimensional scalar space.

Chebyshev() is the space spanned by the Chebyshev polynomials

    T_0(x),T_1(x),T_2(x),…

where T_k(x) = cos(k*acos(x)). This is the default space as there exists a fast transform and general smooth functions on [-1,1] can be easily resolved.

Hemite(L) represents H_k(sqrt(L) * x) where H_k are Hermite polynomials. Hermite() is equivalent to Hermite(1.0).

Jacobi(b,a) represents the space spanned by Jacobi polynomials P_k^{(a,b)}, which are orthogonal with respect to the weight (1+x)^β*(1-x)^α

Laguerre(α) is a space spanned by generalized Laguerre polynomials Lₙᵅ(x) 's on (0, Inf), which satisfy the differential equations

    xy'' + (α + 1 - x)y' + ny = 0

Laguerre() is equivalent to Laguerre(0) by default.

Ultraspherical(λ) is the space spanned by the ultraspherical polynomials

    C_0^{(λ)}(x),C_1^{(λ)}(x),C_2^{(λ)}(x),…

Note that λ=1 this reduces to Chebyshev polynomials of the second kind: C_k^{(1)}(x) = U_k(x). For λ=1/2 this also reduces to Legendre polynomials: C_k^{(1/2)}(x) = P_k(x).

Taylor() is the space spanned by [1,z,z^2,...]. This is a type alias for Hardy{true}.

Hardy{false}() is the space spanned by [1/z,1/z^2,...]. Hardy{true}() is the space spanned by [1,z,z^2,...].

Fourier() is the space spanned by the trigonemtric polynomials

    1,sin(θ),cos(θ),sin(2θ),cos(2θ),…

See also Laurent.

Laurent() is the space spanned by the complex exponentials

    1,exp(-im*θ),exp(im*θ),exp(-2im*θ),…

See also Fourier.

CosSpace() is the space spanned by [1,cos θ,cos 2θ,...]

SinSpace() is the space spanned by [sin θ,sin 2θ,...]

JacobiWeight(β,α,s::Space)

weights a space s by a Jacobi weight, which on -1..1 is (1+x)^β*(1-x)^α. For other domains, the weight is inferred by mapping to -1..1.

LogWeight(β,α,s::Space)

represents a function on -1..1 weighted by log((1+x)^β*(1-x)^α). For other domains, the weight is inferred by mapping to -1..1.

ArraySpace(s::Space,dims...)

is used to represent array-valued expansions in a space s. The coefficients are of each entry are interlaced.

For example,

f = Fun(x->[exp(x),sin(x)],-1..1)
space(f) == ArraySpace(Chebyshev(),2)
TensorSpace(a::Space,b::Space)

represents a tensor product of two 1D spaces a and b. The coefficients are interlaced in lexigraphical order.

For example, consider

Fourier()*Chebyshev()  # returns TensorSpace(Fourier(),Chebyshev())

This represents functions on [-π,π) x [-1,1], using the Fourier basis for the first argument and Chebyshev basis for the second argument, that is, φ_k(x)T_j(y), where

φ_0(x) = 1,
φ_1(x) = sin x,
φ_2(x) = cos x,
φ_3(x) = sin 2x,
φ_4(x) = cos 2x
…

By Choosing (k,j) appropriately, we obtain a single basis:

φ_0(x)T_0(y) (= 1),
φ_0(x)T_1(y) (= y),
φ_1(x)T_0(y) (= sin x),
φ_0(x)T_2(y), …

Accessing information about a Fun

ApproxFunBase.domainFunction.
domain(f::Fun)

returns the domain that f is defined on

source
coefficients(f::Fun) -> Vector

returns the coefficients of f, corresponding to the space space(f).

source
coefficients(f::Fun,s::Space) -> Vector

returns the coefficients of f in the space s, which may not be the same as space(f).

source
coefficients(cfs::AbstractVector,fromspace::Space,tospace::Space) -> Vector

converts coefficients in fromspace to coefficients in tospace

source
extrapolate(f::Fun,x)

returns an extrapolation of f from its domain to x.

source
ncoefficients(f::Fun) -> Integer

returns the number of coefficients of a fun

source
ApproxFunBase.pointsFunction.
points(f::Fun)

returns a grid of points that f can be transformed into values and back

source
points(s::Space,n::Integer)

returns a grid of approximately n points, for which a transform exists from values at the grid to coefficients in the space s.

source
ApproxFunBase.spaceFunction.
space(f::Fun)

returns the space of f

source
Base.valuesFunction.
values(iterator)

For an iterator or collection that has keys and values, return an iterator over the values. This function simply returns its argument by default, since the elements of a general iterator are normally considered its "values".

Examples

julia> d = Dict("a"=>1, "b"=>2);

julia> values(d)
Base.ValueIterator for a Dict{String,Int64} with 2 entries. Values:
  2
  1

julia> values([2])
1-element Array{Int64,1}:
 2
source
values(a::AbstractDict)

Return an iterator over all values in a collection. collect(values(a)) returns an array of values. Since the values are stored internally in a hash table, the order in which they are returned may vary. But keys(a) and values(a) both iterate a and return the elements in the same order.

Examples

julia> D = Dict('a'=>2, 'b'=>3)
Dict{Char,Int64} with 2 entries:
  'a' => 2
  'b' => 3

julia> collect(values(D))
2-element Array{Int64,1}:
 2
 3
source
values(f::Fun)

returns f evaluated at points(f)

source
Base.strideMethod.
stride(f::Fun)

returns the stride of the coefficients, checked numerically

source

Modify a Fun

reverseorientation(f::Fun)

return f on a reversed orientated contour.

source
setdomain(f::Fun,d::Domain)

returns f projected onto domain

source
Base.chopFunction.
chop(s::AbstractString; head::Integer = 0, tail::Integer = 1)

Remove the first head and the last tail characters from s. The call chop(s) removes the last character from s. If it is requested to remove more characters than length(s) then an empty string is returned.

Examples

julia> a = "March"
"March"

julia> chop(a)
"Marc"

julia> chop(a, head = 1, tail = 2)
"ar"

julia> chop(a, head = 5, tail = 5)
""
source
chop(p::Poly{T}; rtol::Real = Base.rtoldefault(real(T)), atol::Real = 0)

Chop off leading values from a polynomial which are approximately zero. The tolerances rtol and atol are passed to isapprox to check for zeros.

chop(f::Fun,tol) -> Fun

reduces the number of coefficients by dropping the tail that is below the specified tolerance.

source

Operators

Operator{T}

is an abstract type to represent linear operators between spaces.

source
bandwidths(op::Operator)

returns the bandwidth of op in the form (l,u), where l ≥ 0 represents the number of subdiagonals and u ≥ 0 represents the number of superdiagonals.

source
domainspace(op::Operator)

gives the domain space of op. That is, op*f will first convert f to a Fun in the space domainspace(op) before applying the operator.

source
rangespace(op::Operator)

gives the range space of op. That is, op*f will return a Fun in the space rangespace(op), provided f can be converted to a Fun in domainspace(op).

source
Base.getindexMethod.
op[f::Fun]

constructs the operator op*Multiplication(f), that is, it multiplies on the right by f first. Note that op*f is different: it applies op to f.

source
LinearAlgebra.qrMethod.
qr(A::Operator)

returns a cached QR factorization of the Operator A. The result QR enables solving of linear equations: if u=QR, then u approximately satisfies A*u = b.

LazyArrays.cacheMethod.
cache(op::Operator)

Caches the entries of an operator, to speed up multiplying a Fun by the operator.

Inbuilt operators

Conversion(fromspace::Space,tospace::Space)

represents a conversion operator between fromspace and tospace, when available.

source

Derivative(sp::Space,k::Int) represents the k-th derivative on sp.

Derivative(sp::Space,k::Vector{Int}) represents a partial derivative on a multivariate space. For example,

Dx = Derivative(Chebyshev()^2,[1,0]) # ∂/∂x
Dy = Derivative(Chebyshev()^2,[0,1]) # ∂/∂y

Derivative(sp::Space) represents the first derivative Derivative(sp,1).

Derivative(k) represents the k-th derivative, acting on an unset space. Spaces will be inferred when applying or manipulating the operator.

Derivative() represents the first derivative on an unset space. Spaces will be inferred when applying or manipulating the operator.

Dirichlet(sp,k) is the operator associated with restricting the k-th derivative on the boundary for the space sp.

Dirichlet(sp) is the operator associated with restricting the the boundary for the space sp.

Dirichlet() is the operator associated with restricting on the the boundary.

Evaluation(sp,x,k) is the functional associated with evaluating the k-th derivative at a point x for the space sp.

Evaluation(sp,x) is the functional associated with evaluating at a point x for the space sp.

Evaluation(x) is the functional associated with evaluating at a point x.

Integral(sp::Space,k::Int) represents a k-th integral on sp. There is no guarantee on normalization.

Integral(sp::Space) represents the first integral Integral(sp,1).

Integral(k)represents thek`-th integral, acting on an unset space. Spaces will be inferred when applying or manipulating the operator.

Intergral() represents the first integral on an unset space. Spaces will be inferred when applying or manipulating the operator.

Laplacian(sp::Space) represents the laplacian on space sp.

Laplacian() represents the laplacian on an unset space. Spaces will be inferred when applying or manipulating the operator.

Multiplication(f::Fun,sp::Space) is the operator representing multiplication by f on functions in the space sp.

Multiplication(f::Fun) is the operator representing multiplication by f on an unset space of functions. Spaces will be inferred when applying or manipulating the operator.

ApproxFunBase.NeumannFunction.

Neumann(sp) is the operator associated with restricting the normal derivative on the boundary for the space sp. At the moment it is implemented as Dirichlet(sp,1).

`Neumann( is the operator associated with restricting the normal derivative on the boundary.

Bivariate

LowRankFun gives an approximation to a bivariate function in low rank form.