How to quickly rescale arrays
In the following code, I use the Julia Optim package to find the optimal matrix relative to an objective function. Unfortunately, the provided optimization function only supports vectors, so I have to convert the matrix to a vector before passing it to the optimization function, and also convert it when used in the objective function.
function opt(A0,X)
I1(A) = sum(maximum(X*A,1))
function transform(A)
# reshape matrix to vector
return reshape(A,prod(size(A)))
end
function transformback(tA)
# reshape vector to matrix
return reshape(tA, size(A0))
end
obj(tA) = -I1(transformback(tA))
result = optimize(obj, transform(A0), method = :nelder_mead)
return transformback(result.minimum)
end
I think Julia allocates a new space for this every time, and it feels slow, so what would be a more efficient way to deal with this problem?
source to share
As long as arrays contain elements that are considered immutable, including all primitives, array elements are contained in 1 large contiguous block of memory. This way, you can break the measurement rules and just treat the two-dimensional array as the one-dimensional array you want to do. This way you don't need to change the shape, but I don't think your problem will change.
Arrays are primary columns and contiguous
Consider the following function
function enumerateArray(a)
for i = 1:*(size(a)...)
print(a[i])
end
end
This function multiplies all sizes of a together and then loops from 1 to that number if a is one-dimensional.
When you define a as next
julia> a = [ 1 2; 3 4; 5 6]
3x2 Array{Int64,2}:
1 2
3 4
5 6
Result
julia> enumerateArray(a)
135246
This illustrates a couple of things.
- Yes it really works
- Matrices are stored in a basic format
Reshape
So the question is, why doesn't this change the use of this fact? Good. Here's the julia source for the change in array.c
a = (jl_array_t*)allocobj((sizeof(jl_array_t) + sizeof(void*) + ndimwords*sizeof(size_t) + 15)&-16);
So yes, a new array is created, but only new size information is created, it points to the original data which is not copied. You can check it simply like this:
b = reshape(a,6);
julia> size(b)
(6,)
julia> size(a)
(3,2)
julia> b[4]=100
100
julia> a
3x2 Array{Int64,2}:
1 100
3 4
5 6
So setting the 4th element to b sets (1,2) to element a .
As for the overall slowness
I1(A) = sum(maximum(X*A,1))
will create a new array.
You can use multiple macros to track this parameter @profile and @time . Time additionally records the amount of memory allocated and can be placed before any expression.
for example
julia> A = rand(1000,1000);
julia> X = rand(1000,1000);
julia> @time sum(maximum(X*A,1))
elapsed time: 0.484229671 seconds (8008640 bytes allocated)
266274.8435928134
The statistics recorded by @profile are printed using Profile.print ()
source to share