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)))

    function transformback(tA)
      # reshape vector to matrix
      return reshape(tA, size(A0))

    obj(tA) = -I1(transformback(tA))
    result = optimize(obj, transform(A0), method = :nelder_mead)
    return transformback(result.minimum)


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

2 answers

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)...)


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



julia> enumerateArray(a)


This illustrates a couple of things.

  • Yes it really works
  • Matrices are stored in a basic format


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)

julia> size(a)

julia> b[4]=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)


The statistics recorded by @profile are printed using Profile.print ()



Also, most of the methods in Optim actually allow you to provide arrays, not just vectors. You can generalize the function nelder_mead

to do the same.



All Articles