General update

- Generation of matrix variates now follows Distributions.jl's abstract style
- Matrix distributions are declared as abstract types
- rand() generates dense matrix variates
- bidrand(), tridrand(), sprand() generate bidiagonal, tridiagonal and general sparse variates
- eigvalrand(), eigvaljpdf() generate random eigenvalue sets and calculate their joint pdf
- Remove old fudges to generate full matrices to work around lack of bidiagonal support in pre-0.2 Julia
- Remove obsolete custom extensions to Distributions.jl
- Retab to 2 space indent
- Make GSL dependency optional

Author: jiahao

src/FastHistogram.jl

@@ -11,22 +11,22 @@ export hist_eig
 #
 #For general matrices, this is slower than hist(eigvals(M), bins) and is NOT recommended
 function hist_eig{GridPoint <: Number}(M::AbstractMatrix, bins::Vector{GridPoint})
-    n = size(M)[1]
-    NumBins = length(bins)
-    histogram = zeros(NumBins)
-    SturmSequence = zeros(n)
-    for BinId=1:NumBins
-        K = M - bins[BinId]*eye(n)
-        #Progression of determinants of lower-right submatrices
-        #SturmSequence = [det(K[end-j+1:end,end-j+1:end]) for j=1:n]
-        #SturmRatioSequence = [SturmSequence[j+1]/SturmSequence[j] for j=1:n]
+  n = size(M)[1]
+  NumBins = length(bins)
+  histogram = zeros(NumBins)
+  SturmSequence = zeros(n)
+  for BinId=1:NumBins
+    K = M - bins[BinId]*eye(n)
+    #Progression of determinants of lower-right submatrices
+    #SturmSequence = [det(K[end-j+1:end,end-j+1:end]) for j=1:n]
+    #SturmRatioSequence = [SturmSequence[j+1]/SturmSequence[j] for j=1:n]
 
-        SturmSequenceSigns = [int(sign(det(K[end-j+1:end,end-j+1:end]))) for j=1:n]
-        SturmSequenceSigns = [1; SturmSequenceSigns]
-        SturmRatioSequence = [SturmSequenceSigns[j+1]/SturmSequenceSigns[j] for j=1:n]
-        histogram[BinId] = sum([r < 0 for r in SturmRatioSequence])
-    end
-    histogram = int(diff([histogram; n]))
+    SturmSequenceSigns = [int(sign(det(K[end-j+1:end,end-j+1:end]))) for j=1:n]
+    SturmSequenceSigns = [1; SturmSequenceSigns]
+    SturmRatioSequence = [SturmSequenceSigns[j+1]/SturmSequenceSigns[j] for j=1:n]
+    histogram[BinId] = sum([r < 0 for r in SturmRatioSequence])
+  end
+  histogram = int(diff([histogram; n]))
 end
 
 
@@ -39,19 +39,19 @@ end
 #   the Beta-Hermite Random Matrix Ensemble", M.Sc. thesis (MIT), 2007
 #
 function hist_eig{GridPoint <: Number}(M::SymTridiagonal, bins::Vector{GridPoint})
-    n = size(M)[1]
-    NumBins = length(bins)
-    histogram = zeros(NumBins)
-    r = zeros(n)
-    for BinId=1:NumBins
-        a, b = M.dv - bins[BinId], M.ev
-        #Formula (2.3)
-        r[1] = a[1]
-        for i=2:n
-            r[i] = a[i] - b[i-1]^2/r[i-1]
-        end
-        histogram[BinId] = sum([SturmRatio < 0 for SturmRatio in r])
+  n = size(M)[1]
+  NumBins = length(bins)
+  histogram = zeros(NumBins)
+  r = zeros(n)
+  for BinId=1:NumBins
+    a, b = M.dv - bins[BinId], M.ev
+    #Formula (2.3)
+    r[1] = a[1]
+    for i=2:n
+      r[i] = a[i] - b[i-1]^2/r[i-1]
     end
-    histogram = int(diff([histogram; n]))
+    histogram[BinId] = sum([SturmRatio < 0 for SturmRatio in r])
+  end
+  histogram = int(diff([histogram; n]))
 end