I am using Julia 1.9.0 and I’m trying to fill an array using views and dot notation where each entry depends on the preceding one(s). Say, I want the first 10 entries of the Fibonacci series in my array, so I initialize the first two entries and want the rest computed without a for-loop.
My approach was:
v = zeros(Int64, 10)
v[1:2] .= 1
@views v[3:10] .= v[1:8] .+ v[2:9]
v
Output:
10-element Vector{Int64}:
1
1
2
1
0
0
0
0
0
0
This doesn’t work, though, because apparently v
doesn’t get updated before the next value is calculated. From reading the docs, I couldn’t find out why, because as far as I understood no arrays are copied and no temporary array is created. This should be equivalent to:
v = zeros(Int64, 10)
v[1:2] .= 1
@views a, b, c = v[1:8], v[2:9], v[3:10]
for i in 1:8
c[i] = a[i] + b[i]
end
v
Output:
10-element Vector{Int64}:
1
1
2
3
5
8
13
21
34
55
Which is exactly the result I expected to get. Why doesn’t the first version work and is it possible to solve this without the for-loop?
Broadcasting does alias detection. Therefore actually your operation copies data, because you have aliases.
To convince yourself that this is what happens check this code:
julia> v = zeros(Int64, 10^7);
julia> v[1:2] .= 1;
julia> @allocated v[3:end] .= @views v[1:end-2] .+ v[2:end-1]
This check is done using the Base.mightalias
function.
Dealiasing is done to protect RHS to be affected by mutation of the LHS. So, on a basic level, the answer is that with broadcasting you cannot implement recursion.
However, you could create a custom type that would always return false
in Base.mightalias
aliasing check in which case what you ask for would work.
Here is an example:
julia> using SentinelArrays
julia> v = zeros(Int64, 10^7);
julia> v[1:2] .= 1;
julia> @allocated v[3:end] .= ChainedVector([@view v[1:end-2]]) .+ ChainedVector([@view v[2:end-1]])
800
julia> v
10000000-element Vector{Int64}:
1
1
2
3
5
8
13
21
34
55
89
144
233
377
610
987
⋮