Skip to content

Commit 027c436

Browse files
dlfivefiftydlfivefifty
committed
Move out singularities
1 parent ac8e454 commit 027c436

26 files changed

+3
-2195
lines changed

‎src/Spaces/IntervalSpace.jl

Lines changed: 0 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -16,16 +16,3 @@ Fun(::typeof(identity), d::IntervalOrSegment{T}) where {T<:Number} =
1616
# the default domain space is higher to avoid negative ultraspherical spaces
1717
Integral(d::IntervalOrSegment,n::Integer) = Integral(Ultraspherical(1,d),n)
1818

19-
for Func in (:DefiniteIntegral,:DefiniteLineIntegral)
20-
@eval begin
21-
#TODO: this may be misleading
22-
$Func(d::IntervalOrSegment) = $Func(JacobiWeight(-0.5,-0.5,Chebyshev(d)))
23-
function $Func(α::Number,β::Number,d::IntervalOrSegment)
24-
@assert α == β
25-
@assert round(Int,α+.5) == α+.5
26-
@assert round(Int,α+.5) >= 0
27-
$Func(JacobiWeight(α,β,Ultraspherical(round(Int,α+.5),d)))
28-
end
29-
$Func(α::Number,β::Number) = $Func(α,β,ChebyshevInterval())
30-
end
31-
end

‎src/Spaces/Jacobi/Jacobi.jl

Lines changed: 1 addition & 60 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
export Jacobi, NormalizedJacobi, Legendre, NormalizedLegendre, WeightedJacobi
1+
export Jacobi, NormalizedJacobi, Legendre, NormalizedLegendre
22

33

