Excellent point, and that's indeed the solution. Thanks!
It works with R = float64(A) instead of R = copy(A).
==========
*julia> **A = [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]*
*5x5 Array{Int64,2}:*
* 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*
*julia> **B=float(A)*
*5x5 Array{Float64,2}:*
* 17.0 24.0 1.0 8.0 15.0*
* 23.0 5.0 7.0 14.0 16.0*
* 4.0 6.0 13.0 20.0 22.0*
* 10.0 12.0 19.0 21.0 3.0*
* 11.0 18.0 25.0 2.0 9.0*
*julia> **bug(A)*
*ERROR: InexactError()*
* in setindex! at multidimensional.jl:76*
* in bug at [...].jl:75*
*julia> **bug(B)*
*(*
*5x5 Array{Float64,2}:*
* -0.523387 0.505754 0.67347 -0.121542 0.0441038*
* -0.708111 -0.696573 -0.0177274 0.0815416 0.0800023*
* -0.12315 0.136742 -0.355751 -0.630744 0.664635 *
* -0.307875 0.191057 -0.412231 -0.424672 -0.720021 *
* -0.338662 0.451439 -0.499631 0.632767 0.177441 ,*
*5x5 Array{Float64,2}:*
* -32.4808 -26.6311 -21.3973 -23.7063 -25.8615 *
* 0.0 19.8943 12.3234 1.94393 4.08558*
* 0.0 0.0 -24.3985 -11.6316 -3.74148*
* 0.0 0.0 0.0 -20.0982 -9.97394*
* 0.0 0.0 0.0 0.0 16.0005 )*
*julia> *
On Saturday, September 20, 2014 5:30:23 PM UTC-7, Andreas Noack wrote:
>
> It is easier to help you if you provide an example that can actually run.
> What is A when you get the error.
>
> A guess could be that it is something like a Matrix{Int} and therefore R =
> copy(A) is also Matrix{Int}. Hence, the line R[k:m, :] = P * R[k:m, :]
> fails as the right hand side consists of floating point values that are not
> integer valued. It works when you use eye() because it is integer valued
> even though it is a Matrix{Float64}.
>
> Med venlig hilsen
>
> Andreas Noack
>
> 2014-09-20 20:24 GMT-04:00 Clemens Heitzinger <[email protected]
> <javascript:>>:
>
>> [I am reposting this from julia-dev, where it was suggested to try here.]
>>
>> Dear all:
>>
>> I am seeing a curious problem when assigning submatrices, e.g., in
>> R[k:m, :] = P * R[k:m, :]
>>
>> I was implementing QR factorization and here is the code.
>>
>> ==========
>> function bug(A)
>> (m, n) = size(A)
>> Q = eye(m)
>> R = copy(A)
>>
>> for k in 1:n
>> z = R[k:m, k]
>> v = [-sign(z[1]) * norm(z) - z[1]; -z[2:end]]
>> P = eye(m-k+1) - 2 * (v * v') / dot(v, v) # does not work
>> ##P = eye(m-k+1) - (v * v') # does not work
>> ##P = eye(m-k+1) # works
>> ##P = rand(m-k+1, m-k+1) # does not work
>> R[k:m, :] = P * R[k:m, :]
>> Q[k:m, :] = P * Q[k:m, :]
>> end
>>
>> (Q', triu(R))
>> end
>> ==========
>>
>> This is the error I get:
>>
>> ==========
>> ERROR: InexactError()
>> in setindex! at multidimensional.jl:76
>> in bug at [...].jl:75
>> in include at ./boot.jl:245
>> in include_from_node1 at ./loading.jl:128
>> in eval at [...]/src/julia-0.3/base/sysimg.jl:7
>> in anonymous at multi.jl:1310
>> in run_work_thunk at multi.jl:621
>> in run_work_thunk at multi.jl:630
>> in anonymous at task.jl:6
>> while loading [...].jl, in expression starting on line 103
>> ==========
>>
>> The offending line seems to be:
>> R[k:m, :] = P * R[k:m, :]
>>
>> Curiously enough, the InexactError does not occur when P is the identity
>> matrix. Here are some data points:
>> P = eye(m-k+1) - 2 * (v * v') / dot(v, v) # does not work
>> P = eye(m-k+1) - (v * v') # does not work
>> P = eye(m-k+1) # works
>> P = rand(m-k+1, m-k+1) # does not work
>>
>> I'd appreciate any insights... It may be a problem with setindex! in
>> multidimensional.jl.
>>
>> Cheers,
>> Clemens
>>
>
>
