Skip to content

Commit 67c7ef6

Browse files
Update basecol.jl (#123)
Simpler crossprod
1 parent e3813e7 commit 67c7ef6

1 file changed

Lines changed: 9 additions & 18 deletions

File tree

src/utils/basecol.jl

Lines changed: 9 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -35,24 +35,15 @@ end
3535
## Crossprod computes [A B C ...]' [A B C ...] without forming it
3636
##
3737
##############################################################################
38-
39-
# Construct [A B C]'[A B C] without generating [A B C]
40-
function crossprod(c::Combination)
41-
out = Array{Float64}(undef, size(c, 2), size(c, 2))
42-
cviews = [view(c, :, i) for i in 1:size(c, 2)]
43-
for j in 1:size(c, 2)
44-
for i in j:size(c, 2)
45-
out[i, j] = dot(cviews[j], cviews[i])
46-
end
47-
end
48-
# make symmetric
49-
for j in 1:size(c, 2), i in 1:(j-1)
50-
out[i, j] = out[j, i]
51-
end
52-
return out
53-
end
5438
crossprod(A::AbstractMatrix) = A'A
55-
crossprod(A::AbstractMatrix...) = crossprod(Combination(A...))
39+
function crossprod(A::AbstractMatrix, B::AbstractMatrix, C::AbstractMatrix)
40+
    u11, u12, u13 = A'A, A'B, A'C
41+
    u22, u23 = B'B, B'C
42+
    u33 = C'C
43+
    hvcat(3, u11,  u12,  u13,
44+
             u12', u22,  u23,
45+
             u13', u23', u33)
46+
end
5647

5748
##############################################################################
5849
##
@@ -76,4 +67,4 @@ end
7667

7768
function getcols(X::AbstractMatrix, basecolX::AbstractVector)
7869
sum(basecolX) == size(X, 2) ? X : X[:, basecolX]
79-
end
70+
end

0 commit comments

Comments
 (0)