Skip to content

Instantly share code, notes, and snippets.

@swannodette
Last active October 3, 2015 00:37
Embed
What would you like to do?
spectral_norm_alioth.js
function A(i,j) {
return 1/((i+j)*(i+j+1)/2+i+1);
}
function Au(u,v) {
for (var i=0; i<u.length; ++i) {
var t = 0;
for (var j=0; j<u.length; ++j)
t += A(i,j) * u[j];
v[i] = t;
}
}
function Atu(u,v) {
for (var i=0; i<u.length; ++i) {
var t = 0;
for (var j=0; j<u.length; ++j)
t += A(j,i) * u[j];
v[i] = t;
}
}
function AtAu(u,v,w) {
Au(u,w);
Atu(w,v);
}
function spectralnorm(n) {
var i, u=[], v=[], w=[], vv=0, vBv=0;
for (i=0; i<n; ++i) {
u[i] = 1; v[i] = w[i] = 0;
}
for (i=0; i<10; ++i) {
AtAu(u,v,w);
AtAu(v,u,w);
}
for (i=0; i<n; ++i) {
vBv += u[i]*v[i];
vv += v[i]*v[i];
}
return Math.sqrt(vBv/vv);
}
console.log(spectralnorm(5500).toFixed(9));
(ns spectral-norm.core)
(defn a [i j]
(/ 1 (+ (/ (* (+ i j) (+ i j 1)) 2) i 1)))
(defn mul-av [n v av]
(dotimes [i n]
(aset av i 0)
(dotimes [j n]
(aset av i (+ (aget av i) (* (a i j) (aget v j)))))))
(defn mul-atv [n v atv]
(dotimes [i n]
(aset atv i 0)
(dotimes [j n]
(aset atv i (+ (aget atv i) (* (a j i) (aget v j)))))))
(defn mul-atav [n v atav]
(let [u (array)]
(mul-av n v u)
(mul-atv n u atav)))
(defn approximate [n]
(let [u (array)
v (array)]
(dotimes [i n] (aset u i 1) (aset v i 0))
(dotimes [i 10]
(mul-atav n u v)
(mul-atav n v u))
(loop [i 0 vbv 0 vv 0]
(if (< i n)
(recur (inc i)
(+ vbv (* (aget u i) (aget v i)))
(+ vv (* (aget v i) (aget v i))))
(Math/sqrt (/ vbv vv))))))
(defn _print [& args]
(.log js/console
(->> args (interpose " ") (apply str))))
(set! *print-fn* _print)
(time
(.log js/console (approximate 5500)))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment