[julia-users] Re: unsafe_wrap from a pointer from a Julia array

[Gunnar Farnebäck]

Is this what you're looking for?

julia> using SIMD

julia> v = rand(2, 5)
2×5 Array{Float64,2}:
 0.391832  0.785861  0.352291  0.874891  0.874593
 0.697768  0.51637   0.711525  0.400077  0.266559

julia> w = reinterpret(Vec{2, Float64}, vec(v))
5-element Array{SIMD.Vec{2,Float64},1}:
 Float64⟨0.391832,0.697768⟩
 Float64⟨0.785861,0.51637⟩ 
 Float64⟨0.352291,0.711525⟩
 Float64⟨0.874891,0.400077⟩
 Float64⟨0.874593,0.266559⟩


Den onsdag 16 november 2016 kl. 05:52:14 UTC+1 skrev Sheehan Olver:
>
>
> In some high performance code I want to convert an 2 x n `Matrix{Float64}` 
> to a n long `Matrix{Vec{2,Float64}}`with no memory allocation (both types 
> will store memory in the same order).  The following works:
>
> *v=rand(2,100)*
> *w=unsafe_wrap(Vector{Vec{2,Float64}},*
> *reinterpret(Ptr{Vec{2,Float64}},pointer(v)),*
> *size(pts,1),true)*
>
> The variable *v* is a throw away variable, however, and as I understand 
> it, the way its written the memory for *w* will be freed as soon as *v* 
> is no longer referenced.  Is there any way to tell the garbage collector to 
> not collect *v*?
>

Re: [julia-users] GC rooting for embedding: what is safe and unsafe?

[Gunnar Farnebäck]

Thanks. That makes things clearer.

