@Douglas Bates I guess you wanted this again for the new matrix A you
generated with randn()?
>> A
A =
-0.1608 -1.1589 -2.0501
-1.4828 0.2281 -0.5360
-0.3844 -1.2645 0.9202
-1.3541 -0.0305 -2.0846
1.0384 -1.7348 2.8537
-1.1627 1.1079 0.2988
-0.0203 -0.0967 0.3023
-0.7476 3.3424 -0.8231
1.3239 1.7039 1.6594
-1.2648 -1.1727 0.8940
>> orth(A)
ans =
-0.1693 -0.4531 0.3329
-0.1941 -0.1162 -0.4222
0.2334 -0.1589 -0.2756
-0.3703 -0.3270 -0.0951
0.6401 0.0629 -0.1205
-0.1513 0.1612 -0.4768
0.0502 0.0141 -0.0595
-0.5121 0.4801 -0.1718
0.1346 0.5890 0.1307
0.1628 -0.2037 -0.5779
On Thursday, April 3, 2014 1:59:22 AM UTC+2, Andreas Noack Jensen wrote:
>
> I have just tried to build Julia with MKL to check this and with MKL I get
> the same result as in MATLAB. However, both results span the same space.
> With the MKL solution named U and the OpenBLAS solution named U2 I get
>
> norm(U*U'-U2*U2')
>
> and the same is true for the pivoted QR.
>
>
> 2014-04-02 17:53 GMT-04:00 Douglas Bates <[email protected] <javascript:>>:
> >
> > On Wednesday, April 2, 2014 10:46:10 AM UTC-5, Andreas Noack Jensen
> wrote:
> >>
> >> Maybe pivoted QR is best solution here.
> >
> >
> > It might be worthwhile trying out some other examples first. I can't
> reproduce the Matlab answer from qrfact(A), qrfact(A;pivot=true) or
> svdfact(A). That matrix does have a singular value of multiplicity 3 so
> the last 3 columns are only defined up to orthogonal transformation.
> >
> > Here's another example
> > julia> srand(1234321)
> >
> > julia> A = randn((10,3))
> > 10x3 Array{Float64,2}:
> > -0.160767 -1.1589 -2.05007
> > -1.48275 0.228077 -0.536018
> > -0.384419 -1.26447 0.920166
> > -1.3541 -0.0305313 -2.08457
> > 1.03842 -1.73483 2.85367
> > -1.16268 1.10785 0.298823
> > -0.0202548 -0.0966725 0.30227
> > -0.747563 3.3424 -0.823113
> > 1.32394 1.70391 1.65939
> > -1.26479 -1.17272 0.893986
> >
> > julia> A = hcat(A, sum(A,2))
> > 10x4 Array{Float64,2}:
> > -0.160767 -1.1589 -2.05007 -3.36974
> > -1.48275 0.228077 -0.536018 -1.79069
> > -0.384419 -1.26447 0.920166 -0.728726
> > -1.3541 -0.0305313 -2.08457 -3.46921
> > 1.03842 -1.73483 2.85367 2.15726
> > -1.16268 1.10785 0.298823 0.243989
> > -0.0202548 -0.0966725 0.30227 0.185343
> > -0.747563 3.3424 -0.823113 1.77173
> > 1.32394 1.70391 1.65939 4.68723
> > -1.26479 -1.17272 0.893986 -1.54353
> >
> > julia> svdfact(A)
> > SVD{Float64,Float64}(10x4 Array{Float64,2}:
> > -0.445203 -0.0352878 -0.381188 0.452273
> > -0.230914 0.114977 0.403526 -0.592758
> > -0.0696868 -0.285529 0.263285 -0.042519
> > -0.46266 0.190456 0.0526681 -0.0314172
> > 0.320501 -0.552839 0.141094 0.068795
> > 0.0305325 0.194927 0.487114 0.633171
> > 0.0282625 -0.0405129 0.0617822 0.0379557
> > 0.191618 0.665113 0.208002 0.0473809
> > 0.599258 0.136042 -0.0668694 0.0420263
> > -0.171123 -0.248158 0.5577 0.173786
> ,[8.96075,5.27466,2.59348,7.29607e-16],4x4 Array{Float64,2}:
> > 0.248925 0.212351 0.39648 0.857757
> > -0.211592 0.824165 -0.516474 0.0960992
> > -0.802038 0.160183 0.571015 -0.0708406
> > -0.5 -0.5 -0.5 0.5 )
> >
> > julia> qq = qrfact(A; pivot=true);
> >
> > julia> qq[:R]
> > 4x4 Array{Float64,2}:
> > 7.70503 2.17425 2.20126 3.32952
> > 0.0 -4.2384 1.47642 2.76199
> > 0.0 0.0 1.87679 -1.87679
> > 0.0 0.0 0.0 1.31664e-15
> >
> > julia> qq[:Q]
> > 10x10 QRPackedQ{Float64}:
> > -0.437343 0.0490778 0.388686 -0.419191
> > -0.232406 -0.173033 -0.381341 -0.115012
> > -0.094578 0.24982 -0.290425 0.665078
> > -0.450252 -0.22377 -0.0173727 -0.14866
> > 0.279981 0.552938 -0.21007 -0.487088
> > 0.0316663 -0.245138 -0.463803 -0.161661
> > 0.0240548 0.0351485 -0.0666562 -0.040157
> > 0.229944 -0.670641 -0.140445 -0.0845015
> > 0.608334 -0.0899486 0.0626806 -0.158882
> > -0.200328 0.173925 -0.575773 -0.221894
> >
> >
> > I don't use Matlab. Could someone who does check what the result of
> orth(A) is?
> >>
> >>
> >> 2014-04-02 11:39 GMT-04:00 Tomas Lycken <[email protected]>:
> >>
> >>> This is what I get from Matlab R2013a
> >>>
> >>> >> A
> >>>
> >>> A =
> >>>
> >>> 1 1 0 0 0
> >>> 1 1 0 0 0
> >>> 1 1 0 0 0
> >>> 1 1 0 0 0
> >>> 1 1 0 0 0
> >>> 1 0 1 0 0
> >>> 1 0 1 0 0
> >>> 1 0 1 0 0
> >>> 1 0 1 0 0
> >>> 1 0 1 0 0
> >>> 1 0 0 1 0
> >>> 1 0 0 1 0
> >>> 1 0 0 1 0
> >>> 1 0 0 1 0
> >>> 1 0 0 1 0
> >>> 1 0 0 0 1
> >>> 1 0 0 0 1
> >>> 1 0 0 0 1
> >>> 1 0 0 0 1
> >>> 1 0 0 0 1
> >>>
> >>> >> orth(A)
> >>>
> >>> ans =
> >>>
> >>> -0.2236 0.1060 0.3725 -0.0000
> >>> -0.2236 0.1060 0.3725 -0.0000
> >>> -0.2236 0.1060 0.3725 -0.0000
> >>> -0.2236 0.1060 0.3725 -0.0000
> >>> -0.2236 0.1060 0.3725 -0.0000
> >>> -0.2236 -0.3495 -0.0348 -0.1633
> >>> -0.2236 -0.3495 -0.0348 -0.1633
> >>> -0.2236 -0.3495 -0.0348 -0.1633
> >>> -0.2236 -0.3495 -0.0348 -0.1633
> >>> -0.2236 -0.3495 -0.0348 -0.1633
> >>> -0.2236 -0.0143 -0.1302 0.3645
> >>> -0.2236 -0.0143 -0.1302 0.3645
> >>> -0.2236 -0.0143 -0.1302 0.3645
> >>> -0.2236 -0.0143 -0.1302 0.3645
> >>> -0.2236 -0.0143 -0.1302 0.3645
> >>> -0.2236 0.2577 -0.2076 -0.2012
> >>> -0.2236 0.2577 -0.2076 -0.2012
> >>> -0.2236 0.2577 -0.2076 -0.2012
> >>> -0.2236 0.2577 -0.2076 -0.2012
> >>> -0.2236 0.2577 -0.2076 -0.2012
> >>>
> >>> On Wednesday, April 2, 2014 5:22:34 PM UTC+2, Douglas Bates wrote:
> >>>>
> >>>> On Tuesday, April 1, 2014 10:44:03 PM UTC-5, Gustavo Lacerda wrote:
> >>>>>
> >>>>> yes, I think that would be useful.
> >>>>
> >>>>
> >>>> Do you know what Matlab's orth function returns for a rank deficient
> matrix? Suppose you ask for an orthogonal basis of the column space of a
> matrix like
> >>>>
> >>>> julia> hcat(fill(1.,(20,)), eye(4)[iceil([1:20] ./ 5),:])
> >>>> 20x5 Array{Float64,2}:
> >>>> 1.0 1.0 0.0 0.0 0.0
> >>>> 1.0 1.0 0.0 0.0 0.0
> >>>> 1.0 1.0 0.0 0.0 0.0
> >>>> 1.0 1.0 0.0 0.0 0.0
> >>>> 1.0 1.0 0.0 0.0 0.0
> >>>> 1.0 0.0 1.0 0.0 0.0
> >>>> 1.0 0.0 1.0 0.0 0.0
> >>>> 1.0 0.0 1.0 0.0 0.0
> >>>> 1.0 0.0 1.0 0.0 0.0
> >>>> 1.0 0.0 1.0 0.0 0.0
> >>>> 1.0 0.0 0.0 1.0 0.0
> >>>> 1.0 0.0 0.0 1.0 0.0
> >>>> 1.0 0.0 0.0 1.0 0.0
> >>>> 1.0 0.0 0.0 1.0 0.0
> >>>> 1.0 0.0 0.0 1.0 0.0
> >>>> 1.0 0.0 0.0 0.0 1.0
> >>>> 1.0 0.0 0.0 0.0 1.0
> >>>> 1.0 0.0 0.0 0.0 1.0
> >>>> 1.0 0.0 0.0 0.0 1.0
> >>>> 1.0 0.0 0.0 0.0 1.0
> >>>>
> >>>>
> >>>> Do you get 4 or 5 columns?
> >>>>
> >>>> The simple calculation gives you 5 columns but there should only be 4
> columns.
> >>
> >>
> >>
> >>
> >> --
> >> Med venlig hilsen
> >>
> >> Andreas Noack Jensen
>
>
>
>
> --
> Med venlig hilsen
>
> Andreas Noack Jensen
>
