Fortran-ähnliche Arrays wie FArray (Float64, -1: 1, -7: 7, -128: 512) in Julia
Im Allgemeinen ist es eine gute Entscheidung, ein 1-basiertes Array für Julia zu haben, aber manchmal ist es wünschenswert, ein Fortran-ähnliches Array mit Indizes zu haben, die einige Unterbereiche von ℤ umfassen:
julia> x = FArray(Float64, -1:1,-7:7,-128:512)
wo es sinnvoll wäre:
in den Codes zum Buch Numerische Lösung hyperbolischer partieller Differentialgleichungen von prof. John A. Trangenstein Diese negativen Indizes werden intensiv für Geisterzellen für Randbedingungen verwendet. Gleiches gilt für Clawpack (steht für „Conservation Laws Package“) von prof. Randall J. LeVequehttp://depts.washington.edu/clawpack/und es gibt viele andere Codes, in denen solche Indizes natürlich wären. Eine solche Hilfsklasse wäre also nützlich, um solche Codes schnell zu übersetzen.
Ich habe gerade damit begonnen, einen solchen Hilfstyp zu implementieren, aber da ich für Julia noch ziemlich neu bin, wäre Ihre Hilfe sehr dankbar.
Ich habe angefangen mit:
type FArray
ranges
array::Array
function FArray(T, r::Range1{Int}...)
dims = map((x) -> length(x), r)
array = Array(T, dims)
new(r, array)
end
end
Ausgabe:
julia> include ("FortranArray.jl")
julia> x = FArray(Float64, -1:1,-7:7,-128:512)
FArray((-1:1,-7:7,-128:512),3x15x641 Array{Float64,3}:
[:, :, 1] =
6.90321e-310 2.6821e-316 1.96042e-316 0.0 0.0 0.0 9.84474e-317 … 1.83233e-316 2.63285e-316 0.0 9.61618e-317 0.0
6.90321e-310 6.32404e-322 2.63285e-316 0.0 0.0 0.0 2.63292e-316 2.67975e-316
...
[:, :, 2] =
...
Da ich für Julia völlig neu bin, wären alle Empfehlungen, die zu mehr Effizienz führen, sehr dankbar.
[Bearbeiten]
Das Thema wurde hier diskutiert:
https://groups.google.com/forum/#!topic/julia-dev/NOF6MA6tb9Y
Während der Diskussion wurden zwei Möglichkeiten ausgearbeitet, Julia-Arrays mit beliebiger Basis zu haben: SubArray-basierte Beispielnutzung mit einer Hilfsfunktion:
function farray(T, r::Range1{Int64}...)
dims = map((x) -> length(x), r)
array = Array(T, dims)
sub_indices = map((x) -> -minimum(x) + 2 : maximum(x), r)
sub(array, sub_indices)
end
julia> y[-1,-7,-128] = 777
777
julia> y[-1,-7,-128] + 33
810.0
julia> y[-2,-7,-128]
ERROR: BoundsError()
in getindex at subarray.jl:200
julia> y[2,-7,-128]
2.3977385e-316
Bitte beachten Sie, dass die Grenzen nicht vollständig überprüft werden. Weitere Details finden Sie hier:https://github.com/JuliaLang/julia/issues/4044
Im Moment hat SubArray Leistungsprobleme, aber möglicherweise wird seine Leistung verbessert, siehe auch:
https://github.com/JuliaLang/julia/issues/5117
https://github.com/JuliaLang/julia/issues/3496
Ein weiterer Ansatz, der im Moment eine bessere Leistung aufweist, überprüft außerdem beide Grenzen:
type FArray{T<:Number, N, A<:AbstractArray} <: AbstractArray
ranges
offsets::NTuple{N,Int}
array::A
function FArray(r::Range1{Int}...)
dims = map((x) -> length(x), r)
array = Array(T, dims)
offs = map((x) -> 1 - minimum(x), r)
new(r, offs, array)
end
end
FArray(T, r::Range1{Int}...) = FArray{T, length(r,), Array{T, length(r,)}}(r...)
getindex{T<:Number}(FA::FArray{T}, i1::Int) = FA.array[i1+FA.offsets[1]]
getindex{T<:Number}(FA::FArray{T}, i1::Int, i2::Int) = FA.array[i1+FA.offsets[1], i2+FA.offsets[2]]
getindex{T<:Number}(FA::FArray{T}, i1::Int, i2::Int, i3::Int) = FA.array[i1+FA.offsets[1], i2+FA.offsets[2], i3+FA.offsets[3]]
setindex!{T}(FA::FArray{T}, x, i1::Int) = arrayset(FA.array, convert(T,x), i1+FA.offsets[1])
setindex!{T}(FA::FArray{T}, x, i1::Int, i2::Int) = arrayset(FA.array, convert(T,x), i1+FA.offsets[1], i2+FA.offsets[2])
setindex!{T}(FA::FArray{T}, x, i1::Int, i2::Int, i3::Int) = arrayset(FA.array, convert(T,x), i1+FA.offsets[1], i2+FA.offsets[2], i3+FA.offsets[3])
getindex und setindex! Methoden für FArray wurden vom Code base / array.jl inspiriert.
Anwendungsfälle:
julia> y = FArray(Float64, -1:1,-7:7,-128:512);
julia> y[-1,-7,-128] = 777
777
julia> y[-1,-7,-128] + 33
810.0
julia> y[-1,2,3]
0.0
julia> y[-2,-7,-128]
ERROR: BoundsError()
in getindex at FortranArray.jl:27
julia> y[2,-7,-128]
ERROR: BoundsError()
in getindex at FortranArray.jl:27