Which notation is best?

[MikaelSlevinsky]

There are multiple notations used to create multivariate operators:

MultivariateOrthogonalPolynomials.jl/src/triangle.jl

Lines 89 to 128 in 04704a8

	function Wx(a,b,c)
	n = mortar(Fill.(oneto(∞),oneto(∞)))
	k = mortar(Base.OneTo.(oneto(∞)))
	dat = BlockVcat(
	((k .+ (c-1)) .* ( k .- n .- 1 ) ./ (2k .+ (b+c-1)))',
	(k .* (k .- n .- a) ./ (2k .+ (b+c-1)))'
	)
	_BandedBlockBandedMatrix(dat, axes(k,1), (1,-1), (1,0))
	end
	function Rx(a,b,c)
	n = mortar(Fill.(oneto(∞),oneto(∞)))
	k = mortar(Base.OneTo.(oneto(∞)))
	dat = BlockHcat(
	BroadcastVector((n,k,bc1,abc) -> (n + k + bc1) / (2n + abc), n, k, b+c-1, a+b+c),
	BroadcastVector((n,k,abc) -> (n + k + abc) / (2n + abc), n, k, a+b+c)
	)
	_BandedBlockBandedMatrix(dat', axes(k,1), (0,1), (0,0))
	end
	function Ry(a,b,c)
	n = mortar(Fill.(oneto(∞),oneto(∞)))
	k = mortar(Base.OneTo.(oneto(∞)))
	dat = PseudoBlockArray(Vcat(
	((k .+ (c-1) ) .* (n .+ k .+ (b+c-1)) ./ ((2n .+ (a+b+c)) .* (2k .+ (b+c-1) )))',
	((k .- n .- a ) .* (k .+ (b+c)) ./ ((2n .+ (a+b+c)) .* (2k .+ (b+c-1) )))',
	((k .+ (c-1) ) .* (k .- n .- 1) ./ ((2n .+ (a+b+c)) .* (2k .+ (b+c-1) )))',
	((n .+ k .+ (a+b+c) ) .* (k .+ (b+c)) ./ ((2n .+ (a+b+c)) .* (2k .+ (b+c-1) )))'
	), (blockedrange(Fill(2,2)), axes(n,1)))
	_BandedBlockBandedMatrix(dat, axes(k,1), (0,1), (0,1))
	end
	function Lx(a,b,c)
	n = mortar(Fill.(oneto(∞),oneto(∞)))
	k = mortar(Base.OneTo.(oneto(∞)))
	dat = PseudoBlockArray(Vcat(
	((n .- k .+ a) ./ (2n .+ (a+b+c)))',
	((n .- k .+ 1) ./ (2n .+ (a+b+c)))'
	), (blockedrange(Fill(1,2)), axes(n,1)))
	_BandedBlockBandedMatrix(dat, axes(k,1), (1,0), (0,0))
	end

When to use the dot notation and when to use BroadcastVector? Is it better to use @. or is that wasteful because it creates some unnecessary broadcasting?

(My impression is that these are due to successive invasions of new notation and not all operators are fully updated, but I could be way off.)

[dlfivefifty]

That's mostly just me trying different things to get building operators to be faster. The big issue with @. is it makes nested broadcasting, which at some point breaks type inference. The manual . versions suffer the same problem, but not as bad.

It's likely that the BroadcastVector(...) version was a specific case where the nested broadcasting was too hard for type inference so by doing this we avoid nested broadcasting completely.

But this would have been a one-off experiment so has not been thoroughly investigated.

[MikaelSlevinsky]

Does type assertion solve issues related to inference?

[dlfivefifty]

No. Here's a rough idea of the problem. Let's say we have a complicated nested BroadcastVector like

v = ((k .+ (c-1)) .* ( k .- n .- 1 ) ./ (2k .+ (b+c-1)))

When building operators the code basically lowers to

for K in Block.(1:N)
   copyto!(view(dest,K), view(v,K))
end

It's smart enough to recognise view(v,K) is itself equivalent to a BroadcastVector thus this becomes equivalent to:

k_K = view(k,K) # returns Base.OneTo(Int(K))
n_K = view(n,K) # returns Fill(Int(K),Int(K))
copyto!(view(dest,K), ((k_K .+ (c-1)) .* ( k_K .- n_K .- 1 ) ./ (2k_K .+ (b+c-1)))) # actually a nested Broadcasted 

When correctly type-inferred the last step then gets reduced to a tight loop, essentially as fast as if done manually. When type-inferrence fails, it doesn't know how to reduce the nested Broadcasted to a simple loop, and hence has to do dynamic dispatch looping over every entry of view(dest,K), so extremely slow.