Propagated Sam's changes to spectralnorm to the generic and unsafe versions.
This commit is contained in:
Vincent St-Amour
parent
cdfbbc5476
commit
ae242e2f88
|
|
@ -1,59 +1,47 @@
|
|||
#lang racket/base
|
||||
|
||||
;; The Computer Language Benchmarks Game
|
||||
;; http://shootout.alioth.debian.org/
|
||||
;; Translated from Mike Pall's Lua version.
|
||||
|
||||
;; Translated directly from the C# version, which was:
|
||||
;; contributed by Isaac Gouy
|
||||
|
||||
(require racket/cmdline)
|
||||
|
||||
(define (Approximate n)
|
||||
(let ([u (make-vector n 1.0)]
|
||||
[v (make-vector n 0.0)])
|
||||
;; 20 steps of the power method
|
||||
(for ([i (in-range 10)])
|
||||
(MultiplyAtAv n u v)
|
||||
(MultiplyAtAv n v u))
|
||||
|
||||
;; B=AtA A multiplied by A transposed
|
||||
;; v.Bv /(v.v) eigenvalue of v
|
||||
(let loop ([i 0][vBv 0][vv 0])
|
||||
(if (= i n)
|
||||
(sqrt (/ vBv vv))
|
||||
(let ([vi (vector-ref v i)])
|
||||
(loop (add1 i)
|
||||
(+ vBv (* (vector-ref u i) vi))
|
||||
(+ vv (* vi vi))))))))
|
||||
|
||||
;; return element i,j of infinite matrix A
|
||||
(define (A i j)
|
||||
(/ 1.0 (+ (* (+ i j) (/ (+ i (+ j 1)) 2.0)) (+ i 1))))
|
||||
|
||||
;; multiply vector v by matrix A
|
||||
(define (MultiplyAv n v Av)
|
||||
(for ([i (in-range n)])
|
||||
(vector-set! Av i
|
||||
(for/fold ([r 0])
|
||||
([j (in-range n)])
|
||||
(+ r (* (A i j) (vector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A transposed
|
||||
(define (MultiplyAtv n v Atv)
|
||||
(for ([i (in-range n)])
|
||||
(vector-set! Atv i
|
||||
(for/fold ([r 0])
|
||||
([j (in-range n)])
|
||||
(+ r (* (A j i) (vector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A and then by matrix A transposed
|
||||
(define (MultiplyAtAv n v AtAv)
|
||||
(let ([u (make-vector n 0.0)])
|
||||
(MultiplyAv n v u)
|
||||
(MultiplyAtv n u AtAv)))
|
||||
|
||||
(printf "~a\n"
|
||||
(real->decimal-string
|
||||
(Approximate (command-line #:args (n) (string->number n)))
|
||||
9))
|
||||
(require racket/cmdline racket/trace racket/contract)
|
||||
|
||||
(let* ([A (lambda (i j)
|
||||
(let ([ij (+ i j)])
|
||||
(/ 1.0 (+ (* (* ij (+ ij 1))
|
||||
0.5)
|
||||
(+ i 1)))))]
|
||||
[Av
|
||||
(lambda (x y N)
|
||||
(for ([i (in-range N)])
|
||||
(vector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (= j N) a
|
||||
(L (+ a (* (vector-ref x j) (A i j)))
|
||||
(+ j 1)))))))]
|
||||
[Atv
|
||||
(lambda (x y N)
|
||||
(for ([i (in-range N)])
|
||||
(vector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (= j N) a
|
||||
(L (+ a (* (vector-ref x j) (A j i)))
|
||||
(+ j 1)))))))]
|
||||
[AtAv (lambda (x y t N) (Av x t N) (Atv t y N))]
|
||||
[N (command-line #:args (n) (string->number n))]
|
||||
[u (make-vector N 1.0)]
|
||||
[v (make-vector N)]
|
||||
[t (make-vector N)])
|
||||
(for ([i (in-range 10)])
|
||||
(AtAv u v t N)
|
||||
(AtAv v u t N))
|
||||
(displayln (real->decimal-string
|
||||
(sqrt
|
||||
(let L ([vBv 0.0] [vv 0.0] [i 0])
|
||||
(if (= i N) (/ vBv vv)
|
||||
(let ([ui (vector-ref u i)] [vi (vector-ref v i)])
|
||||
(L (+ vBv (* ui vi))
|
||||
(+ vv (* vi vi))
|
||||
(+ i 1))))))
|
||||
9)))
|
||||
|
|
|
|||
|
|
@ -1,49 +0,0 @@
|
|||
#lang racket/base
|
||||
;; The Computer Language Benchmarks Game
|
||||
;; http://shootout.alioth.debian.org/
|
||||
;; Translated from Mike Pall's Lua version.
|
||||
|
||||
(require racket/cmdline racket/trace racket/contract
|
||||
racket/unsafe/ops racket/flonum)
|
||||
|
||||
(let* ([A (lambda (i j)
|
||||
(let ([ij (unsafe-fx+ i j)])
|
||||
(unsafe-fl/ 1.0 (unsafe-fl+ (unsafe-fl* (unsafe-fl* (unsafe-fx->fl ij)
|
||||
(unsafe-fx->fl (unsafe-fx+ ij 1)))
|
||||
0.5)
|
||||
(unsafe-fx->fl (unsafe-fx+ i 1))))))]
|
||||
[Av
|
||||
(lambda (x y N)
|
||||
(for ([i (in-range N)])
|
||||
(unsafe-flvector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (unsafe-fx= j N) a
|
||||
(L (unsafe-fl+ a (unsafe-fl* (unsafe-flvector-ref x j) (A i j)))
|
||||
(unsafe-fx+ j 1)))))))]
|
||||
[Atv
|
||||
(lambda (x y N)
|
||||
(for ([i (in-range N)])
|
||||
(unsafe-flvector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (unsafe-fx= j N) a
|
||||
(L (unsafe-fl+ a (unsafe-fl* (unsafe-flvector-ref x j) (A j i)))
|
||||
(unsafe-fx+ j 1)))))))]
|
||||
[AtAv (lambda (x y t N) (Av x t N) (Atv t y N))]
|
||||
[N (command-line #:args (n) (string->number n))]
|
||||
[u (make-flvector N 1.0)]
|
||||
[v (make-flvector N)]
|
||||
[t (make-flvector N)])
|
||||
(for ([i (in-range 10)])
|
||||
(AtAv u v t N)
|
||||
(AtAv v u t N))
|
||||
(displayln (real->decimal-string
|
||||
(unsafe-flsqrt
|
||||
(let L ([vBv 0.0] [vv 0.0] [i 0])
|
||||
(if (unsafe-fx= i N) (unsafe-fl/ vBv vv)
|
||||
(let ([ui (unsafe-flvector-ref u i)] [vi (unsafe-flvector-ref v i)])
|
||||
(L (unsafe-fl+ vBv (unsafe-fl* ui vi))
|
||||
(unsafe-fl+ vv (unsafe-fl* vi vi))
|
||||
(unsafe-fx+ i 1))))))
|
||||
9)))
|
||||
|
|
@ -1,62 +1,49 @@
|
|||
#lang racket/base
|
||||
|
||||
;; The Computer Language Benchmarks Game
|
||||
;; http://shootout.alioth.debian.org/
|
||||
;; Translated from Mike Pall's Lua version.
|
||||
|
||||
;; Translated directly from the C# version, which was:
|
||||
;; contributed by Isaac Gouy
|
||||
|
||||
(require racket/cmdline
|
||||
racket/flonum)
|
||||
|
||||
(define (Approximate n)
|
||||
(let ([u (make-flvector n 1.0)]
|
||||
[v (make-flvector n 0.0)])
|
||||
;; 20 steps of the power method
|
||||
(for ([i (in-range 10)])
|
||||
(MultiplyAtAv n u v)
|
||||
(MultiplyAtAv n v u))
|
||||
|
||||
;; B=AtA A multiplied by A transposed
|
||||
;; v.Bv /(v.v) eigenvalue of v
|
||||
(let loop ([i 0][vBv 0.0][vv 0.0])
|
||||
(if (= i n)
|
||||
(flsqrt (fl/ vBv vv))
|
||||
(let ([vi (flvector-ref v i)])
|
||||
(loop (add1 i)
|
||||
(fl+ vBv (fl* (flvector-ref u i) vi))
|
||||
(fl+ vv (fl* vi vi))))))))
|
||||
|
||||
;; return element i,j of infinite matrix A
|
||||
(define (A i j)
|
||||
(fl/ 1.0 (fl+ (fl* (->fl (+ i j))
|
||||
(fl/ (->fl (+ i (+ j 1))) 2.0))
|
||||
(->fl (+ i 1)))))
|
||||
|
||||
;; multiply vector v by matrix A
|
||||
(define (MultiplyAv n v Av)
|
||||
(for ([i (in-range n)])
|
||||
(flvector-set! Av i
|
||||
(for/fold ([r 0.0])
|
||||
([j (in-range n)])
|
||||
(fl+ r (fl* (A i j) (flvector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A transposed
|
||||
(define (MultiplyAtv n v Atv)
|
||||
(for ([i (in-range n)])
|
||||
(flvector-set! Atv i
|
||||
(for/fold ([r 0.0])
|
||||
([j (in-range n)])
|
||||
(fl+ r (fl* (A j i) (flvector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A and then by matrix A transposed
|
||||
(define (MultiplyAtAv n v AtAv)
|
||||
(let ([u (make-flvector n 0.0)])
|
||||
(MultiplyAv n v u)
|
||||
(MultiplyAtv n u AtAv)))
|
||||
|
||||
(printf "~a\n"
|
||||
(real->decimal-string
|
||||
(Approximate (command-line #:args (n) (string->number n)))
|
||||
9))
|
||||
(require racket/cmdline racket/trace racket/contract
|
||||
racket/unsafe/ops racket/flonum)
|
||||
|
||||
(let* ([A (lambda (i j)
|
||||
(let ([ij (unsafe-fx+ i j)])
|
||||
(unsafe-fl/ 1.0 (unsafe-fl+ (unsafe-fl* (unsafe-fl* (unsafe-fx->fl ij)
|
||||
(unsafe-fx->fl (unsafe-fx+ ij 1)))
|
||||
0.5)
|
||||
(unsafe-fx->fl (unsafe-fx+ i 1))))))]
|
||||
[Av
|
||||
(lambda (x y N)
|
||||
(for ([i (in-range N)])
|
||||
(unsafe-flvector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (unsafe-fx= j N) a
|
||||
(L (unsafe-fl+ a (unsafe-fl* (unsafe-flvector-ref x j) (A i j)))
|
||||
(unsafe-fx+ j 1)))))))]
|
||||
[Atv
|
||||
(lambda (x y N)
|
||||
(for ([i (in-range N)])
|
||||
(unsafe-flvector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (unsafe-fx= j N) a
|
||||
(L (unsafe-fl+ a (unsafe-fl* (unsafe-flvector-ref x j) (A j i)))
|
||||
(unsafe-fx+ j 1)))))))]
|
||||
[AtAv (lambda (x y t N) (Av x t N) (Atv t y N))]
|
||||
[N (command-line #:args (n) (string->number n))]
|
||||
[u (make-flvector N 1.0)]
|
||||
[v (make-flvector N)]
|
||||
[t (make-flvector N)])
|
||||
(for ([i (in-range 10)])
|
||||
(AtAv u v t N)
|
||||
(AtAv v u t N))
|
||||
(displayln (real->decimal-string
|
||||
(unsafe-flsqrt
|
||||
(let L ([vBv 0.0] [vv 0.0] [i 0])
|
||||
(if (unsafe-fx= i N) (unsafe-fl/ vBv vv)
|
||||
(let ([ui (unsafe-flvector-ref u i)] [vi (unsafe-flvector-ref v i)])
|
||||
(L (unsafe-fl+ vBv (unsafe-fl* ui vi))
|
||||
(unsafe-fl+ vv (unsafe-fl* vi vi))
|
||||
(unsafe-fx+ i 1))))))
|
||||
9)))
|
||||
|
|
@ -1,62 +1,47 @@
|
|||
;; The Computer Language Benchmarks Game
|
||||
;; http://shootout.alioth.debian.org/
|
||||
;; Translated from Mike Pall's Lua version.
|
||||
|
||||
;; Translated directly from the C# version, which was:
|
||||
;; contributed by Isaac Gouy
|
||||
|
||||
(require racket/cmdline)
|
||||
|
||||
(: Approximate (Natural -> Float))
|
||||
(define (Approximate n)
|
||||
(let ([u (make-vector n 1.0)]
|
||||
[v (make-vector n 0.0)])
|
||||
;; 20 steps of the power method
|
||||
(for: : Void ([i : Natural (in-range 10)])
|
||||
(MultiplyAtAv n u v)
|
||||
(MultiplyAtAv n v u))
|
||||
|
||||
;; B=AtA A multiplied by A transposed
|
||||
;; v.Bv /(v.v) eigenvalue of v
|
||||
(let: loop : Float ([i : Natural 0][vBv : Float 0.0][vv : Float 0.0])
|
||||
(if (= i n)
|
||||
(assert (sqrt (/ vBv vv)) inexact-real?)
|
||||
(let ([vi (vector-ref v i)])
|
||||
(loop (add1 i)
|
||||
(+ vBv (* (vector-ref u i) vi))
|
||||
(+ vv (* vi vi))))))))
|
||||
|
||||
;; return element i,j of infinite matrix A
|
||||
(: A (Natural Natural -> Float))
|
||||
(define (A i j)
|
||||
(exact->inexact (/ 1.0 (+ (* (+ i j) (/ (+ i (+ j 1)) 2.0)) (+ i 1)))))
|
||||
|
||||
;; multiply vector v by matrix A
|
||||
(: MultiplyAv (Natural (Vectorof Float) (Vectorof Float) -> Void))
|
||||
(define (MultiplyAv n v Av)
|
||||
(for: : Void ([i : Natural (in-range n)])
|
||||
(vector-set! Av i
|
||||
(for/fold: : Float ([r : Float 0.0])
|
||||
([j : Natural (in-range n)])
|
||||
(+ r (* (A i j) (vector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A transposed
|
||||
(: MultiplyAtv (Natural (Vectorof Float) (Vectorof Float) -> Void))
|
||||
(define (MultiplyAtv n v Atv)
|
||||
(for: : Void ([i : Natural (in-range n)])
|
||||
(vector-set! Atv i
|
||||
(for/fold: : Float ([r : Float 0.0])
|
||||
([j : Natural (in-range n)])
|
||||
(+ r (* (A j i) (vector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A and then by matrix A transposed
|
||||
(: MultiplyAtAv (Natural (Vectorof Float) (Vectorof Float) -> Void))
|
||||
(define (MultiplyAtAv n v AtAv)
|
||||
(let ([u (make-vector n 0.0)])
|
||||
(MultiplyAv n v u)
|
||||
(MultiplyAtv n u AtAv)))
|
||||
|
||||
(printf "~a\n"
|
||||
(real->decimal-string
|
||||
(Approximate (command-line #:args (n) (assert (string->number (assert n string?)) exact-nonnegative-integer?)))
|
||||
9))
|
||||
(require racket/cmdline racket/trace racket/contract racket/flonum)
|
||||
|
||||
(let* ([A (lambda: ((i : Natural) (j : Natural))
|
||||
(let ([ij (+ i j)])
|
||||
(/ 1.0 (+ (* (* ij (+ ij 1))
|
||||
0.5)
|
||||
(+ i 1)))))]
|
||||
[Av
|
||||
(lambda: ((x : (Vectorof Float)) (y : (Vectorof Float)) (N : Natural))
|
||||
(for: ([i : Natural (in-range N)])
|
||||
(vector-set!
|
||||
y i
|
||||
(let: L : Float ([a : Float 0.0] [j : Natural 0])
|
||||
(if (= j N) a
|
||||
(L (+ a (* (vector-ref x j) (A i j)))
|
||||
(+ j 1)))))))]
|
||||
[Atv
|
||||
(lambda: ((x : (Vectorof Float)) (y : (Vectorof Float)) (N : Natural))
|
||||
(for: ([i : Natural (in-range N)])
|
||||
(vector-set!
|
||||
y i
|
||||
(let: L : Float ([a : Float 0.0] [j : Natural 0])
|
||||
(if (= j N) a
|
||||
(L (+ a (* (vector-ref x j) (A j i)))
|
||||
(+ j 1)))))))]
|
||||
[AtAv (lambda: ((x : (Vectorof Float)) (y : (Vectorof Float)) (t : (Vectorof Float)) (N : Natural))
|
||||
(Av x t N) (Atv t y N))]
|
||||
[N (command-line #:args (#{n : String}) (assert (string->number n) exact-nonnegative-integer?))]
|
||||
[u (make-vector N 1.0)]
|
||||
[v (make-vector N 0.0)]
|
||||
[t (make-vector N 0.0)])
|
||||
(for: ([i : Natural (in-range 10)])
|
||||
(AtAv u v t N)
|
||||
(AtAv v u t N))
|
||||
(displayln (real->decimal-string
|
||||
(flsqrt
|
||||
(let: L : Float ([vBv : Float 0.0] [vv : Float 0.0] [i : Natural 0])
|
||||
(if (= i N) (/ vBv vv)
|
||||
(let ([ui (vector-ref u i)] [vi (vector-ref v i)])
|
||||
(L (+ vBv (* ui vi))
|
||||
(+ vv (* vi vi))
|
||||
(+ i 1))))))
|
||||
9)))
|
||||
|
|
|
|||
|
|
@ -1 +0,0 @@
|
|||
#lang s-exp "wrapper.rkt"
|
||||
|
|
@ -1 +0,0 @@
|
|||
#lang s-exp "wrapper.rkt"
|
||||
|
|
@ -1,71 +0,0 @@
|
|||
;; The Computer Language Benchmarks Game
|
||||
;; http://shootout.alioth.debian.org/
|
||||
|
||||
;; Translated directly from the C# version, which was:
|
||||
;; contributed by Isaac Gouy
|
||||
|
||||
(require racket/cmdline
|
||||
racket/require (for-syntax racket/base)
|
||||
(rename-in
|
||||
(filtered-in
|
||||
(lambda (name) (regexp-replace #rx"unsafe-" name ""))
|
||||
racket/unsafe/ops)
|
||||
[fx->fl ->fl])
|
||||
(only-in racket/flonum make-flvector))
|
||||
|
||||
|
||||
(: Approximate (Natural -> Float))
|
||||
(define (Approximate n)
|
||||
(let ([u (make-flvector n 1.0)]
|
||||
[v (make-flvector n 0.0)])
|
||||
;; 20 steps of the power method
|
||||
(for ([i (in-range 10)])
|
||||
(MultiplyAtAv n u v)
|
||||
(MultiplyAtAv n v u))
|
||||
|
||||
;; B=AtA A multiplied by A transposed
|
||||
;; v.Bv /(v.v) eigenvalue of v
|
||||
(let: loop : Float ([i : Natural 0][vBv : Float 0.0][vv : Float 0.0])
|
||||
(if (= i n)
|
||||
(flsqrt (fl/ vBv vv))
|
||||
(let ([vi (flvector-ref v i)])
|
||||
(loop (add1 i)
|
||||
(fl+ vBv (fl* (flvector-ref u i) vi))
|
||||
(fl+ vv (fl* vi vi))))))))
|
||||
|
||||
;; return element i,j of infinite matrix A
|
||||
(: A (Natural Natural -> Float))
|
||||
(define (A i j)
|
||||
(fl/ 1.0 (fl+ (fl* (->fl (+ i j))
|
||||
(fl/ (->fl (+ i (+ j 1))) 2.0))
|
||||
(->fl (+ i 1)))))
|
||||
|
||||
;; multiply vector v by matrix A
|
||||
(: MultiplyAv (Natural FlVector FlVector -> Void))
|
||||
(define (MultiplyAv n v Av)
|
||||
(for ([i (in-range n)])
|
||||
(flvector-set! Av i
|
||||
(for/fold ([r 0.0])
|
||||
([j (in-range n)])
|
||||
(fl+ r (fl* (A i j) (flvector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A transposed
|
||||
(: MultiplyAtv (Natural FlVector FlVector -> Void))
|
||||
(define (MultiplyAtv n v Atv)
|
||||
(for ([i (in-range n)])
|
||||
(flvector-set! Atv i
|
||||
(for/fold ([r 0.0])
|
||||
([j (in-range n)])
|
||||
(fl+ r (fl* (A j i) (flvector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A and then by matrix A transposed
|
||||
(: MultiplyAtAv (Natural FlVector FlVector -> Void))
|
||||
(define (MultiplyAtAv n v AtAv)
|
||||
(let ([u (make-flvector n 0.0)])
|
||||
(MultiplyAv n v u)
|
||||
(MultiplyAtv n u AtAv)))
|
||||
|
||||
(printf "~a\n"
|
||||
(real->decimal-string
|
||||
(Approximate (command-line #:args (n) (assert (string->number (assert n string?)) exact-nonnegative-integer?)))
|
||||
9))
|
||||
|
|
@ -1,65 +1,49 @@
|
|||
;; The Computer Language Benchmarks Game
|
||||
;; http://shootout.alioth.debian.org/
|
||||
;; Translated from Mike Pall's Lua version.
|
||||
|
||||
;; Translated directly from the C# version, which was:
|
||||
;; contributed by Isaac Gouy
|
||||
|
||||
(require racket/cmdline
|
||||
racket/flonum)
|
||||
|
||||
(: Approximate (Natural -> Float))
|
||||
(define (Approximate n)
|
||||
(let ([u (make-flvector n 1.0)]
|
||||
[v (make-flvector n 0.0)])
|
||||
;; 20 steps of the power method
|
||||
(for: : Void ([i : Natural (in-range 10)])
|
||||
(MultiplyAtAv n u v)
|
||||
(MultiplyAtAv n v u))
|
||||
|
||||
;; B=AtA A multiplied by A transposed
|
||||
;; v.Bv /(v.v) eigenvalue of v
|
||||
(let loop ([i 0][vBv 0.0][vv 0.0])
|
||||
(if (= i n)
|
||||
(flsqrt (fl/ vBv vv))
|
||||
(let ([vi (flvector-ref v i)])
|
||||
(loop (add1 i)
|
||||
(fl+ vBv (fl* (flvector-ref u i) vi))
|
||||
(fl+ vv (fl* vi vi))))))))
|
||||
|
||||
;; return element i,j of infinite matrix A
|
||||
(: A (Natural Natural -> Float))
|
||||
(define (A i j)
|
||||
(fl/ 1.0 (fl+ (fl* (->fl (+ i j))
|
||||
(fl/ (->fl (+ i (+ j 1))) 2.0))
|
||||
(->fl (+ i 1)))))
|
||||
|
||||
;; multiply vector v by matrix A
|
||||
(: MultiplyAv (Natural FlVector FlVector -> Void))
|
||||
(define (MultiplyAv n v Av)
|
||||
(for: : Void ([i : Natural (in-range n)])
|
||||
(flvector-set! Av i
|
||||
(for/fold ([r 0.0])
|
||||
([j (in-range n)])
|
||||
(fl+ r (fl* (A i j) (flvector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A transposed
|
||||
(: MultiplyAtv (Natural FlVector FlVector -> Void))
|
||||
(define (MultiplyAtv n v Atv)
|
||||
(for: : Void ([i : Natural (in-range n)])
|
||||
(flvector-set! Atv i
|
||||
(for/fold ([r 0.0])
|
||||
([j (in-range n)])
|
||||
(fl+ r (fl* (A j i) (flvector-ref v j)))))))
|
||||
|
||||
;; multiply vector v by matrix A and then by matrix A transposed
|
||||
(: MultiplyAtAv (Natural FlVector FlVector -> Void))
|
||||
(define (MultiplyAtAv n v AtAv)
|
||||
(let ([u (make-flvector n 0.0)])
|
||||
(MultiplyAv n v u)
|
||||
(MultiplyAtv n u AtAv)))
|
||||
|
||||
(printf "~a\n"
|
||||
(real->decimal-string
|
||||
(Approximate (command-line #:args (n) (assert (string->number (assert n string?)) exact-nonnegative-integer?)))
|
||||
9))
|
||||
(require racket/cmdline racket/trace racket/contract
|
||||
racket/unsafe/ops racket/flonum)
|
||||
|
||||
(let* ([A (lambda: ((i : Natural) (j : Natural))
|
||||
(let ([ij (unsafe-fx+ i j)])
|
||||
(unsafe-fl/ 1.0 (unsafe-fl+ (unsafe-fl* (unsafe-fl* (unsafe-fx->fl ij)
|
||||
(unsafe-fx->fl (unsafe-fx+ ij 1)))
|
||||
0.5)
|
||||
(unsafe-fx->fl (unsafe-fx+ i 1))))))]
|
||||
[Av
|
||||
(lambda: ((x : FlVector) (y : FlVector) (N : Natural))
|
||||
(for: ([i : Natural (in-range N)])
|
||||
(unsafe-flvector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (unsafe-fx= j N) a
|
||||
(L (unsafe-fl+ a (unsafe-fl* (unsafe-flvector-ref x j) (A i j)))
|
||||
(unsafe-fx+ j 1)))))))]
|
||||
[Atv
|
||||
(lambda: ((x : FlVector) (y : FlVector) (N : Natural))
|
||||
(for: ([i : Natural (in-range N)])
|
||||
(unsafe-flvector-set!
|
||||
y i
|
||||
(let L ([a 0.0] [j 0])
|
||||
(if (unsafe-fx= j N) a
|
||||
(L (unsafe-fl+ a (unsafe-fl* (unsafe-flvector-ref x j) (A j i)))
|
||||
(unsafe-fx+ j 1)))))))]
|
||||
[AtAv (lambda: ((x : FlVector) (y : FlVector) (t : FlVector) (N : Natural))
|
||||
(Av x t N) (Atv t y N))]
|
||||
[N (command-line #:args (#{n : String}) (assert (string->number n) exact-nonnegative-integer?))]
|
||||
[u (make-flvector N 1.0)]
|
||||
[v (make-flvector N)]
|
||||
[t (make-flvector N)])
|
||||
(for: ([i : Natural (in-range 10)])
|
||||
(AtAv u v t N)
|
||||
(AtAv v u t N))
|
||||
(displayln (real->decimal-string
|
||||
(unsafe-flsqrt
|
||||
(let: L : Float ([vBv : Float 0.0] [vv : Float 0.0] [i : Natural 0])
|
||||
(if (unsafe-fx= i N) (unsafe-fl/ vBv vv)
|
||||
(let ([ui (unsafe-flvector-ref u i)] [vi (unsafe-flvector-ref v i)])
|
||||
(L (unsafe-fl+ vBv (unsafe-fl* ui vi))
|
||||
(unsafe-fl+ vv (unsafe-fl* vi vi))
|
||||
(unsafe-fx+ i 1))))))
|
||||
9)))
|
||||
Loading…
Reference in New Issue
Block a user