Re: matrix library?

[Damien Mattei <damien.mattei@gmail.com>]

and the rest of code for multiplying a matrix and a vector:

(define (vector->matrix-column v)
  (matrix (vector-map (lambda (i x)
            (make-vector 1 x))
              v)))

(define (matrix-column->vector Mc)
  {v <+ (matrix-v Mc)}
  (vector-map (lambda (i v2) {v2[0]})
          v))

;; scheme@(guile-user)> (define M1 (create-matrix-by-function (lambda
(l c) (+ l c)) 2 3))
;; scheme@(guile-user)> (matrix-v M1)
;; $1 = #(#(0 1 2) #(1 2 3))
;; scheme@(guile-user)> {M1 * #(1 2 3)}
;; $2 = #(8 14)
(define-method (* (M <matrix>) (v <vector>)) ;; args: matrix ,vector ;
 return vector
  {Mc <+ (vector->matrix-column v)}
  ;(display Mc)
  (matrix-column->vector {M * Mc}))


;; define getter,setter
(define (matrix-ref M lin col)
  {v <+ (matrix-v M)}
  {v[lin][col]})

(define (matrix-set! M lin col x)
  {v <+ (matrix-v M)}
  {v[lin][col] <- x})


;; >  (overload-square-brackets matrix-ref matrix-set!  (matrix?
number? number?))
;; >  (overload-square-brackets matrix-line-ref     (lambda (x) '())
(matrix? number?))
;; > $ovrld-square-brackets-lst$
;; '(((#<procedure:matrix?> #<procedure:number?>)
(#<procedure:matrix-line-ref> . #<procedure>))
;;   ((#<procedure:matrix?> #<procedure:number?> #<procedure:number?>)
(#<procedure:matrix-ref> . #<procedure:matrix-set!>)))
;; > (define Mv (matrix #(#(1 2 3) #(4 5 6))))
;; > Mv
;; #<matrix>
;; > (find-getter-for-overloaded-square-brackets (list Mv 1))
;; #<procedure:matrix-line-ref>
;; > {Mv[1 0]}
;; 4
;; > {Mv[1]}
;; '#(4 5 6)
;; > {Mv[1][0]}
;; 4
;; >
(define (matrix-line-ref M lin)
  {v <+ (matrix-v M)}
  {v[lin]})


(define (matrix-line-set! M lin vect-line)
  {v <+ (matrix-v M)}
  {v[lin] <- vect-line})

this allow also to indexes matrix as it is possible in Python and C++
by overloading the [ ] operator for the new type <matrix> :

;; overload [ ]
(overload-square-brackets matrix-ref matrix-set!  (matrix? number? number?))
(overload-square-brackets matrix-line-ref matrix-line-set! (matrix? number?))

scheme@(guile-user)> (define M (create-matrix-by-function (lambda (l
c) (+ l c)) 2 3))
scheme@(guile-user)>  (matrix-v M)
$1 = #(#(0 1 2) #(1 2 3))
scheme@(guile-user)> {M[1 2]}
$2 = 3
scheme@(guile-user)> {M[1][2]}
$3 = 3
scheme@(guile-user)> {M[1]}
$4 = #(1 2 3)

On Thu, Oct 5, 2023 at 4:47 PM Damien Mattei <damien.mattei@gmail.com> wrote:
>
> that is:
>
> (define-method (* (M1 <matrix>) (M2 <matrix>))
>
>   {(n1 p1) <+ (dim-matrix M1)}
>   {(n2 p2) <+ (dim-matrix M2)}
>
>   (when {p1 ≠ n2} (error "matrix.* : matrix product impossible,
> incompatible dimensions"))
>
>   {v1 <+ (matrix-v M1)}
>   {v2 <+ (matrix-v M2)}
>
>   (define (res i j)
>     {sum <+ 0}
>     (for ({k <+ 0} {k < p1} {k <- k + 1})
>      {sum <- sum + v1[i][k] * v2[k][j]}))
>
>   {v <+ (create-vector-2d res n1 p2)}
>   ;(display "v=") (display v) (newline)
>
>   (matrix v))
>
> and the full code:
> a matrix multiplication (could add + - but no needs in my code...)
> with GOOPS and Scheme+ , without using the overload features of
> scheme+ but the native ones of guile scheme:
>
> (define-module (matrix)
>
>   #:use-module (oop goops)
>
>   #:use-module ((guile)
>         #:select (*)
>         #:prefix orig:)
>
>   ;;#:use-module ((srfi srfi-43)) ;; vector
>
>   #:export (<matrix>
>         matrix?
>         matrix
>         matrix-v
>         create-matrix-by-function
>         dim-matrix)
>
>   )
>
>
> (include "/Users/mattei/Dropbox/git/Scheme-PLUS-for-Guile/Scheme+.scm")
>
> (include "/Users/mattei/Dropbox/git/Scheme-PLUS-for-Guile/scheme-infix.scm")
>
> (include "/Users/mattei/Dropbox/git/Scheme-PLUS-for-Guile/assignment.scm")
> (include "/Users/mattei/Dropbox/git/Scheme-PLUS-for-Guile/apply-square-brackets.scm")
>
> (include "/Users/mattei/Dropbox/git/Scheme-PLUS-for-Guile/array.scm")
> ;; must know <- from assignment
>
>
> (define-class <matrix> ()
>
>   (v #:accessor matrix-v
>      #:getter matrix-get-v
>      #:setter matrix-set-v!
>      #:init-keyword #:v)
>
>   )
>
>
> (define-method (matrix vinit)
>   (make <matrix> #:v vinit))
>
> (define-method (matrix? obj)
>   (equal? (class-name (class-of obj)) (class-name <matrix>)))
>
> ;; (define M (create-matrix-by-function (lambda (l c) (+ l c)) 2 3))
> ;; $1 = #<<matrix> 1078a0a70>
> (define (create-matrix-by-function fct lin col)
>   (matrix (create-vector-2d fct lin col)))
>
> ;; return the line and column values of dimension
> ;; scheme@(guile-user)> (dim-matrix M)
> ;; $1 = 2
> ;; $2 = 3
> (define (dim-matrix M)
>
>   (when (not (matrix? M))
>     (error "argument is not of type matrix"))
>
>   {v <+ (matrix-v M)}
>   {lin <+ (vector-length v)}
>   {col <+ (vector-length {v[0]})}
>   (values lin col))
>
>
> ;; scheme@(guile-user)> (define M1 (create-matrix-by-function (lambda
> (l c) (+ l c)) 2 3))
> ;; scheme@(guile-user)> (define M2 (create-matrix-by-function (lambda
> (l c) (- l c)) 3 2))
> ;; scheme@(guile-user)> (matrix-v M1)
> ;; $2 = #(#(0 1 2) #(1 2 3))
> ;; scheme@(guile-user)> (matrix-v M2)
> ;; $3 = #(#(0 -1) #(1 0) #(2 1))
> ;; (define M1*M2 {M1 * M2})
> ;; scheme@(guile-user)> (matrix-v M1*M2)
> ;; $2 = #(#(5 2) #(8 2))
>
> (define-method (* (M1 <matrix>) (M2 <matrix>))
>
>   {(n1 p1) <+ (dim-matrix M1)}
>   {(n2 p2) <+ (dim-matrix M2)}
>
>   (when {p1 ≠ n2} (error "matrix.* : matrix product impossible,
> incompatible dimensions"))
>
>   {v1 <+ (matrix-v M1)}
>   {v2 <+ (matrix-v M2)}
>
>   (define (res i j)
>     {sum <+ 0}
>     (for ({k <+ 0} {k < p1} {k <- k + 1})
>      {sum <- sum + v1[i][k] * v2[k][j]}))
>
>   {v <+ (create-vector-2d res n1 p2)}
>   ;(display "v=") (display v) (newline)
>
>   (matrix v))
>
>
> and the test:
>
>  (define M1 (create-matrix-by-function (lambda (l c) (+ l c)) 2 3))
> ;(define M2 (create-matrix-by-function (lambda (l c) (- l c)) 3 2))
> (matrix-v M1)
>  #(#(0 1 2) #(1 2 3))
> (matrix-v M2)
>  #(#(0 -1) #(1 0) #(2 1))
> (define M1*M2 {M1 * M2})
>  (matrix-v M1*M2)
>  #(#(5 2) #(8 2))
>
> additional function used:
>
> ;; create a vector of line and column with a function
> (define (create-vector-2d fct lin col)
>   {v <+ (make-vector lin)}
>   ;;(display "ok") (newline)
>   (for ({l <+ 0} {l < lin} {l <- l + 1})
>        {v[l] <- (make-vector col)}
>        (for ({c <+ 0} {c < col} {c <- c + 1})
>         {v[l][c] <- (fct l c)}))
>   v)
>
> and a corrected version of :
>
> (define-syntax for
>    (lambda (stx)
>      (syntax-case stx ()
>        ((kwd (init test incrmt) body ...)
>         (with-syntax ((BREAK (datum->syntax (syntax kwd) 'break))
>                       (CONTINUE (datum->syntax (syntax kwd) 'continue)))
>              (syntax
>               (call/cc
>                (lambda (escape)
>              (let-syntax ((BREAK (identifier-syntax (escape))))
>                init
>                (let loop ((res 0)) ;; now we will return a result at
> the end if no break but if we continue? what happens?
>                  (if test
>                  (begin
>                    (call/cc
>                     (lambda (next)
>                       (set! res (let-syntax ((CONTINUE
> (identifier-syntax (next))))
>                           (let () ;; allow definitions
>                             body ...)))))
>                    incrmt
>                    (loop res))
>                  res)
>                  ))))
>               ) ;; close syntax
>
>              )))))
>
> note that expression
>
> {sum <- sum + v1[i][k] * v2[k][j]}
>
> is expanded as:
>
> ($nfx$ sum <- sum + ($bracket-apply$ ($bracket-apply$ v1 i) k) *
> ($bracket-apply$ ($bracket-apply$ v2 k) j))
>
> which is not fast to evaluate,and slow down calculus on  intensive
> computation, as the infix precedence procedure is evaluated each time
> ,even if it never change in the computation loop or never in program,
> i'm working on an optimiser (written in scheme) of code before
> evaluation , that modify the code (transform the infix expression with
> precedence in prefix one definitively ) before it is passed to the
> guile compiler.But it is another piece of code...
>
> On Wed, Oct 4, 2023 at 10:29 PM Damien Mattei <damien.mattei@gmail.com> wrote:
> >
> > yes i will re-implement that from scratch because i can not find
> > something easy in libraries to use.
> > I already made it for Racket/Scheme+:
> > (struct matrix-vect (v)) ;; matrix based on vector of vectors
> >
> > (define (create-matrix-vect-by-function fct lin col)
> >   (matrix-vect (create-vector-2d fct lin col)))
> >
> >
> > ;; return the line and column values of dimension
> > (define (dim-matrix-vect M)
> >
> >   (when (not (matrix-vect? M))
> >     (error "argument is not of type matrix"))
> >
> >   {v <+ (matrix-vect-v M)}
> >   {lin <+ (vector-length v)}
> >   {col <+ (vector-length {v[0]})}
> >   (values lin col))
> >
> >
> > (define (multiply-matrix-matrix M1 M2)
> >
> >   {(n1 p1) <+ (dim-matrix-vect M1)}
> >   {(n2 p2) <+ (dim-matrix-vect M2)}
> >
> >   (when {p1 ≠ n2} (error "matrix-by-vectors.* : matrix product
> > impossible, incompatible dimensions"))
> >
> >   {v1 <+ (matrix-vect-v M1)}
> >   {v2 <+ (matrix-vect-v M2)}
> >
> >   (define (res i j)
> >     (for/sum ([k (in-range p1)])
> >          {v1[i][k] * v2[k][j]}))
> >
> >   {v <+ (create-vector-2d res n1 p2)}
> >   ;(display "v=") (display v) (newline)
> >
> >   (matrix-vect v))
> >
> >
> > (define (vector->matrix-column v)
> >   (matrix-vect (vector-map (lambda (x) (make-vector 1 x))
> >                v)))
> >
> > (define (matrix-column->vector Mc)
> >   {v <+ (matrix-vect-v Mc)}
> >   (vector-map (lambda (v2) {v2[0]})
> >           v))
> >
> > (define (multiply-matrix-vector M v) ;; args: matrix ,vector ;  return vector
> >   {Mc <+ (vector->matrix-column v)}
> >   ;(display Mc)
> >   (matrix-column->vector (multiply-matrix-matrix M Mc)))
> >
> >
> >  there should be no difference for computation in Guile ,just that i
> > will use GOOPS instead of struct and that it should be a guile module
> > written in scheme+ syntax....
> > see it when finished....
> >
> >
> > On Wed, Oct 4, 2023 at 8:28 PM <tomas@tuxteam.de> wrote:
> > >
> > > On Wed, Oct 04, 2023 at 05:42:10PM +0200, Damien Mattei wrote:
> > > > thank you. i will try to find a simple example in.
> > > >
> > > > On Wed, Oct 4, 2023 at 1:36 PM Nala Ginrut <nalaginrut@gmail.com> wrote:
> > > > >
> > > > > And I'd mention AIScm which bound tensorflow APIs. It needs more love.
> > > > >
> > > > > https://wedesoft.github.io/aiscm/
> > > > >
> > > > >
> > > > > On Wed, Oct 4, 2023, 17:12 <tomas@tuxteam.de> wrote:
> > > > >>
> > > > >> On Wed, Oct 04, 2023 at 08:46:02AM +0200, Damien Mattei wrote:
> > > > >> > hello,
> > > > >> > does anyone know if it exists a basic matrix library for Guile?
> > > > >> > just need to multiply a matrix and a vector.
> > > > >>
> > > > >> Perhaps this thread from guile-user [1] has something for you
> > > > >>
> > > > >> Cheers
> > > > >>
> > > > >> [1] https://lists.gnu.org/archive/html/guile-user/2018-12/threads.html#00117
> > >
> > > That all said, if your case doesn't get much hairier,
> > > it is fairly easy to navigate with Guile's arrays:
> > >
> > > (define (dot v1 v2)
> > >   "The dot product of v1 and v2"
> > >   (let ((sum 0))
> > >     (array-for-each
> > >       (lambda (x1 x2)
> > >         (set! sum (+ sum (* x1 x2))))
> > >       v1 v2)
> > >     sum))
> > >
> > > (define (row m i)
> > >   "The i-th row of m, as a vector"
> > >   (make-shared-array m (lambda (j) (list i j)) (car (array-dimensions m))))
> > >
> > > (define (mat* m v)
> > >   "Left multiply matrix m and vector v"
> > >   (let* ((height (cadr (array-dimensions m)))
> > >          (res (make-typed-array 'f64 0 height)))
> > >     (do ((i 0 (1+ i)))
> > >         ((>= i height))
> > >         (array-set! res (dot (row m i) v) i))
> > >     res))
> > >
> > >
> > > (define vec #1f64(1 2 1))
> > > (define mat #2f64((0.5 0.5 0) (0 0.5 0.5) (-0.5 0 0.5)))
> > >
> > > (mat* mat vec)
> > >
> > >  => (1.5 1.5 0)
> > >
> > > CAVEAT: only tested with this one example. Input value checking
> > > left as an...
> > >
> > > Cheers
> > > --
> > > t
> > >