44
"""
@@ -41,10 +41,6 @@ Base.promote_rule(::Type{Jacobi{D,R1}},::Type{Jacobi{D,R2}}) where {D,R1,R2} =
4141
convert(::Type{Jacobi{D,R1}},J::Jacobi{D,R2}) where {D,R1,R2} =
4242
Jacobi{D,R1}(J.b,J.a,J.domain)
4343

44-
const WeightedJacobi{D,R} = JacobiWeight{Jacobi{D,R},D,R,R}
45-
46-
WeightedJacobi(β,α,d::Domain) = JacobiWeight(β,α,Jacobi(β,α,d))
47-
WeightedJacobi(β,α) = JacobiWeight(β,α,Jacobi(β,α))
4844

4945
spacescompatible(a::Jacobi,b::Jacobi) = a.a ≈ b.a && a.b ≈ b.b && domainscompatible(a,b)
5046

@@ -179,34 +175,6 @@ function bilinearform(f::Fun{J},g::Fun{J}) where {J<:Jacobi}
179175
end
180176
end
181177

182-
function bilinearform(f::Fun{JacobiWeight{J,DD,RR,TT}},g::Fun{J}) where {J<:Jacobi,DD<:IntervalOrSegment,RR,TT}
183-
@assert domain(f) == domain(g)
184-
if f.space.β == f.space.space.a == g.space.a && f.space.α == f.space.space.b == g.space.b
185-
return complexlength(domain(f))/2*conjugatedinnerproduct(g.space,f.coefficients,g.coefficients)
186-
else
187-
return defaultbilinearform(f,g)
188-
end
189-
end
190-
191-
function bilinearform(f::Fun{J},
192-
g::Fun{JacobiWeight{J,DD,RR,TT}}) where {J<:Jacobi,DD<:IntervalOrSegment,RR,TT}
193-
@assert domain(f) == domain(g)
194-
if g.space.β == g.space.space.a == f.space.a && g.space.α == g.space.space.b == f.space.b
195-
return complexlength(domain(f))/2*conjugatedinnerproduct(f.space,f.coefficients,g.coefficients)
196-
else
197-
return defaultbilinearform(f,g)
198-
end
199-
end
200-
201-
function bilinearform(f::Fun{JacobiWeight{J,DD,RR,TT}},
202-
g::Fun{JacobiWeight{J,DD,RR,TT}}) where {J<:Jacobi,DD<:IntervalOrSegment,RR,TT}
203-
@assert domain(f) == domain(g)
204-
if f.space.β + g.space.β == f.space.space.a == g.space.space.a && f.space.α + g.space.α == f.space.space.b == g.space.space.b
205-
return complexlength(domain(f))/2*conjugatedinnerproduct(f.space.space,f.coefficients,g.coefficients)
206-
else
207-
return defaultbilinearform(f,g)
208-
end
209-
end
210178

211179
function linebilinearform(f::Fun{J},g::Fun{J}) where {J<:Jacobi}
212180
@assert domain(f) == domain(g)
@@ -216,30 +184,3 @@ function linebilinearform(f::Fun{J},g::Fun{J}) where {J<:Jacobi}
216184
return defaultlinebilinearform(f,g)
217185
end
218186
end
219-
220-
function linebilinearform(f::Fun{JacobiWeight{J,DD,RR,TT}},g::Fun{J}) where {J<:Jacobi,DD<:IntervalOrSegment,RR,TT}
221-
@assert domain(f) == domain(g)
222-
if f.space.β == f.space.space.a == g.space.a && f.space.α == f.space.space.b == g.space.b
223-
return arclength(domain(f))/2*conjugatedinnerproduct(g.space,f.coefficients,g.coefficients)
224-
else
225-
return defaultlinebilinearform(f,g)
226-
end
227-
end
228-
229-
function linebilinearform(f::Fun{J},g::Fun{JacobiWeight{J,DD,RR,TT}}) where {J<:Jacobi,DD<:IntervalOrSegment,RR,TT}
230-
@assert domain(f) == domain(g)
231-
if g.space.β == g.space.space.a == f.space.a && g.space.α == g.space.space.b == f.space.b
232-
return arclength(domain(f))/2*conjugatedinnerproduct(f.space,f.coefficients,g.coefficients)
233-
else
234-
return defaultlinebilinearform(f,g)
235-
end
236-
end
237-
238-
function linebilinearform(f::Fun{JacobiWeight{J,DD,RR,TT}}, g::Fun{JacobiWeight{J,DD,RR,TT}}) where {J<:Jacobi,DD<:IntervalOrSegment,RR,TT}
239-
@assert domain(f) == domain(g)
240-
if f.space.β + g.space.β == f.space.space.a == g.space.space.a && f.space.α + g.space.α == f.space.space.b == g.space.space.b
241-
return arclength(domain(f))/2*conjugatedinnerproduct(f.space.space,f.coefficients,g.coefficients)
242-
else
243-
return defaultlinebilinearform(f,g)
244-
end
245-
end

‎src/Spaces/Jacobi/JacobiOperators.jl

Lines changed: 0 additions & 125 deletions
Original file line numberDiff line numberDiff line change
@@ -108,22 +108,6 @@ getindex(T::ConcreteDerivative{J},k::Integer,j::Integer) where {J<:Jacobi} =
108108

109109

110110

111-
function Derivative(S::WeightedJacobi{DDD,RR}) where {DDD<:IntervalOrSegment,RR}
112-
if S.β>0 && S.β>0 && S.β==S.space.b && S.α==S.space.a
113-
ConcreteDerivative(S,1)
114-
else
115-
jacobiweightDerivative(S)
116-
end
117-
end
118-
119-
bandwidths(D::ConcreteDerivative{WeightedJacobi{DDD,RR}}) where {DDD<:IntervalOrSegment,RR} = 1,0
120-
rangespace(D::ConcreteDerivative{WeightedJacobi{DDD,RR}}) where {DDD<:IntervalOrSegment,RR} =
121-
WeightedJacobi(domainspace(D).β-1,domainspace(D).α-1,domain(D))
122-
123-
getindex(D::ConcreteDerivative{WeightedJacobi{DDD,RR}},k::Integer,j::Integer) where {DDD<:IntervalOrSegment,RR} =
124-
j==k-1 ? eltype(D)(-4(k-1)./complexlength(domain(D))) : zero(eltype(D))
125-
126-
127111

128112
## Integral
129113

@@ -201,27 +185,6 @@ for (Func,Len,Sum) in ((:DefiniteIntegral,:complexlength,:sum),(:DefiniteLineInt
201185
end
202186
end
203187

204-
205-
function getindex(Σ::$ConcFunc{JacobiWeight{Jacobi{D,R},D,R,TT},T},k::Integer) where {D<:IntervalOrSegment,R,T,TT}
206-
dsp = domainspace(Σ)
207-
208-
if dsp.β == dsp.space.b && dsp.α == dsp.space.a
209-
# TODO: copy and paste
210-
k == 1 ? convert(T,$Sum(Fun(dsp,[one(T)]))) : zero(T)
211-
else
212-
convert(T,$Sum(Fun(dsp,[zeros(T,k-1);1])))
213-
end
214-
end
215-
216-
function bandwidths(Σ::$ConcFunc{JacobiWeight{Jacobi{D,R},D,R,TT}}) where {D<:IntervalOrSegment,R,TT}
217-
β,α = domainspace(Σ).β,domainspace(Σ).α
218-
if domainspace(Σ).β == domainspace(Σ).space.b && domainspace(Σ).α == domainspace(Σ).space.a
219-
0,0 # first entry
220-
else
221-
0,∞
222-
end
223-
end
224-
225188
function bandwidths(Σ::$ConcFunc{Jacobi{D,R}}) where {D<:IntervalOrSegment,R}
226189
if domainspace(Σ).b == domainspace(Σ).a == 0
227190
0,0 # first entry
@@ -581,94 +544,6 @@ hasconversion(a::Ultraspherical,b::Jacobi) = hasconversion(Jacobi(a),b)
581544
# (1+x) or (1-x) by _decreasing_ the parameter. Thus the
582545

583546

584-
## <: IntervalOrSegment avoids a julia bug
585-
function Multiplication(f::Fun{JacobiWeight{C,DD,RR,TT}}, S::Jacobi) where {C<:ConstantSpace,DD<:IntervalOrSegmentDomain,RR,TT}
586-
# this implements (1+x)*P and (1-x)*P special case
587-
# see DLMF (18.9.6)
588-
d=domain(f)
589-
if ((space(f).β==1 && space(f).α==0 && S.b >0) ||
590-
(space(f).β==0 && space(f).α==1 && S.a >0))
591-
ConcreteMultiplication(f,S)
592-
elseif isapproxinteger(space(f).β) && space(f).β ≥ 1 && S.b >0
593-
# decrement β and multiply again
594-
M=Multiplication(f.coefficients[1]*jacobiweight(1.,0.,d),S)
595-
MultiplicationWrapper(f,Multiplication(jacobiweight(space(f).β-1,space(f).α,d),rangespace(M))*M)
596-
elseif isapproxinteger(space(f).α) && space(f).α ≥ 1 && S.a >0
597-
# decrement α and multiply again
598-
M=Multiplication(f.coefficients[1]*jacobiweight(0.,1.,d),S)
599-
MultiplicationWrapper(f,Multiplication(jacobiweight(space(f).β,space(f).α-1,d),rangespace(M))*M)
600-
else
601-
# default JacobiWeight
602-
M=Multiplication(Fun(space(f).space,f.coefficients),S)
603-
rsp=JacobiWeight(space(f).β,space(f).α,rangespace(M))
604-
MultiplicationWrapper(f,SpaceOperator(M,S,rsp))
605-
end
606-
end
607-
608-
Multiplication(f::Fun{JacobiWeight{C,DD,RR,TT}},
609-
S::Union{Ultraspherical,Chebyshev}) where {C<:ConstantSpace,DD<:IntervalOrSegmentDomain,RR,TT} =
610-
MultiplicationWrapper(f,Multiplication(f,Jacobi(S))*Conversion(S,Jacobi(S)))
611-
612-
function rangespace(M::ConcreteMultiplication{JacobiWeight{C,DD,RR,TT},J}) where {J<:Jacobi,C<:ConstantSpace,DD<:IntervalOrSegmentDomain,RR,TT}
613-
S=domainspace(M)
614-
if space(M.f).β==1
615-
# multiply by (1+x)
616-
Jacobi(S.b-1,S.a,domain(S))
617-
elseif space(M.f).α == 1
618-
# multiply by (1-x)
619-
Jacobi(S.b,S.a-1,domain(S))
620-
else
621-
error("Not implemented")
622-
end
623-
end
624-
625-
bandwidths(::ConcreteMultiplication{JacobiWeight{C,DD,RR,TT},J}) where {J<:Jacobi,C<:ConstantSpace,DD<:IntervalOrSegmentDomain,RR,TT} = 1,0
626-
627-
628-
function getindex(M::ConcreteMultiplication{JacobiWeight{C,DD,RR,TT},J},k::Integer,j::Integer) where {J<:Jacobi,C<:ConstantSpace,DD<:IntervalOrSegmentDomain,RR,TT}
629-
@assert ncoefficients(M.f)==1
630-
a,b=domainspace(M).a,domainspace(M).b
631-
c=M.f.coefficients[1]
632-
if space(M.f).β==1
633-
@assert space(M.f).α==0
634-
# multiply by (1+x)
635-
if j==k
636-
c*2(k+b-1)/(2k+a+b-1)
637-
elseif k > 1 && j==k-1
638-
c*(2k-2)/(2k+a+b-3)
639-
else
640-
zero(eltype(M))
641-
end
642-
elseif space(M.f).α == 1
643-
@assert space(M.f).β==0
644-
# multiply by (1-x)
645-
if j==k
646-
c*2(k+a-1)/(2k+a+b-1)
647-
elseif k > 1 && j==k-1
648-
-c*(2k-2)/(2k+a+b-3)
649-
else
650-
zero(eltype(M))
651-
end
652-
else
653-
error("Not implemented")
654-
end
655-
end
656-
657-
658-
# We can exploit the special multiplication to construct a Conversion
659-
660-
661-
for FUNC in (:maxspace_rule,:union_rule,:hasconversion)
662-
@eval function $FUNC(A::WeightedJacobi{DD},B::Jacobi) where DD<:IntervalOrSegment
663-
if A.β==A.α+1 && A.space.b>0
664-
$FUNC(Jacobi(A.space.b-1,A.space.a,domain(A)),B)
665-
elseif A.α==A.β+1 && A.space.a>0
666-
$FUNC(Jacobi(A.space.b,A.space.a-1,domain(A)),B)
667-
else
668-
$FUNC(A,JacobiWeight(0.,0.,B))
669-
end
670-
end
671-
end
672547

673548

674549

‎src/Spaces/Laguerre/Laguerre.jl

Lines changed: 1 addition & 48 deletions
Original file line numberDiff line numberDiff line change
@@ -302,51 +302,4 @@ function Multiplication(f::Fun{LaguerreWeight{H,T}},S::Laguerre) where {H<:Lague
302302
end
303303

304304

305-
## Integration
306-
function integrate(f::Fun{LW}) where LW <: LaguerreWeight{LL} where LL <: Laguerre
307-
α = space(f).α;
308-
n = length(f.coefficients);
309-
if space(f).space.α != α
310-
throw(ArgumentError("`integrate` is applicable if only LaguerreWeight parameter and Laguerre parameter are equal."))
311-
else
312-
if n == 0
313-
Fun(0)
314-
else
315-
if !isinteger(α) || α == 2
316-
if n == 1
317-
f̃ = Fun(f, JacobiWeight(α, 0, Chebyshev(Ray())));
318-
g = integrate(f̃);
319-
g = g - last(g)
320-
else
321-
if f.coefficients[1] == 0
322-
Fun(WeightedLaguerre(α + 1), f.coefficients[2:end] ./ (1:n-1))
323-
else
324-
f₀ = Fun(WeightedLaguerre(α), [f.coefficients[1]]);
325-
f₀ = Fun(f₀, JacobiWeight(α, 0, Chebyshev(Ray())));
326-
g₀ = integrate(f₀);
327-
g₀ = g₀ - last(g₀);
328-
g = Fun(WeightedLaguerre(α + 1), f.coefficients[2:end] ./ (1:n-1));
329-
g₀ + g
330-
end
331-
end
332-
else
333-
if n == 1
334-
f̃ = Fun(f, Chebyshev(Ray()));
335-
g = integrate(f̃);
336-
g = g - last(g)
337-
else
338-
if f.coefficients[1] == 0
339-
Fun(WeightedLaguerre(α + 1), f.coefficients[2:end] ./ (1:n-1))
340-
else
341-
f₀ = Fun(WeightedLaguerre(α), [f.coefficients[1]]);
342-
f₀ = Fun(f₀, Chebyshev(Ray()));
343-
g₀ = integrate(f₀);
344-
g₀ = g₀ - last(g₀);
345-
g = Fun(WeightedLaguerre(α + 1), f.coefficients[2:end] ./ (1:n-1));
346-
g₀ + g
347-
end
348-
end
349-
end
350-
end
351-
end
352-
end
305+

‎src/Spaces/Singularities/ExpWeight.jl

Lines changed: 0 additions & 53 deletions
This file was deleted.

0 commit comments

Comments
 (0)