Den fredag 14 oktober 2016 kl. 14:16:54 UTC+2 skrev Yichao Yu:
>
>
>
> On Fri, Oct 14, 2016 at 7:03 AM, Bart Janssens <ba...@bartjanssens.org 
> > wrote:
>
>> Hi,
>>
>> Replies below, to the best of my understanding of the Julia C interface:
>>
>> On Fri, Oct 14, 2016 at 11:47 AM Gunnar Farnebäck <gun...@lysator.liu.se 
>> > wrote:
>>
>>> Reading through the threads and issues on gc rooting for embedded code, 
>>> as well as the code comments above the JL_GC_PUSH* macros in julia.h, I'm 
>>> still uncertain about how defensive it's necessary to be and best 
>>> practices. I'll structure this into a couple of cases with questions.
>>>
>>> 1. One of the use cases in examples/embedding.c is written:
>>>
>>> jl_function_t *func = jl_get_function(jl_base_module, "sqrt");
>>> jl_value_t* argument = jl_box_float64(2.0);
>>> jl_value_t* ret = jl_call1(func, argument);
>>>
>>> if (jl_is_float64(ret)) {
>>> double retDouble = jl_unbox_float64(ret);
>>> printf("sqrt(2.0) in C: %e\n", retDouble);
>>> }
>>>
>>>  
>>>
>> Is this missing gc rooting for argument during the call to jl_call1 or is 
>>> it safe without it?
>>> Would ret need gc rooting to be safe during the calls to jl_is_float64 
>>> and/or jl_unbox_float64?
>>>
>>
>> The jl_call argument must be rooted since func may allocate. I don't 
>> think the operations on ret allocate, but if you're unsure it's better to 
>> root it. Also, as your code evolves you may decide to perform extra 
>> operations on ret and then it's easy to forget the GC rooting at that 
>> point, so I'd root ret here.
>>
>
> jl_call1 (and other jl_call* functions) are special in the sense that it 
> roots it's argument so you don't have to root `argument`
> The ret doesn't have to be rooted if these are all what you are doing with 
> it.
> The only one that should in principle be rooted is actually `func`. 
> However, since it is known to be a global constant you won't get into any 
> trouble without. (If it's a non-const global then you have to)
>  
>
>>  
>>
>>>
>>> 2. 
>>> jl_value_t *a = jl_box_float64(1.0);
>>> jl_value_t *b = jl_box_float64(2.0);
>>> JL_GC_PUSH2(, );
>>>
>>> Is this unsafe since a is not yet rooted during the second call to 
>>> jl_box_float64 and must instead be written like below?
>>>
>>> jl_value_t *a = 0;
>>> jl_value_t *b = 0;
>>> JL_GC_PUSH2(, );
>>> a = jl_box_float64(1.0);
>>> b = jl_box_float64(2.0);
>>>
>>> For a single variable it's just fine to do like this though?
>>> jl_value_t *a = jl_box_float64(1.0);
>>> JL_GC_PUSH1();
>>>
>>>
>> Yes, since jl_box will allocate.
>>
>
> This is correct.
>  
>
>>  
>>
>>> 3. Are
>>> jl_function_t *func = jl_get_function(jl_base_module, "println");
>>> jl_value_t *a = 0;
>>> jl_value_t *b = 0;
>>> JL_GC_PUSH2(, );
>>> a = jl_box_float64(1.0);
>>> b = jl_box_float64(2.0);
>>> jl_call2(func, a, b);
>>> JL_GC_POP();
>>>
>>> and
>>>
>>> jl_function_t *func = jl_get_function(jl_base_module, "println");
>>> jl_value_t **args;
>>> JL_GC_PUSHARGS(args, 2);
>>> args[0] = jl_box_float64(1.0);
>>> args[1] = jl_box_float64(2.0);
>>> jl_call(func, args, 2);
>>> JL_GC_POP();
>>>
>>> equivalent and both safe? Are there any reasons to choose one over the 
>>> other, apart from personal preferences?
>>>
>>
>> They are equivalent, looking at the code for the macro it seems that the 
>> JL_GC_PUSHARGS variant heap-allocates the array of pointers to root, so 
>> that might be slightly slower. I'd only use the JL_GC_PUSHARGS version if 
>> the number of arguments comes from a variable or a parameter or similar.
>>
>
> No JL_GC_PUSHARGS does **NOT** heap allocate the array.
> The difference between the two is pretty small and the performance should 
> be almost not noticeable unless you are doing a lot of things with the 
> boxed variables.
>  
>
>>  
>>
>>>
>>> 4. Can any kind of exception checking be done safely without rooting the 
>>> return value?
>>> jl_value_t *ret = jl_call1(...);
>>> if (jl_exception_occured())
>>> printf("%s \

[julia-users] GC rooting for embedding: what is safe and unsafe?

[Gunnar Farnebäck]

Reading through the threads and issues on gc rooting for embedded code, as 
well as the code comments above the JL_GC_PUSH* macros in julia.h, I'm 
still uncertain about how defensive it's necessary to be and best 
practices. I'll structure this into a couple of cases with questions.

1. One of the use cases in examples/embedding.c is written:

jl_function_t *func = jl_get_function(jl_base_module, "sqrt");
jl_value_t* argument = jl_box_float64(2.0);
jl_value_t* ret = jl_call1(func, argument);

if (jl_is_float64(ret)) {
double retDouble = jl_unbox_float64(ret);
printf("sqrt(2.0) in C: %e\n", retDouble);
}

Is this missing gc rooting for argument during the call to jl_call1 or is 
it safe without it?
Would ret need gc rooting to be safe during the calls to jl_is_float64 
and/or jl_unbox_float64?

2. 
jl_value_t *a = jl_box_float64(1.0);
jl_value_t *b = jl_box_float64(2.0);
JL_GC_PUSH2(, );

Is this unsafe since a is not yet rooted during the second call to 
jl_box_float64 and must instead be written like below?

jl_value_t *a = 0;
jl_value_t *b = 0;
JL_GC_PUSH2(, );
a = jl_box_float64(1.0);
b = jl_box_float64(2.0);

For a single variable it's just fine to do like this though?
jl_value_t *a = jl_box_float64(1.0);
JL_GC_PUSH1();

3. Are
jl_function_t *func = jl_get_function(jl_base_module, "println");
jl_value_t *a = 0;
jl_value_t *b = 0;
JL_GC_PUSH2(, );
a = jl_box_float64(1.0);
b = jl_box_float64(2.0);
jl_call2(func, a, b);
JL_GC_POP();

and

jl_function_t *func = jl_get_function(jl_base_module, "println");
jl_value_t **args;
JL_GC_PUSHARGS(args, 2);
args[0] = jl_box_float64(1.0);
args[1] = jl_box_float64(2.0);
jl_call(func, args, 2);
JL_GC_POP();

equivalent and both safe? Are there any reasons to choose one over the 
other, apart from personal preferences?

4. Can any kind of exception checking be done safely without rooting the 
return value?
jl_value_t *ret = jl_call1(...);
if (jl_exception_occured())
printf("%s \n", jl_typeof_str(jl_exception_occured()));
else
d = jl_unbox_float64(ret);

5. What kind of costs, other than increased code size, can be expected from 
overzealous gc rooting?

6. Is it always safe not to gc root the function pointer returned from 
jl_get_function?

[julia-users] Re: help understanding different ways of wrapping functions

[Gunnar Farnebäck]

It's normal that manually inlined code of this kind is faster than wrapped 
code unless the compiler manages to see the full inlining potential. In 
this case the huge memory allocations for the wrapped solutions indicates 
that it's nowhere near doing that at all. I doubt it will take you all the 
way but start with modifying your inner M_CPS function to only take 
positional arguments or declaring the type of the keyword argument as 
suggested in the performance tips section of the manual.

Den onsdag 28 september 2016 kl. 06:29:37 UTC+2 skrev K leo:
>
> I tested a few different ways of wrapping functions.  It looks different 
> ways of wrapping has slightly different costs.  But the most confusing to 
> me is that putting everything inline looks much faster than wrapping things 
> up.  I would understand this in other languages, but I thought Julia 
> advocates simple wrapping.  Can anyone help explain what is happening 
> below, and how I can do most efficient wrapping in the demo code?
>
> Demo code is included below.
>
> julia> versioninfo()
> Julia Version 0.5.0
> Commit 3c9d753 (2016-09-19 18:14 UTC)
> Platform Info:
>   System: Linux (x86_64-pc-linux-gnu)
>   CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
>   WORD_SIZE: 64
>   BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY Haswell)
>   LAPACK: libopenblas64_
>   LIBM: libopenlibm
>   LLVM: libLLVM-3.7.1 (ORCJIT, broadwell)
>
> julia> testFunc()
> calling LoopCP (everything inline)
>   0.097556 seconds (2.10 k allocations: 290.625 KB)
> elapsed time (ns): 97555896
> bytes allocated:   297600
> pool allocs:   2100
> [0.0,4200.0,0.0,0.0,4200.0,4200.0,4200.0,4200.0,0.0,4200.0,4200.0]
>
> calling LoopCP0 (slightly wrapped)
>   4.173830 seconds (49.78 M allocations: 2.232 GB, 5.83% gc time)
> elapsed time (ns): 4173830495
> gc time (ns):  243516584
> bytes allocated:   2396838538
> pool allocs:   49783357
> GC pauses: 104
> full collections:  1
> [4200.0,0.0,4200.0,4200.0,0.0,0.0,0.0,0.0,4200.0,0.0,0.0]
>
> calling LoopCP1 (wrapped one way)
>   5.274723 seconds (59.59 M allocations: 2.378 GB, 3.62% gc time)
> elapsed time (ns): 5274722983
> gc time (ns):  191036337
> bytes allocated:   2553752638
> pool allocs:   59585834
> GC pauses: 112
> [8400.0,0.0,8400.0,8400.0,0.0,0.0,0.0,0.0,8400.0,0.0,0.0]
>
> calling LoopCP2 (wrapped another way)
>   5.212895 seconds (59.58 M allocations: 2.378 GB, 3.60% gc time)
> elapsed time (ns): 5212894550
> gc time (ns):  187696529
> bytes allocated:   2553577600
> pool allocs:   59582100
> GC pauses: 111
> [0.0,8400.0,0.0,0.0,8400.0,8400.0,8400.0,8400.0,0.0,8400.0,8400.0]
>
> const dim=1000
>>
>>
>>> type Tech
>>
>> a::Array{Float64,1}
>>
>> c::Array{Int,1}
>>
>>
>>> function Tech()
>>
>> this = new()
>>
>> this.a = zeros(Float64, dim)
>>
>> this.c = rand([0,1;], dim)
>>
>> this
>>
>> end
>>
>> end
>>
>>
>>> function LoopCP(csign::Int, tech::Tech)
>>
>> for j=1:10
>>
>> for xRat in [1.:20.;]
>>
>> @inbounds for i = 1:dim
>>
>> if csign == tech.c[i]
>>
>> tech.a[i] += 2.*xRat
>>
>> else
>>
>> tech.a[i] = 0.
>>
>> end
>>
>> end
>>
>> end #
>>
>> end
>>
>> nothing
>>
>> end
>>
>>
>>> function M_CPS(i::Int, csign::Int, tech::Tech; xRat=0.)
>>
>> if csign == tech.c[i]
>>
>> tech.a[i] += 2.*xRat
>>
>> else
>>
>> tech.a[i] = 0.
>>
>> end
>>
>> nothing
>>
>> end
>>
>>
>>> function LoopCP0(csign::Int, tech::Tech)
>>
>> for j=1:10
>>
>> for xRat in [1.:20.;]
>>
>> @inbounds for i = 1:dim
>>
>> M_CPS(i, csign, tech, xRat=xRat)
>>
>> end
>>
>> end #
>>
>> end
>>
>> nothing
>>
>> end
>>
>>
>>> function MoleculeWrapS(csign::Int, tech::Tech, molecule::Function, 
>>> xRat=0.)
>>
>> @inbounds for i = 1:dim
>>
>> molecule(i, csign, tech; xRat=xRat)
>>
>> end
>>
>> nothing
>>
>> end
>>
>>
>>> function LoopRunnerM1(csign::Int, tech::Tech, molecule::Function)
>>
>> for j=1:10
>>
>> for xRat in [1.:20.;]
>>
>> MoleculeWrapS(csign, tech, molecule, xRat)
>>
>> end #
>>
>> end
>>
>> nothing
>>
>> end
>>
>>
>>> LoopCP1(csign::Int, tech::Tech) = LoopRunnerM1(csign, tech, M_CPS)
>>
>>
>>> WrapCPS(csign::Int, tech::Tech, xRat=0.) = MoleculeWrapS(csign, tech, 
>>> M_CPS, xRat)
>>
>>
>>> function LoopRunnerM2(csign::Int, tech::Tech, loop::Function)
>>
>> for j=1:10
>>
>> for xRat in [1.:20.;]
>>
>> loop(csign, tech, xRat)
>>
>> end #
>>
>> end
>>
>> nothing
>>
>> end
>>
>>
>>> LoopCP2(csign::Int, tech::Tech) = LoopRunnerM2(csign, tech, WrapCPS)
>>
>>
>>> function testFunc()
>>
>> tech = Tech()
>>
>> nloops = 100
>>
>>
>>> println("calling LoopCP (everything inline)")
>>
>> tech.a = zeros(tech.a)

[julia-users] Re: map(mean, zip(v...)) doesn't scale

[Gunnar Farnebäck]

I realize the question is about the performance of the zip approach but may 
I suggest just running mean(v)?

On Wednesday, August 10, 2016 at 11:34:08 AM UTC+2, Tamas Papp wrote:
>
> I was benchmarking code for calculating means of vectors, eg 
>
> v = [rand(100) for i in 1:100] 
>
> m1(v) = map(mean, zip(v...)) 
> m2(v) = mapslices(mean, hcat(v...), 2) 
>
> @elapsed m1(v) 
> @elapsed m2(v) 
>
> m2 is faster, but for 1000 vectors, 
>
> v = [rand(100) for i in 1:1000] 
>
> m1 with zip takes "forever" (CPU keeps running at 100%, I interrupted 
> after a 10 minutes). Is this supposed to happen, or should I report it 
> as a bug? 
>
> Tamas 
>

[julia-users] Re: Strange performance issue in filling in a matrix column

[Gunnar Farnebäck]

fill_W1! allocates memory because it makes copies when constructing the 
right hand sides. fill_W2 allocates memory in order to construct the 
comprehensions (that you then discard). In both cases memory allocation 
could plausibly be avoided by a sufficiently smart compiler, but until 
Julia becomes that smart, have a look at the sub function to provide views 
instead of copies for the right hand sides of fill_W1!.

On Thursday, July 21, 2016 at 5:07:34 PM UTC+2, Michael Prange wrote:
>
> I'm a new user, so have mercy in your responses. 
>
> I've written a method that takes a matrix and vector as input and then 
> fills in column icol of that matrix with the vector of given values that 
> have been shifted upward by ishift indices with periodic boundary 
> conditions. To make this clear, given the matrix
>
> W = [1  2
> 3  4
> 5  6]
>
> the vector w = [7  8  9], icol = 2 and ishift = 1, the new value of W is 
> given by
>
> W = [1  8
> 3  9
> 5  7]
>
> I need a fast way of doing this for large matrices. I wrote three methods 
> that should (In my naive mind) give the same performance results, but @time 
> reports otherwise.  The method definitions and the performance results are 
> given below. Can someone teach me why the results are so different? The 
> method fill_W! is too wordy for my tastes, but the more compact notation in 
> fill_W1! and fill_W2! achieve poorer results. Any why do these latter two 
> methods allocate so much memory when the whole point of these methods is to 
> use already-allocated memory.
>
> Michael
>
> ### Definitions
>
>
> function fill_W1!{TF}(W::Matrix{TF}, icol::Int, w::Vector{TF}, 
> ishift::Int)
> @assert(size(W,1) == length(w), "Dimension mismatch between W and w")
> W[1:(end-ishift),icol] = w[(ishift+1):end]
> W[(end-(ishift-1)):end,icol] = w[1:ishift]
> return
> end
>
>
> function fill_W2!{TF}(W::Matrix{TF}, icol::Int, w::Vector{TF}, 
> ishift::Int)
> @assert(size(W,1) == length(w), "Dimension mismatch between W and w")
> [W[i,icol] = w[i+ishift] for i in 1:(length(w)-ishift)]
> [W[end-ishift+i,icol] = w[i] for i in 1:ishift]
> return
> end
>
>
> function fill_W!{TF}(W::Matrix{TF}, icol::Int, w::Vector{TF}, 
> ishift::Int)
> @assert(size(W,1) == length(w), "Dimension mismatch between W and w")
> n = length(w)
> for j in 1:(n-ishift)
> W[j,icol] = w[j+ishift]
> end
> for j in (n-(ishift-1)):n
> W[j,icol] = w[j-(n-ishift)]
> end
> end
>
>
> # Performance Results
> julia>
> W = rand(100,2)
> w = rand(100)
> println("fill_W!:")
> println(@time fill_W!(W, 2, w, 2))
> println("fill_W1!:")
> println(@time fill_W1!(W, 2, w, 2))
> println("fill_W2!:")
> println(@time fill_W2!(W, 2, w, 2))
>
>
> Out>
> fill_W!:
>  0.002801 seconds (4 allocations: 160 bytes)
> nothing
> fill_W1!:
>  0.007427 seconds (9 allocations: 7.630 MB)
> [0.152463397611579,0.6314166578356002]
> fill_W2!:
>  0.005587 seconds (7 allocations: 7.630 MB)
> [0.152463397611579,0.6314166578356002]
>
>
>

Re: [julia-users] why does Julia have both map and broadcast?

[Gunnar Farnebäck]

In case the accidental matrix happens to be huge you can also get into all 
sorts of memory shortage related inconveniences.

Den tisdag 24 maj 2016 kl. 07:27:24 UTC+2 skrev Tony Kelman:
>
> I can see some situations where getting a matrix instead of a vector 
> because one of the inputs was unexpectedly transposed, then sending the 
> output to backslash or similar wouldn't trigger an error, just give a 
> meaningless and unexpected result that would be really difficult to debug.
>
>
> On Sunday, May 22, 2016 at 3:30:11 PM UTC-7, Stefan Karpinski wrote:
>>
>> The only reason I can think of for map to not just do what broadcast does 
>> is map's more strict requirements could potentially catch errors. However, 
>> I suspect that as soon as you used the result, you would rapidly get an 
>> error anyway, so it's possible that we should merge these. There may be 
>> some other reason. Probably a good question for julia-dev or a GitHub issue.
>>
>> On Saturday, May 21, 2016, Tony Kelman  wrote:
>>
>>> I don't think the behavior of repeating scalars, or expanding outer 
>>> product dimensions, is really map mathematically, it's something a little 
>>> different. Dimension mismatches are likely a bug when you really want map.
>>>
>>>
>>> On Saturday, May 21, 2016 at 10:23:36 AM UTC-7, Steven G. Johnson wrote:



 On Saturday, May 21, 2016 at 12:36:31 PM UTC-4, Tony Kelman wrote:
>
> I can't think of any aside from the behavior when dimensions don't 
> match, which falls under your exceptions... exception.
>
> Maybe map could have a flag to set whether to be picky about matching 
> dimensions, or behave like broadcast if not?
>

 Why not just have map behave all the time like broadcast? 

>>>

[julia-users] Re: Declaration of types in function construction...

[Gunnar Farnebäck]

Your problem is probably that you try to call your function with four 
positional arguments, although the last two should have been named. Compare

julia> f(x, y; z=0) = x + y + z
f (generic function with 1 method)

julia> f(2, 3)
5

julia> f(2, 3, 4)
ERROR: MethodError: `f` has no method matching f(::Int64, ::Int64, ::Int64)
Closest candidates are:
  f(::Any, ::Any)

julia> f(2, 3, z = 4)
9

Den torsdag 12 maj 2016 kl. 12:51:13 UTC+2 skrev Charles Ll:
>
> Dear all,
>
> There is something I did not understood well in Julia regarding the type 
> fo the variables and this creates difficulties with the functions I am 
> writing.
>
> For one of the functions of Spectra.jl package, I wrote for instance: 
>
> function gaussianarea(Amplitude::Array{Float64},HWHM::Array{Float64}; 
> eseAmplitude::Array{Float64} = 0, eseHWHM::Array{Float64} = 0)
>
> I was thinking that using the ::Array{Float64} will allow users to enter 
> either vectors or arrays. However, when I try to enter a vectors for 
> Amplitude or HWHM for instance, I get the following error:
>
> LoadError: MethodError: `gaussianarea` has no method matching 
> gaussianarea(::Array{Float64,1}, ::Array{Float64,1}, ::Array{Float64,1}, 
> ::Array{Float64,1})
> Closest candidates are:
>   gaussianarea(::Array{Float64,N}, ::Array{Float64,N})
> while loading In[46], in expression starting on line 9
>
>
> I'm quite disappointed by that and don't know how to avoid it... Because 
> for me a vector is an Array {Float64,1}... But the above message says no?
>
> What did I miss here? 
>

[julia-users] Re: compact syntax for mapping vector to tuple

[Gunnar Farnebäck]

This is at least shorter.

julia> (v.^2...)
(1,4,81)

Den torsdag 31 mars 2016 kl. 10:28:09 UTC+2 skrev Tamas Papp:
>
> Hi, 
>
> Is there a preferred syntax for mapping a vector to a tuple? Eg an 
> alternative for (contrived example) 
>
> v = [1,2,9] 
> ntuple(i -> v[i]^2, length(v)) 
> tuple([e^2 for e in v]...) 
>
> Best, 
>
> Tamas 
>

[julia-users] Re: Concurrently install two versions of Julia

[Gunnar Farnebäck]

The environment variable JULIA_PKGDIR may be of use.

Den måndag 7 mars 2016 kl. 05:44:29 UTC+1 skrev Vishnu Raj:
>
> In my system I have Juno with 0.4.2 and julia with 0.4.3. Every time I 
> start one, all my packages are getting recompiled. They are all under same 
> v0.4 folder and I think that's why this is happening.
> Is there a way to avoid this?
>
> On Sunday, March 6, 2016 at 3:28:40 PM UTC+5:30, Andreas Lobinger wrote:
>>
>> Depends on your definition of installed and what system you use. I use 
>> since 0.2 (and the 0.3dev) a local build -on a linux system- and this quite 
>> nicely encapsulated in a single directory inside my home so they live in 
>> parallel. The package directory which is inside .julia is has versioning 
>> (v0.3/v0.4/v0.5), too.
>>
>> On Saturday, March 5, 2016 at 11:48:37 PM UTC+1, Pulkit Agarwal wrote:
>>>
>>> Hi,
>>>
>>> Is there a way to have the stable version of Julia as the global Julia 
>>> (i.e. something which can be accessed by calling `julia` from the command 
>>> line) and the development version of Julia (which will be stored in some 
>>> other folder).
>>>
>>> Thanks,
>>> Pulkit
>>>
>>

[julia-users] Re: if .. else on a whole array

[Gunnar Farnebäck]

If you're not concerned about performing unneeded calls to u2_0 you can get 
away with (x.>Lbox/2) .* (0.5*(u2_0(x)+u2_0(Lbox-x))-u2_
0(Lbox/2))


Den onsdag 10 februari 2016 kl. 11:58:58 UTC+1 skrev Ferran Mazzanti:
>
> Hi folks,
>
> probably a stupid question but can't find the answer, so please help if 
> you can :)
> I would like to evaluate a if.. else.. statement on a whole array. 
> Actually it's a bit more complicated, as I have a function that
> previously was
>
> function u2(x)
> return 0.5*(u2_0(x)+u2_0(Lbox-x))-u2_0(Lbox/2)
> end;
>
> for some other defined function u2(x) and constant L_box. The thing is 
> that I could directly evaluate that on a scalar x and on an array.
> Now I have to change it and check if x is smaller than Lbox/2, returning 
> the same as above if it is, or 0 otherwise.
> I tried something of the form 
>
> function u2(x)
> return x.>Lbox/2 ? 0 : 0.5*(u2_0(x)+u2_0(Lbox-x))-u2_0(Lbox/2)
> end;
>
> but that complains when x is an array. What would be the easiest way to 
> achieve this? It should work for both scalars and arrays...
>
> Best regards and thanks,
>
> Ferran.
>
>

