| 2008-05-18 Michael Gogins
* Opcodes/chua/ChuaOscillator.cpp
Csound opcode implementation of Chua's Oscillator,
a nonlinear dynamical system capable of many kinds of periodic and chaotic
behavior, including multiple simultaneous attractors. See the
reference manual entry for "chuap" for documentation. See
examples/chuap.csd and examples/chuas_oscillator/ChuasOscillatorApp.py for
sample Csound usage.
* Opcodes/linear_algebra.cpp
Csound i-rate and k-rate opcodes for much standard linear algebra over
real and complex vectors and matrices: elementwise arithmetic, norms,
transpose and conjugate, inner products, matrix inverse, LU decomposition,
QR decomposition, and QR-based eigenvalue decomposition. Includes
copying vectors to and from a-rate signals, function tables, and f-signals.
See the source code for preliminary documentation and tests/test24.csd for
sample Csound usage.
See below for a listing of the linear algebra opcodes with some comments and notes.
I would appreciate it if Csound developers would review these linear algebra opcodes. I would particularly appreciate it if anything jumps out at you as being missing or wanted.
I have not exhaustively tested these opcodes, and in particular I have not yet tested them in use. I hope to develop a sample csd with some use cases in it soon. Suggestions for examples also would be appreciated.
There are a lot of these opcodes. They are not polymorphic. This is to make them run faster. Also, I did not define a new "signal" data type. That is because the vectors and matrices, although they can carry signals and operate upon signals, really are not signals and could be used for other purposes as well.
Regards,
Mike
/**
* L I N E A R A L G E B R A O P C O D E S F O R C S O U N D
* Michael Gogins
*
* These opcodes implement many linear algebra operations,
* from scalar, vector, and matrix arithmetic up to
* and including QR based eigenvalue decompositions.
* The opcodes are designed for digital signal processing,
* and of course other mathematical operations,
* in the Csound orchestra language.
*
* The numerical implementation uses the gmm++ library
* from http://home.gna.org/getfem/gmm_intro.
*
* NOTE: SDFT must be #defined in order to build and use these opcodes.
* For application with f-sig variables, arithmetic must be
* performed only when the f-sig is "current," because f-rate
* is some fraction of k-rate; currency can be determined with the
* la_k_current_f opcode.
*
* Linear Algebra Data Types
* -------------------------
*
* Mathematical Type Code Corresponding Csound Type or Types
* ----------------- ---- -------------------------------------------------
* real scalar r i-rate or k-rate variable
* complex scalar c pair of i-rate or k-rate variables, e.g. "kr, ki"
* real vector vr i-rate variable holding address of array
* a a-rate variable
* t function table number
* complex vector vc i-rate variable holding address of array
* f fsig variable
* real matrix mr i-rate variable holding address of array
* complex matrix mc i-rate variable holding address of array
*
* All arrays are 0-based; the first index iterates rows to give columns,
* the second index iterates columns to give elements.
*
* All arrays are general and dense; banded, Hermitian, symmetric
* and sparse routines are not implemented.
*
* An array can be of type code vr, vc, mr, or mc
* and is stored in an i-rate object.
* In orchestra code, an array is passed as a
* MYFLT i-rate variable that contains the address
* of the array object, which is actually stored
* in the allocator opcode instance.
* Although array variables are i-rate, of course
* their values and even shapes may change at i-rate or k-rate.
*
* All operands must be pre-allocated; except for the creation
* opcodes, no opcode ever allocates any arrays.
* This is true even if the array appears on the
* left-hand side of an opcode! However, some operations
* may reshape arrays to hold results.
*
* Arrays are automatically deallocated when their instrument
* is deallocated.
*
* Not only for more efficient performance,
* but also to make it easier to remember opcode names,
* the performance rate, output value types, operation names,
* and input value types are deterministically
* encoded into the opcode name:
* 1. "la" for "linear algebra opcode family".
* 2. "i" or "k" for performance rate.
* 3. Type code(s) (see above table) for output value(s),
* but only if the type is not implicit from the input values.
* 4. Operation name: common mathematical name
* (preferred) or abbreviation.
* 5. Type code(s) for input values, if not implicit.
*
* For details, see the gmm++ documentation at
* http://download.gna.org/getfem/doc/gmmuser.pdf.
*
* Array Creation
* --------------
*
* ivr la_i_vr_create irows
* ivc la_i_vc_create irows
* imr la_i_mr_create irows, icolumns [, odiagonal]
* imc la_i_mc_create irows, icolumns [, odiagonal_r, odiagonal_i]
*
* Array Introspection
* -------------------
*
* irows la_i_size_vr ivr
* irows la_i_size_vc ivc
* irows, icolumns la_i_size_mr imr
* irows, icolumns la_i_size_mc imc
*
* kfiscurrent la_k_current_f fsig
*
* la_i_print_vr ivr
* la_i_print_vc ivc
* la_i_print_mr imr
* la_i_print_mc imc
*
* Array Assignment and Conversion
* -------------------------------
*
* ivr la_i_assign_vr ivr
* ivr la_k_assign_vr ivr
* ivc la_i_assign_vc ivc
* ivc la_k_assign_vc ivr
* imr la_i_assign_mr imr
* imr la_k_assign_mr imr
* imc la_i_assign_mc imc
* imc la_k_assign_mc imr
*
* ivr la_k_assign_a asig
* ivr la_i_assign_t itablenumber
* ivr la_k_assign_t itablenumber
* ivc la_k_assign_f fsig
*
* asig la_k_a_assign ivr
* itablenum la_i_t_assign ivr
* itablenum la_k_t_assign ivr
* fsig la_k_f_assign ivc
*
* Fill Arrays with Random Elements
* --------------------------------
*
* ivr la_i_random_vr [ifill_fraction]
* ivr la_k_random_vr [kfill_fraction]
* ivc la_i_random_vc [ifill_fraction]
* ivc la_k_random_vc [kfill_fraction]
* imr la_i_random_mr [ifill_fraction]
* imr la_k_random_mr [kfill_fraction]
* imc la_i_random_mc [ifill_fraction]
* imc la_k_random_mc [kfill_fraction]
*
* Array Element Access
* --------------------
*
* ivr la_i_vr_set irow, ivalue
* kvr la_k_vr_set krow, kvalue
* ivc la_i_vc_set irow, ivalue_r, ivalue_i
* kvc la_k_vc_set krow, kvalue_r, kvalue_i
* imr la_i mr_set irow, icolumn, ivalue
* kmr la_k mr_set krow, kcolumn, ivalue
* imc la_i_mc_set irow, icolumn, ivalue_r, ivalue_i
* kmc la_k_mc_set krow, kcolumn, kvalue_r, kvalue_i
*
* ivalue la_i_get_vr ivr, irow
* kvalue la_k_get_vr ivr, krow,
* ivalue_r, ivalue_i la_i_get_vc ivc, irow
* kvalue_r, kvalue_i la_k_get_vc ivc, krow
* ivalue la_i_get_mr imr, irow, icolumn
* kvalue la_k_get_mr imr, krow, kcolumn
* ivalue_r, ivalue_i la_i_get_mc imc, irow, icolumn
* kvalue_r, kvalue_i la_k_get_mc imc, krow, kcolumn
*
* Single Array Operations
* -----------------------
*
* imr la_i_transpose_mr imr
* imr la_k_transpose_mr imr
* imc la_i_transpose_mc imc
* imc la_k_transpose_mc imc
* ivr la_i_conjugate_vr ivr
* ivr la_k_conjugate_vr ivr
* ivc la_i_conjugate_vc ivc
* ivc la_k_conjugate_vc ivc
* imr la_i_conjugate_mr imr
* imr la_k_conjugate_mr imr
* imc la_i_conjugate_mc imc
* imc la_k_conjugate_mc imc
*
* Scalar Operations
* -----------------
*
* ir la_i_norm1_vr ivr
* kr la_k_norm1_vr ivc
* ir la_i_norm1_vc ivc
* kr la_k_norm1_vc ivc
* ir la_i_norm1_mr imr
* kr la_k_norm1_mr imr
* ir la_i_norm1_mc imc
* kr la_k_norm1_mc imc
*
* ir la_i_norm_euclid_vr ivr
* kr la_k_norm_euclid_vr ivr
* ir la_i_norm_euclid_vc ivc
* kr la_k_norm_euclid_vc ivc
* ir la_i_norm_euclid_mr mvr
* kr la_k_norm_euclid_mr mvr
* ir la_i_norm_euclid_mc mvc
* kr la_k_norm_euclid_mc mvc
*
* ir la_i_distance_vr ivr
* kr la_k_distance_vr ivr
* ir la_i_distance_vc ivc
* kr la_k_distance_vc ivc
*
* ir la_i_norm_max imr
* kr la_k_norm_max imc
* ir la_i_norm_max imr
* kr la_k_norm_max imc
*
* ir la_i_norm_inf_vr ivr
* kr la_k_norm_inf_vr ivr
* ir la_i_norm_inf_vc ivc
* kr la_k_norm_inf_vc ivc
* ir la_i_norm_inf_mr imr
* kr la_k_norm_inf_mr imr
* ir la_i_norm_inf_mc imc
* kr la_k_norm_inf_mc imc
*
* ir la_i_trace_mr imr
* kr la_k_trace_mr imr
* ir, ii la_i_trace_mc imc
* kr, ki la_k_trace_mc imc
*
* ir la_i_lu_det imr
* kr la_k_lu_det imr
* ir la_i_lu_det imc
* kr la_k_lu_det imc
*
* Elementwise Array-Array Operations
* ----------------------------------
*
* ivr la_i_add_vr ivr_a, ivr_b
* ivc la_k_add_vc ivc_a, ivc_b
* imr la_i_add_mr imr_a, imr_b
* imc la_k_add_mc imc_a, imc_b
*
* ivr la_i_subtract_vr ivr_a, ivr_b
* ivc la_k_subtract_vc ivc_a, ivc_b
* imr la_i_subtract_mr imr_a, imr_b
* imc la_k_subtract_mc imc_a, imc_b
*
* ivr la_i_multiply_vr ivr_a, ivr_b
* ivc la_k_multiply_vc ivc_a, ivc_b
* imr la_i_multiply_mr imr_a, imr_b
* imc la_k_multiply_mc imc_a, imc_b
*
* ivr la_i_divide_vr ivr_a, ivr_b
* ivc la_k_divide_vc ivc_a, ivc_b
* imr la_i_divide_mr imr_a, imr_b
* imc la_k_divide_mc imc_a, imc_b
*
* Inner Products
* --------------
*
* ir la_i_dot_vr ivr_a, ivr_b
* kr la_k_dot_vr ivr_a, ivr_b
* ir, ii la_i_dot_vc ivc_a, ivc_b
* kr, ki la_k_dot_vc ivc_a, ivc_b
*
* imr la_i_dot_mr imr_a, imr_b
* imr la_k_dot_mr imr_a, imr_b
* imc la_i_dot_mc imc_a, imc_b
* imc la_k_dot_mc imc_a, imc_b
*
* ivr la_i_dot_mr_vr imr_a, ivr_b
* ivr la_k_dot_mr_vr imr_a, ivr_b
* ivc la_i_dot_mc_vc imc_a, ivc_b
* ivc la_k_dot_mc_vc imc_a, ivc_b
*
* Matrix Inversion
* ----------------
*
* imr, icondition la_i_invert_mr imr
* imr, kcondition la_k_invert_mr imr
* imc, icondition la_i_invert_mc imc
* imc, kcondition la_k_invert_mc imc
*
* Matrix Decompositions and Solvers
* ---------------------------------
*
* ivr la_i_upper_solve_mr imr [, j_1_diagonal]
* ivr la_k_upper_solve_mr imr [, j_1_diagonal]
* ivc la_i_upper_solve_mc imc [, j_1_diagonal]
* ivc la_k_upper_solve_mc imc [, j_1_diagonal]
*
* ivr la_i_lower_solve_mr imr [, j_1_diagonal]
* ivr la_k_lower_solve_mr imr [, j_1_diagonal]
* ivc la_i_lower_solve_mc imc [, j_1_diagonal]
* ivc la_k_lower_solve_mc imc [, j_1_diagonal]
*
* imr, ivr_pivot, isize la_i_lu_factor_mr imr
* imr, ivr_pivot, ksize la_k_lu_factor_mr imr
* imc, ivr_pivot, isize la_i_lu_factor_mc imc
* imc, ivr_pivot, ksize la_k_lu_factor_mc imc
*
* ivr_x la_i_lu_solve_mr imr, ivr_b
* ivr_x la_k_lu_solve_mr imr, ivr_b
* ivc_x la_i_lu_solve_mc imc, ivc_b
* ivc_x la_k_lu_solve_mc imc, ivc_b
*
* imr_q, imr_r la_i_qr_factor_mr imr
* imr_q, imr_r la_k_qr_factor_mr imr
* imc_q, imc_r la_i_qr_factor_mc imc
* imc_q, imc_r la_k_qr_factor_mc imc
*
* ivr_eig_vals la_i_qr_eigen_mr imr, i_tolerance
* ivr_eig_vals la_k_qr_eigen_mr imr, k_tolerance
* ivr_eig_vals la_i_qr_eigen_mc imc, i_tolerance
* ivr_eig_vals la_k_qr_eigen_mc imc, k_tolerance
*
* NOTE: Matrix must be Hermitian in order to compute eigenvectors.
*
* ivr_eig_vals, imr_eig_vecs la_i_qr_sym_eigen_mr imr, i_tolerance
* ivr_eig_vals, imr_eig_vecs la_k_qr_sym_eigen_mr imr, k_tolerance
* ivc_eig_vals, imc_eig_vecs la_i_qr_sym_eigen_mc imc, i_tolerance
* ivc_eig_vals, imc_eig_vecs la_k_qr_sym_eigen_mc imc, k_tolerance
*
*/
-------------------------------------------------------------------------
This SF.net email is sponsored by: Microsoft
Defy all challenges. Microsoft(R) Visual Studio 2008.
http://clk.atdmt.com/MRT/go/vse0120000070mrt/direct/01/
_______________________________________________
Csound-devel mailing list
Csound-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/csound-devel |