
From:  Vincent Belaïche 
Subject:  Submission of Matrix Kronecker Product for calc 
Date:  Sun, 2 Mar 2008 22:37:40 +0100 
Hello, I am new to this list (I have just subscribed). I have just written a function for CALC that computes the Kronecker product of two matrices. Maybe this can be useful to somebody else. I am not sure that this is the correct forum to submit that. Sorry if I bothered anybody. BR, Vincent. PS1 : I think that this is perfectible to the following extent a) I use reverse instead of nreverse, which is more memory consuming. Not knowing the internals of Emacs I was unsure how the garbage collector free the memory if you nreverse a list not passing to nreverse the pointer of the first element (ie the point to the list itself). This is why I used reverse. b) I started from a template created by the "ZP" key for creating a permanent user function. The starting point was a simple multiplication. Then I modified the function to make a Kronecker product. However I don't know how to set the 'calcuserdefn property (what is there is still what there was originally after "ZP" keystroke in calc mode. PS2 : The function works on two dimensional matrices. I had no need to make it generic for ndimensional (and maybe I was not proficient enough in Lisp to do it... ;/ ). PS3 : Code below the  line :  ;;; Definition of Kronecker matrix product (put 'calcdefine 'calckron '(progn (defun calckron nil (interactive) (calcwrapper (calcenterresult 2 "kron" (cons (quote calcFunckron) (calctoplistn 2))))) (put 'calckron 'calcuserdefn 't) (defun calcFunckron (x y) "Kronecker product of matrices x and y. If x and y are not matrices then they are first embedded into matrices before computation. After computation the result may be desembedded from matrix so that the Kronecker product of two scalars is a scalar, of one scalar and a vector or vice verse is a vector, or of two vectors is a vector. " ;; (matrixbacketlevel x) returns 0 if x is not a matrix, ;; 1 if x is a vector that is not a matrix ;; 2 if x is a matrix (defun matrixbacketlevel (x) (let ((returnedvalue 0) (xsub x)) (while (and (Mathvectorp xsub) (< returnedvalue 2)) (setq xsub (elt xsub 1)) (setq returnedvalue (1+ returnedvalue)) ) returnedvalue )) ;; (matrixembed x n) return x embeded in n level of brackets ;; (matrixembed x 0) returns x ;; (matrixembed x 1) returns (vec x) ;; (matrixembed x 2) returns (vec (vec x)) (defun matrixembed (x n) (let ((returnedvalue x) (nctr 0)) (while (< nctr n) (setq returnedvalue (list 'vec returnedvalue)) (setq nctr (1+ nctr)) ) returnedvalue )) ;; evaluate level of bracketting of arguments x and y (let ( (blx (matrixbacketlevel x)) (bly (matrixbacketlevel y)) prodval ) (let ( (embx (matrixembed x ( 2 blx))) (emby (matrixembed y ( 2 bly))) prodnextrow ) (dolist (xrow (reverse (cdr embx)) prodval) (dolist (yrow (reverse (cdr emby))) (setq prodnextrow nil) (dolist (xcol (reverse (cdr xrow)) prodnextrow) (dolist (ycol (reverse (cdr yrow))) (setq prodnextrow (cons (mathmul xcol ycol) prodnextrow)) ) ) (setq prodnextrow (cons 'vec prodnextrow)) (setq prodval (cons prodnextrow prodval)) ) ) (setq prodval (cons 'vec prodval)) ) (cond ((and (eq blx 0) (eq bly 0)) (elt (elt prodval 1) 1)) ((and (< blx 2) (< bly 2)) (elt prodval 1)) (t prodval) ) ) ) ;; I don't know how to set this property ... this statement does not implement the Kronecker product (put 'calcFunckron 'calcuserdefn '(* (var x varx) (var y vary))) (definekey calcmodemap "zk" 'calckron) )) Windows Live Messenger 2008 vient de sortir, encore plus de fun ! Téléchargez gratuitement Messenger 2008 
[Prev in Thread]  Current Thread  [Next in Thread] 