Re: [julia-users] style question about Float initialization

2014-05-01 Thread Stefan Karpinski
If it's in a function, you can just make the function generic, but a bigger
issue is that that push! pattern is bound to be quite slow. It's *much*
better if you can allocate the entire array up front. I'm having a hard
time coming up with examples based on what you wrote – it's a little too
vague. Do you have a more specific example?




On Thu, May 1, 2014 at 12:19 PM, Ethan Anderes ethanande...@gmail.comwrote:

 I have a question which is somewhat related to this question (
 https://groups.google.com/forum/m/?fromgroups#!topic/julia-users/bldp27KJCjY)
 but mine is more of a style question rather than a language design question.

 I find myself often  filling an array of Array{Float64,1}'s in julia.
 I usually do this by initializing an empty array then using push!

 x = Array{Float64,1}[]
 while someconditional
 nextvec = getvec(...)
 push!(x,nextvec)
 end

 I'm wondering how to make the Float64 specification more flexible so that
 the code works when used by someone on a machine using Float32, for
 example, but still keep the performance of having the entries of x all a
 homogeneous type. I would be tempted to initialize
 x=Array{FloatingPoint,1}[] but that allows some entries of x to be
 Array{Float64,1} and others to be Array{Float32,1} for example.






Re: [julia-users] style question about Float initialization

2014-05-01 Thread Ethan Anderes
Ok, here is quickest example that I came up with. 


theta_init = [2.0, 3.0, 4.0]
chain = Array{Float64,1}[init]
for k=1:10_000
   theta_prop = chain[end]
   theta_prop += rand(3)
   likelihoodratio = computelr(theta_prop, theta_curr)
   if rand()  minimum([1,likelihoodratio])
   push!(chain, theta_prop)
   else
   push!(chain, theta_curr)
   end
end

I'm generating a markov chain of vectors. The code proposes a new state vector, 
theta_prop, and pushes it to chain if it satisfies a criterion. The problem is 
that I may run it for 10_000 steps, look at the results and then decide to run 
it for longer... so preinitializing isn't as natural. 

BTW: thanks a ton...I can't believe I'm having Stefan look at my code!




Re: [julia-users] style question about Float initialization

2014-05-01 Thread Cameron McBride
Would linked lists be a more natural storage for these events?

On Thursday, May 1, 2014, Ethan Anderes ethanande...@gmail.com wrote:

 Oop the second line of the code block should read

 chain = Array{Float64,1}[theta_init]



Re: [julia-users] style question about Float initialization

2014-05-01 Thread Ethan Anderes
Stefan:


Ok so maybe that wasn't the best example to illustrate the Float64 issue. I 
guess using the parametric types in functions works so solve this (but not the 
push! issue)

function chaininit{T}(theta_init::Vector{T}, n)
chain = Array{T,1}[theta_init]
for k=1:n
   theta_prop = chain[end]
   theta_prop += rand(3)
   likelihoodratio = computelr(theta_prop, theta_curr)
   if rand()  minimum([1,likelihoodratio])
   push!(chain, theta_prop)
   else
   push!(chain, theta_curr)
   end  
end
chain
end
function chainadd!{T}(chain::Vector{Vector{T}}, n)
for k=1:n
   theta_prop = chain[end]
   theta_prop += rand(3)
   likelihoodratio = computelr(theta_prop, theta_curr)
   if rand()  minimum([1,likelihoodratio])
   push!(chain, theta_prop)
   else
   push!(chain, theta_curr)
   end  
end
end




Re: [julia-users] style question about Float initialization

2014-05-01 Thread Ethan Anderes
Cameron:

I'm not sure what linked lists in Julia are but I will look into it. Thanks!


Re: [julia-users] style question about Float initialization

2014-05-01 Thread Jameson Nash
I don't know of a machine that can run Julia for which Float64 is not fast.
However, you could make a typealias and use that for embedded processors
without an fpu for the eventual day when someone ports Julia.

On Thursday, May 1, 2014, Ethan Anderes ethanande...@gmail.com wrote:

 Oop the second line of the code block should read

 chain = Array{Float64,1}[theta_init]



Re: [julia-users] style question about Float initialization

2014-05-01 Thread Ethan Anderes
Jameson:

My response is going to expose my ignorance here ... but I wasn't worried that 
Float64 wouldn't be fast, but rather that on someones machine typeof(2.0) would 
say Float32 instead of Float64. In which case my code wouldn't run without 
changing all the Float64 type initialization to Float32.  


Re: [julia-users] style question about Float initialization

2014-05-01 Thread Ethan Anderes
Ok, so I'm getting the message that my push! stuff isn't a good  default goto. 
I just felt so natural but I can see that John's code is much more performance 
critical. I'll wean myself off push! :) 


Re: [julia-users] style question about Float initialization

2014-05-01 Thread Stefan Karpinski
2.0 is a Float64 everywhere – it's only integer literals whose type depend on 
the system. Even in old 32-bit machines, the FPUs had 64-bit floats, so there 
was never a place where that was an issue; but you really don't want to be 
using 64-bit integers on a 32-bit machine.


 On May 1, 2014, at 1:01 PM, Ethan Anderes ethanande...@gmail.com wrote:
 
 Jameson:
 
 My response is going to expose my ignorance here ... but I wasn't worried 
 that Float64 wouldn't be fast, but rather that on someones machine 
 typeof(2.0) would say Float32 instead of Float64. In which case my code 
 wouldn't run without changing all the Float64 type initialization to Float32. 
  


Re: [julia-users] style question about Float initialization

2014-05-01 Thread Stefan Karpinski
It's not terrible, and you can help it out with sizehint, but it's just 
possible that it will have to move the data around a bunch.

 On May 1, 2014, at 1:06 PM, Ethan Anderes ethanande...@gmail.com wrote:
 
 Ok, so I'm getting the message that my push! stuff isn't a good  default 
 goto. I just felt so natural but I can see that John's code is much more 
 performance critical. I'll wean myself off push! :) 


Re: [julia-users] style question about Float initialization

2014-05-01 Thread Cameron McBride
Is there an idiomatic Julian way to deal with these cases without a large
pre-allocation or push? I understand why your suggestion here would be
faster than the original. This application aside, I see a lot of push!()
usage in some of the libs (eg saving state in Optim.jl, readdlm in base,
etc).

Linked lists (they only live in DataStructures.jl?) have all the obvious
problems that Stefan points out (scattered in memory), but prevent big
re-allocations if you can take the performance hit on any calculations.

Happy to be schooled here.

Cameron

On Thursday, May 1, 2014, John Myles White johnmyleswh...@gmail.com wrote:

 For this kind of problem, you could allocate a big matrix in advance and
 fill it in column-by-column, rather than push onto an array. Something
 roughly like the following would work:

 theta = Array(Float64, 3, 10_000)
 theta[:, 1] = [2.0, 3.0, 4.0]
 for k= 2:10_000
   theta_cur = theta[:, k - 1]
   theta_prop = theta_cur + rand(3)
   likelihoodratio = computelr(theta_prop, theta_curr)
   if rand()  min(1,likelihoodratio)
   theta[:, k] = theta_prop
   else
   theta[:, k] = theta_curr
   end
 end

 You could make something much faster still with devectorization. The code
 for simulated annealing in Optim.jl might be helpful in this regard,
 although it doesn't retain the full history of the chain.

  -- John

 On May 1, 2014, at 9:46 AM, Ethan Anderes 
 ethanande...@gmail.comjavascript:;
 wrote:

  Ok, here is quickest example that I came up with.
 
 
  theta_init = [2.0, 3.0, 4.0]
  chain = Array{Float64,1}[init]
  for k=1:10_000
theta_prop = chain[end]
theta_prop += rand(3)
likelihoodratio = computelr(theta_prop, theta_curr)
if rand()  minimum([1,likelihoodratio])
push!(chain, theta_prop)
else
push!(chain, theta_curr)
end
  end
 
  I'm generating a markov chain of vectors. The code proposes a new state
 vector, theta_prop, and pushes it to chain if it satisfies a criterion. The
 problem is that I may run it for 10_000 steps, look at the results and then
 decide to run it for longer... so preinitializing isn't as natural.
 
  BTW: thanks a ton...I can't believe I'm having Stefan look at my code!