[racket-dev] case-> and for/sum:

I solved it by not using `for/sum' and writing this ridiculous function 
for the recursive base case and the initial values in `x':

   (: zero-of (case-> (Real -> Real)
                      (Number -> Number)))
   (define (zero-of x) 0)

Fortunately, it should get inlined. I also renamed `U' to `V', because 
it was overshadowing the type name.

(define matrix-solve-upper
   ; solve the equation Ux=b
   ; using back substitution
   (case-lambda
     [(V b)
      (define (default-fail)
        (raise-argument-error
         'matrix-solve-upper
         "The upper triangular matrix is not invertible." 0 V))
      (matrix-solve-upper V b default-fail)]
     [(V b fail)
      (define zero (zero-of (matrix-ref V 0 0)))
      (define m (matrix-num-rows V))
      (define x (make-vector m zero))
      (for ([i (in-range (- m 1) -1 -1)])
        (define bi (matrix-ref b i 0))
        (define Vii (matrix-ref V i i))
        (when (zero? Vii) (fail))
        (let loop ([j  (+ i 1)] [s  zero])
          (cond [(j . < . m)
                 (loop (+ j 1) (+ s (* (matrix-ref V i j)
                                       (vector-ref x j))))]
                [else
                 (vector-set! x i (/ (- bi s) Vii))])))
      (vector->matrix m 1 x)]))

Neil ⊥

On 01/03/2013 10:47 AM, Jens Axel Søgaard wrote:
> Ignore the previous example. Here is the example again, now
> with correct usage of case-lambda. The for/sum problem remains.
>
> /Jens Axel
>
>
> #lang typed/racket
> (require math)
>
> (: matrix-solve-upper
>     (All (A) (case->
>               ((Matrix Real)   (Matrix Real)          -> (Matrix Real))
>               ((Matrix Real)   (Matrix Real)   (-> A) -> (U A (Matrix Real)))
>               ((Matrix Number) (Matrix Number)        -> (Matrix Number))
>               ((Matrix Number) (Matrix Number) (-> A) -> (U A (Matrix
> Number))))))
> (define matrix-solve-upper
>    ; solve the equation Ux=b
>    ; using back substitution
>    (case-lambda
>      [(U b)
>       (define (default-fail)
>         (raise-argument-error
>          'matrix-solve-upper
>          "The upper triangular matrix is not invertible." 0 U))
>       (matrix-solve-upper U b default-fail)]
>      [(U b fail)
>        (define m (matrix-num-rows U))
>        (define x (make-vector m 0))
>        (for ([i (in-range (- m 1) -1 -1)])
>          (define bi (matrix-ref b i 0))
>          (define Uii (matrix-ref U i i))
>          (when (zero? Uii) (fail))
>          (define x.Ui (for/sum ([j (in-range (+ i 1) m)])
>                         (* (matrix-ref U i j) (vector-ref x j))))
>          (vector-set! x i (/ (- bi x.Ui) Uii)))
>        (vector->matrix m 1 x)]))
>
>   (define U (matrix [[4 -1 2  3]
>                      [0 -2 7 -4]
>                      [0  0 6  5]
>                      [0  0 0  3]]))
>   (define b (col-matrix [20 -7 4 6]))
>   b ; expected
>   (define x (solve-upper-triangular U b))
>   (matrix* U x)
> _________________________
>    Racket Developers list:
>    http://lists.racket-lang.org/dev
>