(defun svdcmp (a m n mp np w v &key (nmax 100))
 (declare (type (simple-array double-float (* *)) a)) (declare (type fixnum m))
 (declare (type fixnum n)) (declare (type fixnum mp))
 (declare (type fixnum np)) (declare (type (simple-array double-float (*)) w))
 (declare (type (simple-array double-float (* *)) v))
 (declare (type fixnum nmax))
 (prog
  ((rv1 (make-array '(nmax) :element-type 'double-float)) (x 0.0d0) (z 0.0d0)
   (y 0.0d0) (c 0.0d0) (nm 0) (its 0) (j 0) (h 0.0d0) (f 0.0d0) (k 0) (s 0.0d0)
   (l 0) (i 0) (anorm 0.0d0) (scale 0.0d0) (g 0.0d0)
  )
  (declare (type (simple-array double-float (*)) rv1))
  (declare (type double-float x)) (declare (type double-float z))
  (declare (type double-float y)) (declare (type double-float c))
  (declare (type fixnum nm)) (declare (type fixnum its))
  (declare (type fixnum j)) (declare (type double-float h))
  (declare (type double-float f)) (declare (type fixnum k))
  (declare (type double-float s)) (declare (type fixnum l))
  (declare (type fixnum i)) (declare (type double-float anorm))
  (declare (type double-float scale)) (declare (type double-float g))
  (setf g 0.0) (setf scale 0.0) (setf anorm 0.0)
  (fdo ((i 1 (+ i 1))) ((> i n) nil)
   (tagbody (setf l (+ i 1)) (fset (fref rv1 i) (* scale g)) (setf g 0.0)
    (setf s 0.0) (setf scale 0.0)
    (cond
     ((<= i m)
      (fdo ((k i (+ k 1))) ((> k m) nil)
       (tagbody (setf scale (+ scale (abs (fref a k i)))))
      )
      (cond
       ((/= scale 0.0)
        (fdo ((k i (+ k 1))) ((> k m) nil)
         (tagbody (fset (fref a k i) (/ (fref a k i) scale))
          (setf s (+ s (* (fref a k i) (fref a k i))))
        ))
        (setf f (fref a i i)) (setf g (- (sign (sqrt s) f)))
        (setf h (+ (* f g) (- s))) (fset (fref a i i) (+ f (- g)))
        (cond
         ((/= i n)
          (fdo ((j l (+ j 1))) ((> j n) nil)
           (tagbody (setf s 0.0)
            (fdo ((k i (+ k 1))) ((> k m) nil)
             (tagbody (setf s (+ s (* (fref a k i) (fref a k j)))))
            )
            (setf f (/ s h))
            (fdo ((k i (+ k 1))) ((> k m) nil)
             (tagbody (fset (fref a k j) (+ (fref a k j) (* f (fref a k i)))))
        )))))
        (fdo ((k i (+ k 1))) ((> k m) nil)
         (tagbody (fset (fref a k i) (* scale (fref a k i))))
    )))))
    (fset (fref w i) (* scale g)) (setf g 0.0) (setf s 0.0) (setf scale 0.0)
    (cond
     ((and (<= i m) (/= i n))
      (fdo ((k l (+ k 1))) ((> k n) nil)
       (tagbody (setf scale (+ scale (abs (fref a i k)))))
      )
      (cond
       ((/= scale 0.0)
        (fdo ((k l (+ k 1))) ((> k n) nil)
         (tagbody (fset (fref a i k) (/ (fref a i k) scale))
          (setf s (+ s (* (fref a i k) (fref a i k))))
        ))
        (setf f (fref a i l)) (setf g (- (sign (sqrt s) f)))
        (setf h (+ (* f g) (- s))) (fset (fref a i l) (+ f (- g)))
        (fdo ((k l (+ k 1))) ((> k n) nil)
         (tagbody (fset (fref rv1 k) (/ (fref a i k) h)))
        )
        (cond
         ((/= i m)
          (fdo ((j l (+ j 1))) ((> j m) nil)
           (tagbody (setf s 0.0)
            (fdo ((k l (+ k 1))) ((> k n) nil)
             (tagbody (setf s (+ s (* (fref a j k) (fref a i k)))))
            )
            (fdo ((k l (+ k 1))) ((> k n) nil)
             (tagbody (fset (fref a j k) (+ (fref a j k) (* s (fref rv1 k)))))
        )))))
        (fdo ((k l (+ k 1))) ((> k n) nil)
         (tagbody (fset (fref a i k) (* scale (fref a i k))))
    )))))
    (setf anorm (max anorm (+ (abs (fref w i)) (abs (fref rv1 i)))))
  ))
  (fdo ((i n (+ i (- 1)))) ((> i 1) nil)
   (tagbody
    (cond
     ((< i n)
      (cond
       ((/= g 0.0)
        (fdo ((j l (+ j 1))) ((> j n) nil)
         (tagbody (fset (fref v j i) (/ (/ (fref a i j) (fref a i l)) g)))
        )
        (fdo ((j l (+ j 1))) ((> j n) nil)
         (tagbody (setf s 0.0)
          (fdo ((k l (+ k 1))) ((> k n) nil)
           (tagbody (setf s (+ s (* (fref a i k) (fref v k j)))))
          )
          (fdo ((k l (+ k 1))) ((> k n) nil)
           (tagbody (fset (fref v k j) (+ (fref v k j) (* s (fref v k i)))))
      )))))
      (fdo ((j l (+ j 1))) ((> j n) nil)
       (tagbody (fset (fref v i j) 0.0) (fset (fref v j i) 0.0))
    )))
    (fset (fref v i i) 1.0) (setf g (fref rv1 i)) (setf l i)
  ))
  (fdo ((i n (+ i (- 1)))) ((> i 1) nil)
   (tagbody (setf l (+ i 1)) (setf g (fref w i))
    (cond
     ((< i n)
      (fdo ((j l (+ j 1))) ((> j n) nil) (tagbody (fset (fref a i j) 0.0)))
    ))
    (cond
     ((/= g 0.0) (setf g (/ 1.0 g))
      (cond
       ((/= i n)
        (fdo ((j l (+ j 1))) ((> j n) nil)
         (tagbody (setf s 0.0)
          (fdo ((k l (+ k 1))) ((> k m) nil)
           (tagbody (setf s (+ s (* (fref a k i) (fref a k j)))))
          )
          (setf f (* (/ s (fref a i i)) g))
          (fdo ((k i (+ k 1))) ((> k m) nil)
           (tagbody (fset (fref a k j) (+ (fref a k j) (* f (fref a k i)))))
      )))))
      (fdo ((j i (+ j 1))) ((> j m) nil)
       (tagbody (fset (fref a j i) (* (fref a j i) g)))
     ))
     (t (fdo ((j i (+ j 1))) ((> j m) nil) (tagbody (fset (fref a j i) 0.0))))
    )
    (fset (fref a i i) (+ (fref a i i) 1.0))
  ))
  (fdo ((k n (+ k (- 1)))) ((> k 1) nil)
   (tagbody
    (fdo ((its 1 (+ its 1))) ((> its 30) nil)
     (tagbody
      (fdo ((l k (+ l (- 1)))) ((> l 1) nil)
       (tagbody (setf nm (+ l (- 1)))
        (if (= (+ (abs (fref rv1 l)) anorm) anorm) (go label2))
        (if (= (+ (abs (fref w nm)) anorm) anorm) (go label1))
      ))
      label1 (setf c 0.0) (setf s 1.0)
      (fdo ((i l (+ i 1))) ((> i k) nil)
       (tagbody (setf f (* s (fref rv1 i)))
        (cond
         ((/= (+ (abs f) anorm) anorm) (setf g (fref w i))
          (setf h (sqrt (+ (* f f) (* g g)))) (fset (fref w i) h)
          (setf h (/ 1.0 h)) (setf c (* g h)) (setf s (- (* f h)))
          (fdo ((j 1 (+ j 1))) ((> j m) nil)
           (tagbody (setf y (fref a j nm)) (setf z (fref a j i))
            (fset (fref a j nm) (+ (* y c) (* z s)))
            (fset (fref a j i) (+ (- (* y s)) (* z c)))
      ))))))
      label2 (setf z (fref w k))
      (cond
       ((= l k)
        (cond
         ((< z 0.0) (fset (fref w k) (- z))
          (fdo ((j 1 (+ j 1))) ((> j n) nil)
           (tagbody (fset (fref v j k) (- (fref v j k))))
        )))
        (go label3)
      ))
      (if (= its 30) (error "No convergence in 30 iterations"))
      (setf x (fref w l)) (setf nm (+ k (- 1))) (setf y (fref w nm))
      (setf g (fref rv1 nm)) (setf h (fref rv1 k))
      (setf f
       (/ (+ (* (+ y (- z)) (+ y z)) (* (+ g (- h)) (+ g h))) (* (* 2.0 h) y))
      )
      (setf g (sqrt (+ (* f f) 1.0)))
      (setf f
       (/ (+ (* (+ x (- z)) (+ x z)) (* h (+ (/ y (+ f (sign g f))) (- h)))) x)
      )
      (setf c 1.0) (setf s 1.0)
      (fdo ((j l (+ j 1))) ((> j nm) nil)
       (tagbody (setf i (+ j 1)) (setf g (fref rv1 i)) (setf y (fref w i))
        (setf h (* s g)) (setf g (* c g)) (setf z (sqrt (+ (* f f) (* h h))))
        (fset (fref rv1 j) z) (setf c (/ f z)) (setf s (/ h z))
        (setf f (+ (* x c) (* g s))) (setf g (+ (- (* x s)) (* g c)))
        (setf h (* y s)) (setf y (* y c))
        (fdo ((nm 1 (+ nm 1))) ((> nm n) nil)
         (tagbody (setf x (fref v nm j)) (setf z (fref v nm i))
          (fset (fref v nm j) (+ (* x c) (* z s)))
          (fset (fref v nm i) (+ (- (* x s)) (* z c)))
        ))
        (setf z (sqrt (+ (* f f) (* h h)))) (fset (fref w j) z)
        (cond ((/= z 0.0) (setf z (/ 1.0 z)) (setf c (* f z)) (setf s (* h z)))
        )
        (setf f (+ (* c g) (* s y))) (setf x (+ (- (* s g)) (* c y)))
        (fdo ((nm 1 (+ nm 1))) ((> nm m) nil)
         (tagbody (setf y (fref a nm j)) (setf z (fref a nm i))
          (fset (fref a nm j) (+ (* y c) (* z s)))
          (fset (fref a nm i) (+ (- (* y s)) (* z c)))
      ))))
      (fset (fref rv1 l) 0.0) (fset (fref rv1 k) f) (fset (fref w k) x)
    ))
    label3
  ))
  (return (values a m n mp np w v))
))