[julia-users] Re: Google releases TensorFlow as open source

[Gunnar Farnebäck]

Den torsdag 12 november 2015 kl. 06:36:28 UTC+1 skrev Alireza Nejati
>
> Anyway, the problem I'm facing right now is that even though TensorFlow's 
> python interface works fine, I can't get TensorFlow's C library to build! 
> Has anyone else had any luck with this? I've had to update java AND gcc 
> just to make some progress in building (they use c++11 features, don't 
> ask). Plus I had to install google's own bizarre and buggy build manager 
> (bazel). TensorFlow.jl would be kind of pointless if everyone faced the 
> same build issues...
>

I managed to build from source on an Ubuntu machine at work, following the 
instructions at 
https://github.com/tensorflow/tensorflow/blob/master/tensorflow/g3doc/get_started/os_setup.md.
 
The most difficult part was getting the Python package installation 
dependencies in working order, which was not covered by the instructions. 
(The low mark being "pip install --no-use-wheel --upgrade distribute", 
which a google search suggested to get past the final hurdle and actually 
did work. Please don't emulate this in Julia.)

Whether it actually built anything useful beyond what the Python package 
needed I have no idea though. There's a forest of bazel generated 
directories that I'm not particularly tempted to try to navigate unless 
there's something specific I should look into.

[julia-users] Re: obtain T from Ptr{T} at runtime

[Gunnar Farnebäck]

julia> eltype(Ptr{Float64})
Float64

Den tisdag 3 november 2015 kl. 07:51:43 UTC+1 skrev Petr Hlavenka:
>
> Hi,
> is there a way to obtain a type T from a Ptr{T}, e.g. get UInt32 from 
> Ptr{UInt32}.
>
> I'm trying to fix NIDAQ for 0.4 (tupplecalypse)
>
> This is needed if I have a signature of a function (demonstration example 
> from base, in real NIDAQ package there is huge C header processed by 
> Clang.jl)
>
> cfunction = getfield(Base.dSFMT,:dsfmt_fill_array_close1_open2!)
> signature = cfunction.env.defs.sig
> pt = signature.types[2]
> out> Ptr{Float64}
>
> how to obtain the Float64 base type from Ptr{T} to construct the 
> appropriate array to pass into the function?
>
> Petr
>
>
>

[julia-users] Re: Code runs 500 times slower under 0.4.0

[Gunnar Farnebäck]

If you don't have deprecation warnings I would suspect some change in 0.4 
has introduced type instabilities. If you are using typed concatenations 
you could be hit by https://github.com/JuliaLang/julia/issues/13254.

Den torsdag 22 oktober 2015 kl. 23:03:00 UTC+2 skrev Kris De Meyer:
>
> Are there any general style guidelines for moving code from 0.3.11 to 
> 0.4.0? Running the unit and functionality tests for a module that I 
> developed under 0.3.11 in 0.4, I experience a 500 times slowdown of blocks 
> of code that I time with @time. 
>
> Can't even imagine where I have to start looking, and find it 
> flabbergasting that perfectly valid julia code under 0.3.11 (not generating 
> a single warning) can show such a performance degradation under 0.4.0.
>
> Anyone seen anything similar? Is there some fundamental difference in how 
> code is JIT-compiled under 0.4.0?
>
> Thanks,
>
> Kris
>
>
>
>

Re: [julia-users] Is there a direct way to get binary representation of a float?

[Gunnar Farnebäck]

It might help if you explain something about what you want to do with the 
binary representation.

For example,

julia> f(x::Float64) = reinterpret(Float64, reinterpret(UInt64, x) & 
0x7fff)
f (generic function with 1 method)

generates very simple LLVM code

julia> code_llvm(f, (Float64,))

define double @julia_f_21398(double) {
top:   
  %1 = bitcast double %0 to i64
  %2 = and i64 %1, 9223372036854775807
  %3 = bitcast i64 %2 to double
  ret double %3
}

which is compiled by LLVM into

julia> code_native(f, (Float64,))
.text
Filename: none
Source line: 0
pushq   %rbp
movq%rsp, %rbp
Source line: 1
vmovq   %xmm0, %rax
movabsq $9223372036854775807, %rcx # imm = 0x7FFF
andq%rax, %rcx
vmovq   %rcx, %xmm0
popq%rbp
retq
nopl(%rax)
pushq   %rbp
movq%rsp, %rbp
subq$16, %rsp
movq(%rsi), %rax
vmovsd  (%rax), %xmm0
movabsq $f, %rax
callq   *%rax
vmovsd  %xmm0, -8(%rbp)
movabsq $__pool_alloc, %rax
callq   *%rax
movabsq $140629191800336, %rcx  # imm = 0x7FE6C905B610
movq%rcx, -8(%rax)
vmovsd  -8(%rbp), %xmm0
vmovsd  %xmm0, (%rax)
addq$16, %rsp
popq%rbp
retq

The lines up to the first retq look reasonable to me but I have no idea 
what's up with the rest of the output. Possibly a bug in code_native.

Den måndag 19 oktober 2015 kl. 06:59:29 UTC+2 skrev Uliano Guerrini:

> I see, however my question is if there is a way to perform the conversion 
> shown (or the like) without any (or minimal) run time penalty that is the 
> reason why when I found that reinterpret was calling intrinsics I focussed 
> on them. In C time penalty is 0 If Julia proposes itself as solving the two 
> language problem I should not think to write parts of the code in C 
>
>

Re: [julia-users] Is there a direct way to get binary representation of a float?

[Gunnar Farnebäck]

Den måndag 19 oktober 2015 kl. 10:28:15 UTC+2 skrev Uliano Guerrini:
>
>
>
> Il giorno lunedì 19 ottobre 2015 08:56:32 UTC+2, Gunnar Farnebäck ha 
> scritto:
>>
>> It might help if you explain something about what you want to do with the 
>> binary representation.
>>
>
> I would like to have an efficient way to access bytes, 16 and 32 bit words 
> of larger types (in this case 8 bytes, not necessarily floating point). 
> Possibly as efficient and as easy as it could be in C/C++
>
> struct Words32 {
>
> int32_t lo;
>
> int32_t hi;
>
> };
>
> union Both {
>
> struct Words32 w;
>
> int64_t i;
>
> };
>
> union Both b;
>
> b.i=0x7fff3fff;
>
> printf("lo:0x%08x\nhi:0x%08x\n",b.w.lo,b.w.hi);
>
> gives
>
> *lo:0x3fff*
>
> *hi:0x7fff*
>
> *type Words32*
>
>   *hi::Int32*
>
>   *lo::Int32*
>
> *end*
>
> *a=0x7fff3fff*
>
> *reinterpret(Words32,a)*
>
>
> *gives*
>
>
> *ERROR: reinterpret: expected bits type as first argument*
>
> * in reinterpret at essentials.jl:115*
>
>
>  
>
>>
>> For example,
>>
>> julia> f(x::Float64) = reinterpret(Float64, reinterpret(UInt64, x) & 
>> 0x7fff)
>> f (generic function with 1 method)
>>
>
>
> I have seen that I can reinterpret a bits type into another bits type of 
> the SAME size. unfortunately my Words32 isn't, even if it size is 
> completely specified and should be known ahead of time by the compiler. I'm 
> looking for ways to access its smaller chunks without passing through a 
> conversion into an hexadecimal string.
>
 
Is this good enough?

julia> immutable Words32
 hi::Int32
 lo::Int32
   end

julia> a = 0x7fff3fff
0x7fff3fff

julia> reinterpret(Words32, [a])
1-element Array{Words32,1}:
 Words32(1073741823,2147483647)

[julia-users] Re: non-modifying push (and unshift, etc.)

[Gunnar Farnebäck]

Array{T, N} is a non-concrete type and thus an inference failure when a 
concrete type could have been inferred. Most of the concatenation inference 
failures existed in 0.3 but there are also some regressions, most 
annoyingly for some constructions referred to by deprecation messages.

Den måndag 19 oktober 2015 kl. 17:17:29 UTC+2 skrev Andras Niedermayer:
>
> Yes, https://github.com/JuliaLang/julia/issues/13254 and 
> https://github.com/JuliaLang/julia/issues/13665 is what I meant. I'm not 
> sure any more whether "type instability" is the right expression, since the 
> type inferencer indeed consistently gives you Array{T,N} for [[1,2];3] for 
> all values. But you'd get more efficient code if the type inferencer gave 
> you Array{Int64,1} instead of Array{T,N} for vcat.
>
> It doesn't look like a regression to me, I get this already on Julia 0.3:
>
> julia> g() = [[1,2];3]
> g (generic function with 1 method)
>
>
> julia> @code_typed g()
> 1-element Array{Any,1}:
>  :($(Expr(:lambda, {}, {{:_var0,:_var1},{{:_var0,(Array{Int64,1},Int64),0
> },{:_var1,Array{Int64,1},18}},{}}, :(begin  # none, line 1:
> _var1 = vcat(1,2)::Array{Int64,1}
> return cat(1,_var1::Array{Int64,1},3)::Array{T,N}
> end::Array{T,N}
>
>
>
>
> On Monday, October 19, 2015 at 5:02:57 PM UTC+2, Kristoffer Carlsson wrote:
>>
>> [a;b] is a bit regressed now when it comes to type stability, see: 
>> https://github.com/JuliaLang/julia/issues/13254
>>
>> On Monday, October 19, 2015 at 4:56:40 PM UTC+2, Josh Langsfeld wrote:
>>>
>>> [a; b] / vcat(...) is type stable. It may produce different output 
>>> depending on the types of a and b but it won't change behavior depending on 
>>> their values.
>>>
>>> On Monday, October 19, 2015 at 10:20:06 AM UTC-4, Andras Niedermayer 
>>> wrote:

 In light of the recent discussions (
 https://groups.google.com/forum/#!topic/julia-users/xJ7GpKAa16E and 
 https://groups.google.com/forum/#!topic/julia-users/_lIVpV0e_WI) I got 
 curious, whether there is a non-modifying version of push!, since 
 push!(copy(a),b) doesn't feel right (I try to avoid modifying functions, 
 except if I really need performance) and [a; b] is not type stable:

 f(a) = [a; 3]
 f2(a) = push!(copy(a),3)


 julia> Base.return_types(f,(Array{Int64,1},))
 1-element Array{Any,1}:
  Array{T,N}


 julia> Base.return_types(f2,(Array{Int64,1},))
 1-element Array{Any,1}:
  Array{Int64,1}


 I couldn't find anything in the Julia docs. Of course, I could just 
 define my own

 push(a,vars...) = push!(copy(a),vars...)


 but if there is a standard way to do that, I'd prefer that. Maybe 
 there's also some clever way to avoid making a copy if `a` is a literal 
 (e.g. 
 push([1,2])
 ).

 The same applies to unshift, etc.

 (Non-modifying `shift!` and `pop!` already exist with the names `first` 
 and `last`. I'd find it easier to remember `shift` and `pop` -- it would 
 also be more similar to `merge` vs `merge!` for dictionaries, but I guess 
 that's just my taste.)

>>>

[julia-users] Re: Syntax for slicing with an array of ranges?

[Gunnar Farnebäck]

julia [foo[r] for r in some_ranges]
2-element Array{Any,1}:
 [6,2]
 [4,9]

julia vcat([foo[r] for r in some_ranges]...)
4-element Array{Int64,1}:
 6
 2
 4
 9


Den onsdag 19 augusti 2015 kl. 04:53:55 UTC+2 skrev John Brock:

 Is there an easy way to slice an array using an array of ranges? I'm 
 looking for something like the following:

 foo = [2,6,2,8,4,9]
 some_ranges = UnitRange{Int64}[2:3, 5:6]
 foo[some_ranges] # gives error; desired output is [6,2,4,9]

[julia-users] Re: Convert Array{Tuple} to Matrix

[Gunnar Farnebäck]

This is reasonably clean:

[y[i] for y in x, i in 1:3]

If performance matters you are probably better off with the other 
suggestions.

Den måndag 6 juli 2015 kl. 22:10:01 UTC+2 skrev Júlio Hoffimann:

 Hi,

 How to convert:

 1000-element Array{Tuple{Integer,Integer,Integer},1}:
  (10,2,1) 
  (5,7,10) 
  (5,7,4)  
  (1,1,6)  
  (2,3,6)  
  (8,6,4)  
  (10,2,4) 
  (1,3,9)  
  (9,3,7)  
  (5,2,4)  
  ⋮
  (1,6,8)  
  (4,6,6)  
  (3,9,5)  
  (10,4,10)
  (8,7,4)  
  (4,8,9)  
  (2,6,10) 
  (3,6,5)  
  (1,7,10) 

 into the corresponding 1000x3 matrix in a clean fashion?

 -Júlio

[julia-users] Re: Efficient way to split an array/dataframe?

[Gunnar Farnebäck]

It doesn't behave the same with respect to duplicates in A but I would do

N = 100
A = rand(N)
n = iceil(0.7 * N)
testindex = sample(1:size(A,1), replace=false, n)
testA = A[testindex]
trainindex = setdiff(1:N, testindex)
trainA = A[trainindex]

Den torsdag 26 mars 2015 kl. 06:21:42 UTC+1 skrev verylu...@gmail.com:

 Hi,
 I have an array of 100 elements. I want to split the array to 70 (test 
 set) and 30 (train set) randomly.

 N=100
 A = rand(N);
 n = convert(Int, ceil(N*0.7))
 testindex = sample(1:size(A,1), replace=false,n)
 testA = A[testindex];

 How can I get the train set?

 I could loop through testA and A to get trainA as below

 trainA = Array(eltype(testA), N-n);
 k=1
 for elem in A
 if !(elem in testA)
 trainA[k] = elem
 k=k+1
 end
 end

 Is there a more efficient or elegant way to do this?

 Thanks!

Re: [julia-users] Some simple use cases for multi-threading

[Gunnar Farnebäck]

Den fredag 20 mars 2015 kl. 05:36:15 UTC+1 skrev Kiran Pamnany:

 Hey Stefan,
  

 Depending on how well Julia-generated code performs, some important 
 concurrent data structures might be better implemented in C and wrapped.


 I would argue that until we can implement efficient concurrent data 
 structures in pure Julia, the work isn't complete. It's ok to use C code as 
 an interim solution while we're developing this, but in the longer term, we 
 should not be forced to rely on C to write libraries.


 Think of the implementations in C as though we're dropping into inline 
 assembly for maximizing performance. In actual fact, using alignment 
 directives, padding, SIMD intrinsics, or specialized load/store 
 instructions _is_ pretty much assembly. I don't see value in adding this 
 sort of architecture-specific explicit control to Julia--it negates the 
 high-level, high-productivity aspect if you have to throw in pragmas like 
 assume_aligned in your code for performance reasons. And if you go down 
 that road, then you're competing with C.

 On the contrary I see an enormous value in being able to do what's 
effectively inline assembly within an environment with the metaprogramming 
capabilities of Julia and seamless integration with the parts of the code 
that are not quite as performance critical.

[julia-users] Re: How to quickly count the moving average with one vector and save in the other vec? Whitout while

[Gunnar Farnebäck]

If speed is what's important you almost certainly want some kind of loop.

If it's more important not to have a loop you can satisfy your 
specification with

b = reshape(repmat(mean(reshape(a,10,10),1),10),100)

although that's not what I would call a moving average.

Den tisdag 10 februari 2015 kl. 10:00:36 UTC+1 skrev paul analyst:

 How to quickly count the moving average with one vector and save the other 
 vector?
 I have vec a and need to fill vec b moving average of a, in this same range.
 Is posible to do it whitout while ?

 julia a=rand(100);
 julia b=zeros(100);

 julia b[1:10]=mean(a[1:10])
 0.6312220153427996

 julia b[11:20]=mean(a[11:20])
 0.6356771620528772

 julia b
 100-element Array{Float64,1}:
  0.631222
  0.631222
  0.631222
  0.631222
  0.631222
  0.631222
  0.631222
  0.631222
  0.631222
  0.631222
  0.635677
  0.635677
  0.635677
  0.635677
  0.635677
  0.635677
  0.635677
  0.635677
  0.635677
  ?
  0.0


 Paul
  0.0

[julia-users] Re: Please help, simple double while

[Gunnar Farnebäck]

Syntactically you're still missing a comma between the loop variables. 
Semantically that's not how the double for loop works.

Try either

julia for i = 1:length(x)
   println(x[i], y[i])
   end

or

julia for (i, j) in zip(x, y)
   println(i, j)
   end


Den tisdag 27 januari 2015 kl. 11:46:03 UTC+1 skrev paul analyst:

 ;) it is 100 steps.

 I need 10 steps driven by 2 vectors x and y

 julia l=193
 193

 julia x=int(rand(10)*l)
 10-element Array{Int64,1}:
  126
   70
  132
   51
   77
  138
   94
  150
  142
   93

 julia y=int(rand(10)*l)
 10-element Array{Int64,1}:
  170
  111
   93
  126
   77
  135
   24
  182
   45
  179

 julia for i in [x] j in [y]
println(x,y)
end
 ERROR: j not defined
  in anonymous at no file:1

 What wrong here ?
 Pau;

Re: [julia-users] simplex in Julia is very slow

[Gunnar Farnebäck]

It does have issues. I tracked down a part of the story to 
https://github.com/JuliaLang/julia/pull/9086

Den onsdagen den 26:e november 2014 kl. 20:52:57 UTC+1 skrev Stefan 
Karpinski:

 I wonder if \ might be type unstable.

[julia-users] Re: ternary operator with return() gives (sometimes) errors

[Gunnar Farnebäck]

Most likely Julia tries to parse it as x  10 ? return ((0):nothing)

Try x  10 ? (return 0) : nothing

instead.

Or the more popular construction x = 10 || return 0

Den tisdagen den 25:e november 2014 kl. 11:21:39 UTC+1 skrev andreas:


 Dear list,

 I observed that this

 function foo(x)
 x = 10 ? nothing : return(0)
 y = 2*x
 return(y)
 end

 works as expected, while that

 function bar(x)
 x  10 ? return(0) : nothing
 y = 2*x
 return(y)
 end

 throws an error: ERROR: syntax: colon expected in ? expression.

 Is this an intended behavior?

 Cheers,
 Andreas

 PS: I'm using Julia v 0.3.3 on Windows 7

[julia-users] Re: Getting Julia to emit an unrolled inner loop with SIMD instructions

[Gunnar Farnebäck]

It might be possible to get further doing llvmcall with LLVM inline 
assembly but my feeble attempts at that have only resulted in crashes 
(well, sort of, getting dumped back to shell prompt without ceremony). 
Getting good load/store instructions by default would be a lot better of 
course but I would still be interested in knowing whether llvmcall with 
inline assembly can be made to work, in order to access a wider range of 
SIMD instructions.

Den söndagen den 2:e november 2014 kl. 09:28:58 UTC+1 skrev Toivo 
Henningsson:

 I was able to get further by using llvmcall in Julia 0.4 as per 
 https://github.com/JuliaLang/julia/issues/8786#issuecomment-60511713 
 (thanks, @GunnarFarneback!). By manually unrolling the loop in can get an 
 unrolled loop with SIMD instructions.

 But llvm creates a lot pack/unpack instructions to only load/store 64 bits 
 from memory at a time. (See this gist 
 https://gist.github.com/toivoh/436cb54b6d23cdfb9888.) Arrays in Julia 
 seem to be always 16 byte aligned. Is it possible to make llvm generate 
 code for the 16-byte aligned case?

[julia-users] Re: Ported a CSPRNG (Isaac64) from C, run times ~100x slower

[Gunnar Farnebäck]

I leave it to someone with more insights to comment on the efficiency of 
nested functions, but this seems at least to be substantially faster:

type State
a::Uint64
b::Uint64
aa::Uint64
bb::Uint64
cc::Uint64
RANDSIZL::Int
RANDSIZ::Int
end

function ind(s::State, mm, x)
return mm[1 + (x  ((s.RANDSIZ-1)3))3]
end

function rngstep(s::State, randrsl, mm, a_mix, m1, m2, r, i)
x = mm[m1 + i]
s.a = a_mix + mm[m2 + i]
mm[m1 + i] = y = ind(s, mm, x) + s.a + s.b
randrsl[r + i] = s.b = ind(s, mm, y  s.RANDSIZL) + x
end

function isaac64(s::State, randrsl, mm)
s.a = s.aa
s.cc += 1
s.b = s.bb + s.cc
m1 = r = 0
mid = m2 = s.RANDSIZ  1
for i = 1:4:mid-2
rngstep(s, randrsl, mm, ~(s.a $ (s.a  21)), m1, m2, r, i)
rngstep(s, randrsl, mm,   s.a $ (s.a  5),   m1, m2, r, i + 1)
rngstep(s, randrsl, mm,   s.a $ (s.a  12),  m1, m2, r, i + 2)
rngstep(s, randrsl, mm,   s.a $ (s.a  33),  m1, m2, r, i + 3)
end
m2 = 0
r = m1 = s.RANDSIZ  1
for i = 1:4:mid-2
rngstep(s, randrsl, mm, ~(s.a $ (s.a  21)), m1, m2, r, i)
rngstep(s, randrsl, mm,   s.a $ (s.a  5),   m1, m2, r, i + 1)
rngstep(s, randrsl, mm,   s.a $ (s.a  12),  m1, m2, r, i + 2)
rngstep(s, randrsl, mm,   s.a $ (s.a  33),  m1, m2, r, i + 3)
end
s.bb = s.b
s.aa = s.a;
end

macro mix()
quote
a -= e; f $= h  9;  h += a
b -= f; g $= a  9;  a += b
c -= g; h $= b  23; b += c
d -= h; a $= c  15; c += d
e -= a; b $= d  14; d += e
f -= b; c $= e  20; e += f
g -= c; d $= f  17; f += g
h -= d; e $= g  14; g += h
end
end

function randinit(randrsl, mm, RANDSIZL, RANDSIZ, flag=false)
a = b = c = d = e = f = g = h = 0x9e3779b97f4a7c13
for i = 1:4
@mix
end
for i = 1:8:RANDSIZ
if flag
a += randrsl[i  ]
b += randrsl[i+1]
c += randrsl[i+2]
d += randrsl[i+3]
e += randrsl[i+4]
f += randrsl[i+5]
g += randrsl[i+6]
h += randrsl[i+7]
end
@mix
mm[i  ] = a
mm[i+1] = b
mm[i+2] = c
mm[i+3] = d
mm[i+4] = e
mm[i+5] = f
mm[i+6] = g
mm[i+7] = h
end
if flag
for i = 1:8:RANDSIZ
a += mm[i  ]
b += mm[i+1]
c += mm[i+2]
d += mm[i+3]
e += mm[i+4]
f += mm[i+5]
g += mm[i+6]
h += mm[i+7]
@mix
mm[i  ] = a
mm[i+1] = b
mm[i+2] = c
mm[i+3] = d
mm[i+4] = e
mm[i+5] = f
mm[i+6] = g
mm[i+7] = h
end
end
z = 0x
s2 = State(z, z, z, z, z, RANDSIZL, RANDSIZ)
isaac64(s2, randrsl, mm)
return s2
end

function main(randrsl::Vector{Uint64}=Uint64[])
RANDSIZL = 8
RANDSIZ = length(randrsl)
if RANDSIZ == 0
RANDSIZ = 256
resize!(randrsl, RANDSIZ)
fill!(randrsl, uint64(0))
elseif RANDSIZ % 8 != 0
error(dimension of seeding array must be a factor of 8)
end
mm = zeros(Uint64, RANDSIZ)
s = randinit(randrsl, mm, RANDSIZL, RANDSIZ, true)
for i = 1:2
isaac64(s, randrsl, mm)
end
end

function timeit(f, n=10^4)
f()
t = time()
for i = 1:n
f()
end
@printf(avg running time: %llf \n, ((time() -t) / n ))
end


Den fredagen den 17:e oktober 2014 kl. 22:58:35 UTC+2 skrev alexander 
maznev:

 I was trying out Julia for the first time and did not seem to find a 
 CSPRNG available. Below is a port of Isaac64 (
 http://burtleburtle.net/bob/rand/isaacafa.html). My runtime meassurements 
 seem to be very poor - reading a little bit about Julia prior to trying the 
 language, and especially looking at the benchmarks I was under the 
 impression that run times should be relatively comparable with C. Would 
 appreciate some feedback, thanks. 


 function main(randrsl::Array{Uint64}=Uint64[]) #takes Uint64 array as a 
 seed
 #wrap everything to pass variables into sub-scopes implicitly
 function isaac64()
 function ind(mm,x)#nesting ind() into rngstep 
 makes it another 2-3 times slower -!
 return mm[int((x  ((RANDSIZ-1)3))/8)+1]
 end
 function rngstep(a_mix,m1,m2,r)
 x = mm[m1+i]
 a = (a_mix) + mm[m2+i]
 mm[m1+i] = y = (ind(mm, x) + a + b)
 randrsl[r+i] = b = (ind(mm, yRANDSIZL) + x)
 i+=1
 end
 a=aa; cc+=1; b=bb+cc; x=0x00::Uint64; 
 m1=r=0;mid=m2=int(RANDSIZ/2); i=1;
 while (i  (mid-2)) #do it without pointers
 rngstep(~(a$(a21)), m1, m2, r) #sub-indexing an array 
 creates a new object, pass indexes as variables
 rngstep(a$(a5), m1, m2, r)
 rngstep(a$(a12), m1, m2, r)
  

[julia-users] Re: Loading data just once

[Gunnar Farnebäck]

I wrote this when I wanted to cache the results of matread. Your problem 
sounds like it could be similar.

let matread_cache = (String = Any)[]
global caching_matread
function caching_matread(filename)
if !haskey(matread_cache, filename)
matread_cache[filename] = matread(filename)
end
return matread_cache[filename]
end
end


Den torsdagen den 9:e oktober 2014 kl. 10:08:23 UTC+2 skrev 
cormu...@mac.com:

 A beginner's question...

 I'm writing a function that wants to load a set of data from a file, 
 depending on an argument passed to the function (so a different argument 
 requires a different set of data to be loaded). I'd like each set of data 
 to be stored somehow in separate variables so that when the function is 
 called again with the same argument, it doesn't have to load that 
 particular data again (because it takes 5 seconds to load).

 I'm not sure whether this requires the use of global variables? I looked 
 through the documents (
 http://julia.readthedocs.org/en/latest/search/?q=global) but didn't gain 
 enlightenment. :)

 I think I can test for the existence of a previously-defined variable 
 using:

 if !isdefined(symbol(string(dataset)))
dataset = include($(dataset).jl)
 end

 but I'm not convinced this works correctly, because this creates a 
 variable inside the function...

Re: [julia-users] Re: Article on `@simd`

[Gunnar Farnebäck]

There are still three arguments to max in the last of those examples. 
Actually it's not clear that you can make an equivalent expression with min 
and max. Functionally (with intended use)
x[i] = max(a, min(b, x[i]))
does the same thing as the earlier examples but it expands to
x[i] = ifelse(ifelse(b  x[i], b, x[i])  a, a, ifelse(b  x[i], b, x[i]))
which should be hard for a compiler to optimize to the earlier examples 
since they don't give the same result in the degenerate case of a  b.

A closer correspondence is given by the clamp function which is implemented 
as a nested ifelse in the same way as example 2 (although in the opposite 
order, so it also differs for ab).

Den onsdagen den 17:e september 2014 kl. 16:28:45 UTC+2 skrev Arch Robison:

 Thanks.  Now fixed.

 On Wed, Sep 17, 2014 at 4:14 AM, Gunnar Farnebäck gun...@lysator.liu.se 
 javascript: wrote:

 In the section The Loop Body Should Be Straight-Line Code, the first 
 and second code example look identical with ifelse constructions. I assume 
 the first one should use ? instead. Also the third code example has a stray 
 x[i]a argument to the max function.

Re: [julia-users] unexpected domain error for ^(float,float)

[Gunnar Farnebäck]

It's not like Julia is doing anything strange or uncommon here. Most people 
would be really surprised if -10² meant positive 100.

Den torsdagen den 18:e september 2014 kl. 15:01:44 UTC+2 skrev Jutho:

 because it is not recognized/parsed as literal but as the application of a 
 unary minus, which has lower precedence than ^

 I guess it is not possible to give binary minus a lower precedence than ^ 
 and unary minus of higher precedence, since these are just different 
 methods of the same function/operator.

 Op donderdag 18 september 2014 14:54:26 UTC+2 schreef Florian Oswald:

 yes - not sure why -0.4 and (-0.4) are any different.

 On 18 September 2014 13:52, Patrick O'Leary patrick...@gmail.com wrote:

 Seems like the literal -0.4^2.5 should throw the same error, though?


 On Thursday, September 18, 2014 6:42:56 AM UTC-5, Tim Holy wrote:

 http://docs.julialang.org/en/latest/manual/faq/#why-does-
 julia-give-a-domainerror-for-certain-seemingly-sensible-operations 

 On Thursday, September 18, 2014 03:24:00 AM Florian Oswald wrote: 
  # define a variable gamma: 
  
  gamma = 1.4 
  mgamma = 1.0-gamma 
  
  julia mgamma 
  -0.3999 
  
  # this works: 
  
  julia -0.3999^2.5 
  -0.10119288512475567 
  
  # this doesn't: 
  
  julia mgamma^2.5 
  ERROR: DomainError 
  in ^ at math.jl:252 

[julia-users] Re: Article on `@simd`

[Gunnar Farnebäck]

In the section The Loop Body Should Be Straight-Line Code, the first and 
second code example look identical with ifelse constructions. I assume the 
first one should use ? instead. Also the third code example has a stray 
x[i]a argument to the max function.

Den måndagen den 15:e september 2014 kl. 23:39:20 UTC+2 skrev Arch Robison:

 I've posted an article on the @simd feature to 
 https://software.intel.com/en-us/articles/vectorization-in-julia .   
 @simd is an experimental feature 
 http://julia.readthedocs.org/en/release-0.3/manual/performance-tips/#performance-annotations
  
 in Julia 0.3 that gives the compiler more latitude to .vectorize loops.   
 Corrections/suggestions appreciated.

 - Arch D. Robison
   Intel Corporation

Re: [julia-users] slow julia version of c code

[Gunnar Farnebäck]

This function helped me find a couple of type problems in a somewhat 
complex code base. The idea is to first run your code, then use this 
function to dig up the type signatures that were actually called and 
identify type problematic variables from code_typed.

I strongly suspect that the inspection of f.env.cache_arg1 is not the 
proper way to find the used type signatures, so if someone who understands 
the internals could clean this up I think it could become a valuable tool.

function inspect(modulename, filename)
for name in names(modulename)
f = getfield(modulename, name)
if isgeneric(f)
if basename(string(f.env.defs.func.code.file)) == filename
println(f)
c = f.env.cache_arg1
for k = 1:length(c)
if isdefined(c, k)
sig = c[k].sig
println(  , sig)
t = code_typed(f, sig)
vars = t[1].args[2][2]
for var in vars
if !isleaftype(var[2])
println($(var[1]) $(var[2]))
end
end
end
end
end
end
end
end




Den tisdagen den 16:e september 2014 kl. 06:02:54 UTC+2 skrev Simon 
Kornblith:

 There is most certainly a type problem. You're not getting type 
 information for sparse_grid.lvl_l which deoptimizes a lot of things. In 
 your code, you have:

 type sparse_gridd::Int64q::Int64n::Int64
 grid::Array{Float64}ind::Array{Int64}lvl::Array{Int64}   
  lvl_l::Array{Int64}lvl_s::Array{Float64}bounds::Array
 {Float64}gtypeBBfubend

 Try specifying the number of dimensions for the arrays, e.g. 
 grid::Array{Float64,2}. Also using any of the untyped fields will slow 
 code down, but I don't see them in that function.

 A hackish way to find type issues is:

 filter!(x-!isleaftype(x[2]), @code_typed(f(args))[1].args[2][2])

 which returns a list of variables where type inference couldn't infer a 
 concrete/leaf type. In some cases these will be hidden variables so you 
 may still need to inspect the output of code_typed. TypeCheck.jl and 
 Lint.jl may also provide more sophisticated tools for this purpose.

 Simon

 On Monday, September 15, 2014 5:16:36 AM UTC-4, Zac wrote:

 https://gist.github.com/Zac12345/519bd7a503a1fd1b8d98 has the updated 
 function and code_typed output

 Pkg.clone(https://github.com/Zac12345/Sparse;) should be fine for the 
 test as the gist does'nt require functions using the shared lib.

 Are there any specific things i should be looking for in the code_typed 
 results?

Re: [julia-users] Re: How quickly calculate the entropy for long vectors ? -sum(x-x*log(2,x), [count(x-x==c,s)/length(s) for c in unique(s)])

[Gunnar Farnebäck]

You get NaN because the binning has left some bin empty and x * log(2, x) 
for x = 0.0 results in 0.0 * -Inf which is NaN. You can avoid this 
particular problem by replacing x * log(2,x) with something that special 
cases x = 0.0 to return 0.0 rather than NaN.

Note that the result you get from this operation will depend on the number 
of bins and will typically be different from the result of your first 
approach.

Den fredagen den 5:e september 2014 kl. 09:16:18 UTC+2 skrev paul analyst:

  THX, work but:
 julia t=h5read(FMLo2_z_reversem.h5,FMLo2,(:,1));

 julia eltype(t)
 Float64

 julia t
 5932868x1 Array{Float64,2}:
   0.0181719
   0.303473
  -0.526979
   ?
  -0.526979
   0.912295
  -0.0281875

 julia entropy(t)
 NaN

 julia entropy(vec(t))
 NaN
 Why ?
 julia s

 100-element Array{Float64,1}:
  0.737228
  0.162308
  0.688503
  0.108727
  0.40552
  0.654883
  0.618194
  0.137147
  0.934373
  0.077236
  ?
  0.413675
  0.463914
  0.504321
  0.976408
  0.998662
  0.343927
  0.477739
  0.660533
  0.918326
  0.579264

 julia entropy(s)
 13.280438296634793


 julia

 julia entropy(t[:,1])
 NaN

 ???
 Paul


 W dniu 2014-09-05 07:34, David P. Sanders pisze:
  


 El jueves, 4 de septiembre de 2014 23:42:20 UTC-5, paul analyst escribió: 

 julia entropy(s)=-sum(x-x*log(2,x), [count(x-x==c,s)/length(s) for c 
 in unique(s)]);

 julia s=rand(10^3);

 julia @time entropy(s)
 elapsed time: 0.167097546 seconds (20255140 bytes allocated)
 9.965784284662059

 julia s=rand(10^4);

 julia @time entropy(s)
 elapsed time: 3.62008077 seconds (1602061320 bytes allocated, 21.81% gc 
 time)
 13.287712379549843

 julia s=rand(10^5);

 julia @time entropy(s)
 elapsed time: 366.181311932 seconds (160021245832 bytes allocated, 21.89% 
 gc time)
 16.609640474434073

  
  You can see from these last two that the time is multiplied by 100 when 
 the length of the vector is multiplied by 10, i.e. your method has O(n^2) 
 complexity. This is due to the way that you are counting repeats. What you 
 are basically doing is a histogram. 

  If your data are really floats, then in any case some binning is 
 required. If they are ints, you could use a dictionary. I think there's 
 even a counting object already implemented (but I don't remember what it's 
 called).

  How about this:

  function entropy(s)
N = length(s)
num_bins = 1
h = hist(s, num_bins)
p = h[2] ./ N  # probabilities
-sum([x * log(2,x) for x in p])
 end
  
  julia @time entropy(rand(10^6))
 elapsed time: 0.199634039 seconds (79424624 bytes allocated, 31.51% gc 
 time)
 13.28044771568381

  julia @time entropy(rand(10^7))
 elapsed time: 1.710673571 seconds (792084208 bytes allocated, 26.20% gc 
 time)
 13.286992511965552

  julia @time entropy(rand(10^8))
 elapsed time: 18.20088571 seconds (7918627344 bytes allocated, 24.03% gc 
 time)
 13.28764216804997
  
  The calculation is now O(n) instead.

   

 julia s=rand(10^6);

 julia @time entropy(s)
 
 After 12 h not yet counted :/

 Paul
  
  
 

[julia-users] Re: Best way to convert Any of vectors into matrix

[Gunnar Farnebäck]

For maximum efficiency Diego's loop solution is likely best but for 
succinctness, try hcat(A...).

Den lördagen den 23:e augusti 2014 kl. 04:05:19 UTC+2 skrev Bradley Setzler:

 Good evening,

 I often have Any types filled with vectors (especially as the return of a 
 pmap), and need them as a matrix or DataFrame. For example, suppose I have,

 julia A
 3-element Array{Any,1}: 
 [1,2] 
 [3,4] 
 [5,6]

 But I want,

 julia B
 2x3 Array{Int64,2}: 
 1 3 5 
 2 4 6

 I came up with the following, which successfully constructs B from A, but 
 it's a bit inefficient/messy:

 B = A[1]
 for i=2:length(A)
 B = hcat(B,A[i])
 end

 Is there a better way to do this? Something like,
 julia B = hcatForAny(A)

 Thanks,
 Bradley

Re: [julia-users] Re: question about array comparison

[Gunnar Farnebäck]

Den fredagen den 8:e augusti 2014 kl. 22:25:01 UTC+2 skrev Steven G. 
Johnson:



 On Friday, August 8, 2014 3:19:00 AM UTC-4, Gunnar Farnebäck wrote:

 1. As in Ethan's response (modulo a trivial bug):

 Matrix[x[k,:] for  k=1:size(x,1)]

 Depending on how well typed the result needs to be, the leading 'Matrix' 
 can be dropped or be made more specific.


 You could also just do

  reshape(x.', length(x))
  

That produces something rather different and I'm not sure how it matches a 
1-dimensional array of row arrays but maybe it's some Pythonic terminology 
I'm not familiar with.

Re: [julia-users] Re: question about array comparison

[Gunnar Farnebäck]

1. As in Ethan's response (modulo a trivial bug):

Matrix[x[k,:] for  k=1:size(x,1)]

Depending on how well typed the result needs to be, the leading 'Matrix' 
can be dropped or be made more specific.

2. I wouldn't say so but opinions vary depending on mathematical and 
programming language biases. There are several issues and mailing list 
threads on this subject.
You can start with https://github.com/JuliaLang/julia/issues/5949 and 
follow links or search further if that doesn't satiate your appetite.

Den fredagen den 8:e augusti 2014 kl. 01:31:23 UTC+2 skrev K leo:

 Thanks to both for the responses. 

 Related questions: how does one turn a 2-dimensional array into a 
 1-dimensional array of row arrays?  Also, the default behavior of julia 
 is that a row of a 2-dimensional array is also a 2-dimensional array.   
 Would anyone comment why this is the case? Should it better be a 
 1-dimensional array? 

 On 2014年08月07日 19:59, Ethan Anderes wrote: 
  
  Here is another one…slightly different: 
  
  |in(x[1,:], Matrix[x[k,:] for  k=1:size(x,2)] ) 
  | 
  
  To be honest, the custom for loop is probably your best bet in terms 
  of performance: 
  
  |function isrow(row, x::Matrix) 
  for k=1:size(x,1) 
 if row == x[k,:] 
return true 
 end 
  end 
  false 
  end 
  | 
  
  On Thursday, August 7, 2014 1:36:08 AM UTC-7, Gunnar Farnebäck wrote: 
  
  Is this easy enough? 
  
  julia x = reshape(1:16, 4, 4) 
  4x4 Array{Int64,2}: 
   1  5   9  13 
   2  6  10  14 
   3  7  11  15 
   4  8  12  16 
  
  julia a = x[1,:] 
  1x4 Array{Int64,2}: 
   1  5  9  13 
  
  julia [a == x[k,:] for k = 1:size(x,1)] 
  4-element Array{Any,1}: 
true 
   false 
   false 
   false 
  
  julia any([a == x[k,:] for k = 1:size(x,1)]) 
  true 
  
  julia find([a == x[k,:] for k = 1:size(x,1)]) 
  1-element Array{Int64,1}: 
   1 
  
  
  Den torsdagen den 7:e augusti 2014 kl. 06:09:23 UTC+2 skrev K leo: 
  
  julia x 
  4x4 Array{Int64,2}: 
1  5   9  13 
2  6  10  14 
3  7  11  15 
4  8  12  16 
  
  julia in(x[1,:], x) 
  false 
  
  julia x[1,:] 
  1x4 Array{Int64,2}: 
1  5  9  13 
  
  How can I check if x[1,:] is in x easily?  And with the row 
  index in x? 
  
  ​ 

[julia-users] Re: question about array comparison

[Gunnar Farnebäck]

Is this easy enough?

julia x = reshape(1:16, 4, 4)
4x4 Array{Int64,2}:
 1  5   9  13
 2  6  10  14
 3  7  11  15
 4  8  12  16

julia a = x[1,:]
1x4 Array{Int64,2}:
 1  5  9  13

julia [a == x[k,:] for k = 1:size(x,1)]
4-element Array{Any,1}:
  true
 false
 false
 false

julia any([a == x[k,:] for k = 1:size(x,1)])
true

julia find([a == x[k,:] for k = 1:size(x,1)])
1-element Array{Int64,1}:
 1


Den torsdagen den 7:e augusti 2014 kl. 06:09:23 UTC+2 skrev K leo:

 julia x 
 4x4 Array{Int64,2}: 
   1  5   9  13 
   2  6  10  14 
   3  7  11  15 
   4  8  12  16 

 julia in(x[1,:], x) 
 false 

 julia x[1,:] 
 1x4 Array{Int64,2}: 
   1  5  9  13 

 How can I check if x[1,:] is in x easily?  And with the row index in x? 

[julia-users] Re: Counting instances from an array that fulfill some inequality statements

[Gunnar Farnebäck]

This is not the most efficient way, neither the clearest, but it's compact. 
In a language like Matlab it would be the preferred approach.

sum((outputarray[1:end-1] . 0)  (outputarray[2:end] .= 0))

Den onsdagen den 30:e juli 2014 kl. 22:03:32 UTC+2 skrev yaois...@gmail.com:

 Hi guys,

 I asked this in a previous thread, but because that diverged off-topic 
 from my existing question, I decided to create a new thread.

 Anyhow, say I have an array

 outputarray = 
 Float64[-1.23423,-3.23423,-2.34234,-2.12342,1.23234,2.23423,-2.23432,5.2341,0.0,1.23423]

 This array lists the output of some function. I want to count the number 
 of times that the function passes by or equals 0 while emerging from a 
 negative f(x). 

 In pseudocode, I want to do:

 function counter(outputarray)
 count = 0
 for i in Int64[1:len(outputarray)]
 if outputarray[i] = 0  outputarray[i-1]  0
 count += 1
 end
 end
 return count
 end

 What would be the most efficient way of doing this in Julia?

 Thanks,
 Wally 

Re: [julia-users] Numerical types and machine word size

[Gunnar Farnebäck]

In medical imaging, and probably most other parts of image processing with 
any kind of real-time expectations, single precision is used because the 
factor of 2 for SIMD operations is an enormously big deal, regardless 
whether the computations are done on a GPU or a CPU. There are pitfalls 
related to the limited precision but in general single precision is good 
enough that speed is much more important. In fact it's not unheard of going 
down to 16-bit fixed precision in specific computations to get even more 
out of the SIMD registers.

So yes, there are substantial (or at least non-trivial) parts of the Julia 
target audience that need single precision floats. As others have noted, 
this is completely unrelated to 32 or 64 bit architectures and the width of 
the Int type.

Den onsdagen den 30:e juli 2014 kl. 05:14:41 UTC+2 skrev Tony Kelman:

 There were also architectural problems related to SIMD width (I think) on 
 most GPU designs until fairly recently, causing operations on 64-bit 
 doubles to be significantly more than 2x slower than 32-bit singles. This 
 is not so much the case any more, but old habits die hard and if single 
 precision is good enough, then you still do get a factor of roughly 2 
 benefit in throughput on modern GPU's. There are examples in games (anyone 
 here ever play Kerbal Space Program?) where ingrained habits of single 
 precision have caused noticeable problems.

[julia-users] Re: Trouble with BoundsError()

[Gunnar Farnebäck]

BoundsError means you're trying to index out of bounds. I see nothing in 
the backtrace indicating the splice function being involved and I have no 
way of knowing which line in your function is numbered 204 but there's a 
good chance Varraycut[i-1] will try to index at 0 in the first turn of your 
loop when i=1.

Den onsdagen den 30:e juli 2014 kl. 01:14:53 UTC+2 skrev yaois...@gmail.com:

 Hi Julia users,


 There are a number of things wrong with this script that I am working on, 
 but for the time being, I am trying to resolve a BoundsError() triggered 
 while slicing an array. Code is attached. 

 The source of the BoundsError() appears to be this function. 

 function spikecounter(Varray)
 #spiketimes = Int64[]
 Varraycut = splice!(Varray,1:spliceend)
 spikecount = 0
 for i in Int64[1:len(Varraycut)]
 if Varraycut[i] = 0  Varraycut[i-1]  0
 #= push!(spiketimes,i)
 other option of counting spikes is to append spiketimes to an 
 array, and then count length of array =#
 spikecount += 1
 end
 end
 return spikecount
 end


 Because of other errors in my code, Varray usually contains NaN entries. 
 However, these NaN entries are still appending into my arrays, and I was 
 thinking that if these entries exist, it should still be possible to splice 
 the array. Do these NaN entries preclude Varray from being spliced? 

 This is the exact error:

 ERROR: BoundsError()
  in getindex at array.jl:246
  in spikecounter at /Users/xieh/Dropbox/spiking/j2test.jl:204
  in moldakarimov at /Users/xieh/Dropbox/spiking/j2test.jl:223
  in calcchi at /Users/xieh/Dropbox/spiking/j2test.jl:306
  in mcmc at /Users/xieh/Dropbox/spiking/j2test.jl:319
 while loading /Users/xieh/Dropbox/spiking/j2test.jl, in expression 
 starting on line 365



 Thanks,
 Wally  

Re: [julia-users] Re: Why is string concatenation done with * not +

[Gunnar Farnebäck]

Here's an example of a language that does a bit more with string operations.

% pike
Pike v7.8 release 700 running Hilfe v3.5 (Incremental Pike Frontend)
 ab + c;
(1) Result: abc
 abc * 2;
(2) Result: abcabc
 abcabc - ca;
(3) Result: abbc
 abcabc / b;
(4) Result: ({ /* 3 elements */
a,
ca,
c
})
 ({a, b, c}) * x;
(5) Result: axbxc

Den fredagen den 4:e juli 2014 kl. 04:31:42 UTC+2 skrev Laszlo Hars:

 Of course, 2str, 2*str and str*2 could be equivalent, as they are for 
 other data types. For me two (copies of) str is just more natural than str1 
 times str2. Similarly, add str2 to str1 for concatenation follows the 
 natural language better. Also, you could define str1-str2 as deleting str2 
 from (the beginning/end of) str1, if it was there. (even though str2/str1 
 works, too). You can find in algebra ample examples for both multiplicative 
 and additive notations for structures, but either way the notation of 
 string concatenation is the least of our difficulties when learning Julia.

 On Thursday, July 3, 2014 10:12:47 AM UTC-6, Stefan Karpinski wrote:

 On Thu, Jul 3, 2014 at 12:09 PM, Patrick O'Leary patrick...@gmail.com 
 wrote:

 I do find the generalization to string exponentiation - repetition is 
 pretty beautiful.


 This is my favorite part about the * choice. From str+str it's not at all 
 obvious whether 2str or str*2 is the right choice for string repetition. I 
 also agree that's fascinating how attached people get to string 
 concatenation operators. Don't get Jeff started on string interpolation, 
 which is essentially a fancy string concatenation syntax.
  

[julia-users] Re: Repeated array indexing

[Gunnar Farnebäck]

You can do B[E[:],:] or better write E as a vector (rather than a 1x10 
array) directly:
E = [5, 5, 5, 4, 5, 5, 5, 5, 1, 5]
Then B[E,:] works as expected.

Den tisdagen den 10:e juni 2014 kl. 15:40:50 UTC+2 skrev Bruno Rodrigues:

 Hi, I posted a first thread but can't find it anywhere, so I probably 
 messed up somewhere. So here is my problem: let's say we have following 
 arrays:

 B = [17 24 1 8 15; 23 5 7 14 16;4 6 13 20 22; 10 12 19 21 3;11 18 25 2 9]

 E =[5 5 5 4 5 5 5 5 1 5]

 If we do this: 

 B[E]

 we get: 

 B[E]
 1x10 Array{Int64,2}:
  11  11  11  10  11  11  11  11  17  11

 No problems. But how would I do that for all columns of B? In matlab, I'd 
 do this: B(E,:) and get:

 ans =

 111825 2 9
 111825 2 9
 111825 2 9
 10121921 3
 111825 2 9
 111825 2 9
 111825 2 9
 111825 2 9
 1724 1 815
 111825 2 9

 In Julia I must do a loop:

 S = zeros(10,5)
 x = [1 2 3 4 5]

 for i in x
S[:,i] = B[:,i][E[end,:]]
end


 Is there a way to avoid this loop and obtain something like in Matlab 
 using simple indexing notation?

[julia-users] Re: logical indexing... broadcasting

[Gunnar Farnebäck]

If nothing else works you can solve it with a manual repmat.

A different approach is to forgo the logical indexing entirely.

species_codes = [setosa = [1 0 2],
 versicolor = [2 1 0],
 virginica = [0 2 1]]

N = length(iris_data.columns[5])
iris_output = zeros(N, 3)
for k = 1:N
iris_output[k,:] = species_codes[iris_data.columns[5][k]]
end

Den onsdagen den 14:e maj 2014 kl. 14:05:50 UTC+2 skrev Héctor Fabio 
Satizábal Mejía:

 Hello

 I am trying to conduct some tests on machine learning algorithms with 
 Julia but given that I am new to the language I am kind of stuck in a very 
 simple thing.

 I am working with the Iris dataset, and I want to translate the species 
 name into a code like this:

 setosa = [1 0 2]
 versicolor = [2 1 0]
 virginica = [0 2 1]

 I would like to have something like:

 using RDatasets
 iris_data = dataset(datasets, iris);

 iris_output = zeros(length(iris_data.columns[5]), 3);
 iris_output[iris_data.columns[5].==setosa,:] = [1 0 2];
 iris_output[iris_data.columns[5].==versicolor,:] = [2 1 0];
 iris_output[iris_data.columns[5].==virginica,:] = [0 2 1];

 Unfortunately this is not working. I get this error message:

 ERROR: DimensionMismatch(tried to assign 1x3 array to 50x3 destination)

 Which I think is a broadcasting problem.

 Any help? Is there a compact way of doing that?

 Thanks in advance

 Héctor


 By the way, I manage to do this in R using the following code:

 iris.data - iris
 iris.output - matrix(0, 3, 150)
 iris.output[,iris.data[,5]==setosa] - c(1,0,0)
 iris.output[,iris.data[,5]==versicolor] - c(0,1,0)
 iris.output[,iris.data[,5]==virginica] - c(0,0,1)
 iris.output - t(iris.output)

Re: [julia-users] line continuation

[Gunnar Farnebäck]

You can also use parentheses to make the first line incomplete.

if (ab
 bc)

Den torsdagen den 17:e april 2014 kl. 00:30:41 UTC+2 skrev K leo:

 Thanks for the answer. 

 I tested previously this but it did not work: 
 if ab  bc  cd 
   de  ef 

 But now this works: 
 if ab  bc  cd  
  de  ef 

 So the  got to be on the end of the first for it to know that the line 
 is not finished. 


 On 04/17/2014 06:25 AM, Stefan Karpinski wrote: 
  If the expression is incomplete, continuation is automatic. There is no 
 explicit continuation syntax. 
  
  On Apr 16, 2014, at 6:01 PM, cnbiz850 cnbi...@gmail.com javascript: 
 wrote: 
  
  Searched but didn't find an answer. 
  
  Does Julia allow statement to continue on the next line? 

[julia-users] Re: Performance expectations and benchmarks

[Gunnar Farnebäck]

I have, in the process of evaluating Julia, ported some proprietary matlab 
image processing code, at the order of 3k lines matlab and 1k lines C mex 
code, into pure Julia. The matlab code is fairly well optimized with most 
computations vectorized (usually comes natural in image processing) and 
some critical hot spots offloaded to C. The Julia code is mostly direct 
translations of the vectorized matlab code and the loopy C code. 

For recent versions of Julia this results in about twice the speed of 
Matlab 2007b. With Matlab 2013b the difference is substantially smaller but 
should still be in Julia's favor (not fully investigated yet). However, the 
same algorithms also have highly optimized C++ implementations, which are 
more than 10 times faster than both Julia and Matlab.

My profiling analysis points to two areas that primarily hold back Julia's 
performance on this code:
1. For the tight loop code ported from C, Julia loses on SIMD optimizations 
and OpenMP multithreading.
2. For the vectorized code Julia seems to spend way too much time on memory 
management.

Den torsdagen den 17:e april 2014 kl. 03:24:54 UTC+2 skrev Gilberto Noronha:

 Hi all,

 I was looking at some benchmark results (mainly the ones on julialang.org) 
 and I could not resist the thought that they are a bit misleading. I say 
 this as a Julia fan (I wish I did not think this way!) because 
 realistically I don't think one can expect a real world Julia application 
 to be even close (performance-wise)  to a highly optimized C, C++ or 
 Fortran program (the website suggests that the performance ratio is almost 
 one to one).

 In my (highly subjective) experience it is possible to hope for a Julia 
 program to be competitive with Java and Go (which is already a very good 
 thing, when Python and Matlab can be hundreds of times slower!). 
 But I am curious about other's experiences and whether I could find more 
 realistic benchmarks. In particular, I am curious about benchmarks 
 involving custom objects/classes. My gut feeling is that Julia is not as 
 good with those, but I believe it still beats Matlab very easily (because, 
 among other things, Matlab classes seem to be horribly slow).

 Any thoughts?

 Gilberto  

[julia-users] Re: question on multiple assignment

[Gunnar Farnebäck]

You can do multiple assignment with anything you can iterate over on the 
right hand side, in this case you could do e.g. aa, bb = 12:3:15 with the 
same assignment to aa and bb. One difference is in the value of the 
expression, as printed by the REPL, which coincides with the RHS but would 
usually be discarded. A more interesting difference is how well the code 
generation can optimize away or avoid generating temporary objects for the 
RHS, in your case a tuple and an array respectively, in my example a range. 
I don't have sufficiently deep insights to say something with confidence 
here but in my opinion your first option is the most natural one and I'd be 
surprised if any other solution would be more efficient.

Den torsdagen den 17:e april 2014 kl. 05:40:18 UTC+2 skrev K leo:

 Can anyone explain the difference of the following 2 ways?  Which is 
 more preferred? 

 julia aa, bb = 12, 15 
 (12,15) 

 julia aa, bb = [12, 15] 
 2-element Array{Int64,1}: 
   12 
   15 

[julia-users] Re: preserve column type when building a matrix

[Gunnar Farnebäck]

Ugly workaround:

julia [Array(Any,3,0) [1:3] float([1:3])]
3x2 Array{Any,2}:
 1  1.0
 2  2.0
 3  3.0

Den tisdagen den 1:e april 2014 kl. 10:39:15 UTC+2 skrev harven:

 When I build a matrix from columns using vcat, the type of the columns is 
 sometimes modified, e.g.

 julia [[1:3] float([1:3])]
 3x2 Array{Float64,2}:
  1.0  1.0
  2.0  2.0
  3.0  3.0

 Here the first column of integers has been converted to float. How can I 
 prevent this conversion so as to obtain
 3x2 Array{Any,2}:# or Array{Real,2} since Real is the join type here
  1  1.0
  2  2.0
  3  3.0

 I don't care about the type of the resulting array. The point here is not 
 losing the information concerning the types of the columns. Sincerely, 

[julia-users] Re: Constructing arrays with dim 2

[Gunnar Farnebäck]

Ok, one might be less redundant if reading the thread properly before 
posting. Sorry about the noise.

Den måndagen den 24:e mars 2014 kl. 23:12:03 UTC+1 skrev Gunnar Farnebäck:

 [f(i,j)[k] for k in 1:K, i in 1:n, j in i:m]

 Den måndagen den 24:e mars 2014 kl. 15:07:49 UTC+1 skrev Linus Mellberg:

 Hi!

 I'm trying to construct a 3 dimensional array from a number of 1 
 dimensional arrays. Essentially what i would like to do is

 a = [f(i, j) for i in 1:n, j in 1:m]

 where f(i, j) is a function that returns an array (note, f has to 
 construct the entire array at the same time). The code above creates a 
 2-dimensional array of arrays, but I would like to get a 3-dimensional 
 array with the arrays returned by f in the first dimension with i and j in 
 the second and third dimension, hope you understand

 a[:,:,1] = [f(1,1) f(2,1) ... f(n,1)]
 a[:,:,2] = [f(1,2) f(2,2) ... f(n,2)]
 .
 .
 .
 a[:,:,m] = [f(1,m) f(2,m) ... f(n,m)]

 f(i,j) are column arrays above.

 It can be achieved by first creating the large matrix and then filling it

 a = zeros(Int64, k, n, m)
 for i in 1:n, j in 1:m
   a[:,i,j] = f(i,j)
 end

 Is this the only way? I find it sort of ugly when its usually possible to 
 do nice construction using comprehensions in other cases.

[julia-users] Re: Method redefinition only half works

[Gunnar Farnebäck]

This is a long-standing issue. 
See https://github.com/JuliaLang/julia/issues/265

Den fredagen den 7:e mars 2014 kl. 17:17:49 UTC+1 skrev David Moon:

 In Julia-0.3.0-prerelease-d2d9bd93ed.app on Mac OS X 10.6:

 julia p() = q()
 p (generic function with 1 method)
 julia p2() = q()
 p2 (generic function with 1 method)
 julia q() = 1
 q (generic function with 1 method)
 julia q()
 1
 julia p()
 1
 julia q() = 2
 q (generic function with 1 method)
 julia q()
 2
 julia p2()
 2
 julia p()
 1   # !!!

 Geeze Louse!  p is still calling the old version of q!  No wonder I am 
 having trouble debugging anything.  I thought I was losing my mind when 
 changes to a previously defined method sometimes did not have the expected 
 effect.

 Shouldn't defining a method decache anything that has inlined the same or 
 a less specific method for that function?

 Would I get better results by using an IDE rather than a plain vanilla 
 read-eval-print loop?  Which IDE?

Re: [julia-users] Julia implementation of zlib and gzip decompression

[Gunnar Farnebäck]

That's actually also a gzip decompression implementation, so there's a 
large overlap in that they both implement the inflate algorithm. Otherwise 
the interfaces are rather different and her implementation includes a 
visualization branch whereas my implementation is two orders of magnitude 
faster.

Den lördagen den 4:e januari 2014 kl. 19:03:33 UTC+1 skrev Stefan Karpinski:

 Possibly relevant is Julia Evan's gzip implementation:

 http://jvns.ca/blog/2013/10/24/day-16-gzip-plus-poetry-equals-awesome/


 On Sat, Jan 4, 2014 at 6:54 AM, Gunnar Farnebäck 
 gun...@lysator.liu.sejavascript:
  wrote:

 At some point when I had trouble loading zlib I wrote a pure Julia 
 implementation of zlib and gzip decompression.

 Pros:
 * Written in Julia.
 * No external dependencies.

 Cons:
 * Only supports decompression and only from a buffer (i.e. no streaming).
 * Substantially slower than zlib (about five times when I wrote it, 
 current status unknown).
 * Much less tested than zlib.

 Is this of interest to someone else, e.g. as a package or as a 
 performance test?

 The code can be viewed at

 https://gist.github.com/GunnarFarneback/8254567

[julia-users] Julia implementation of zlib and gzip decompression

[Gunnar Farnebäck]

At some point when I had trouble loading zlib I wrote a pure Julia 
implementation of zlib and gzip decompression.

Pros:
* Written in Julia.
* No external dependencies.

Cons:
* Only supports decompression and only from a buffer (i.e. no streaming).
* Substantially slower than zlib (about five times when I wrote it, current 
status unknown).
* Much less tested than zlib.

Is this of interest to someone else, e.g. as a package or as a performance 
test?

The code can be viewed at

https://gist.github.com/GunnarFarneback/8254567