I would like to know the right way to compute the products of the same n by m matrix A by around a hundred m by 1 column vectors in Nyquist.
In my use case, n << m, the column vectors are sections m samples long taken from the selected track using the extract function. A, when calculated, is an array of arrays of floats.
What I have so far:
Immediately after calculating A, I convert it to a multichannel sound like so:
(defun mat-to-snds (A &optional (n (length A)))
(dotimes (i n A) (setf (aref A i) (snd-from-array 0 *sound-srate* (aref A i)))))
Then,
(setq t-last (/ (float (1- m)) (get-duration 1.0) *sound-srate*))
(defun mv-mult (A b)
(multichan-expand
#'(lambda (sig) (sref sig t-last))
(scale *sound-srate* (integrate (mult A b)))))
My questions are:
What is the value of t-last such that sref will give me the very last sample in a sound with m samples, start time 0, and sample rate sound-srate within a process-type plugin? What’s written above seems to give me the next-to-last sample and I’m not sure how to fix that.
Does mv-mult calculate any samples besides the one I want?
Does mv-mult overwrite A?
Is there any advantage to using sounds here, or can I expect to gain performance by using arrays of arrays, and writing the loops in lisp?
A is a matrix. Usually A would be the transpose of the “Q” component of the QR factorization of some other matrix B. Calculating this factor may also require several intermediate matrix multiplications, so I would likely be passing other matrices to the same function.
The biggest thing I’d like to do with this is wavetable resynthesis, the reconstruction of TRACK using amplitude envelopes applied to a small number of looped wavetables. The least squares solutions x of Bx = b - where B is an m x n matrix whose columns are wavetables length m (or the ffts thereof), and the b are sections (or possibly, fft windows) of the target sound - gives the best possible amplitude envelope for a given set of wavetables.
MULT does exactly what I need, thank you! No clue why I hadn’t tried it.
Hang on, mult isn’t what I’m looking for, sorry. What that does is componentwise multiplication, but I want to do matrix multiplication.
Let’s say the matrix A is represented as an array of n arrays, each length m, and the vector b is an array length n. The result should be an array length m, where the _i_th entry should be the sum of the componentwise multiplication of (aref A i) and b.
For example,
(print (mv-mult
#(#(2 0 6)
#(0 2 9))
#(5 2 1)))
should print #(16 13).
Below is how I’ve done this so far - it works fine, but if there is a Nyquist built-in available that can do it faster, I would prefer to use that.
;;;Row-major matrix A times vector b.
(defun mv-mult (A b &optional (m (length A)) (n (length b)))
(do ((rst (make-array m))
(i (1- m) (1- i)))
((< i 0) rst)
(setf (aref rst i)
(do ((entry 0)
(row (aref A i))
(j (1- n) (1- j)))
((< j 0) entry)
(setq entry (sum entry
(* (aref row j) (aref b j))))))))
As far as I’m aware, there isn’t a matrix multiplication operator in Nyquist.
That looks fine assuming that we only need to multiply by a N x 1 matrix.
I would probably have written it slightly differently, but only because I find this version a little easier to read:
(defun mxmult (A B)
;; Matrix multiplication, where A is matrix MxN and B is Nx1.
(let* ((M (length A))
(result (make-array M)))
(dotimes (i M result)
(let ((AiB (mult (aref A i) B)))
(setf (aref result i)(vsum AiB))))))
(defun vsum (v &aux (total 0))
;;Return the sum of elements in vector v.
(dotimes (i (length v) total)
(setf total (+ total (aref v i)))))
or for a condensed single function version:
(defun mxmult (A B &aux (M (length A))
(result (make-array (length A))))
(flet ((vsum (v &aux (total 0))
(dotimes (i (length v) total)
(setf total (+ total (aref v i))))))
(dotimes (i M result)
(let ((AiB (mult (aref A i) B)))
(setf (aref result i)(vsum AiB))))))
Personally I prefer the first version as I find it easiest to read, but your version and both of these versions work.