Bass funct.ck

//Bass Funct: a 2D Linear Dynamical System
//watson 2008 w.metaph.org (replace first . with @)
//please email me with comments/ideas!

//comment: note this would not be hard at all to extend to higher dimensions.

SinOsc s1 => dac.left;
SinOsc s2 => dac.right;

//initial condition
[1.0,-1.0] @=> float x_this[];

//frequency scaling: scale the state variables by this much before => freq
40 => float freq_scale; 

//time step duration (numerical integration)
.001 => float delta_t;

//time dilation (how to comparitively step through real time)
1 => float time_dil;

//extra parameters
[-10.0] @=> float r[]; //static parameter value
//dynamic periodic parameter
SinOsc param_osc => blackhole; 
.01 => param_osc.freq; 

//A matrix: deterimines the dynamical system behavior
//"Nonlinear Dynamics and Chaos" by Strogatz is a great place to start
[ -0.3 , 0.15 , .4, -.6 ] @=> float A[]; 

while ((delta_t*time_dil)::second => now){

        //<<< x_this[0] >>>;
        //<<< x_this[1] >>>;
        x_next ( x_this, A, r, delta_t) @=> x_this; 
        x_this[0] * freq_scale => s1.freq;
        x_this[1] * freq_scale => s2.freq;
        
        //      <<<r[0]>>> + .003 => r[0];
        param_osc.last() => r[0];
}

//next x using runge-kata numerical integration
fun float[] x_next( float x_last[], float A[], float r[], float delta_t ) {
   f(x_last, A, r) @=> float test[];
   mult_array ( f(x_last, A, r), delta_t) @=> float k1[];
   mult_array ( f( add_array(x_last, div_array(k1,2)), A, r), delta_t ) @=> float k2[];
   mult_array ( f( add_array(x_last, div_array(k2,2)), A, r), delta_t ) @=> float k3[];
   mult_array ( f( add_array(x_last, k3), A, r), delta_t ) @=> float k4[];
   return add_array( x_last, div_array( add_array(k1, add_array( mult_array(k2,2), add_array(
   mult_array(k3,2), k4) ) ), 6 ) ); 
}

//dx/dt
//note: currently assumes a 2D vector from x
fun float[] f( float x[], float A[], float r[]) { 
        float x_new[2];
    A[0] * x[0] + A[1] * x[1] + r[0] => x_new[0];
    A[2] * x[0] + A[3] * x[1] => x_new[1];      
        return x_new;
} 

//ARRAY UTILITIES *****************

fun float[] add_array( float a1[], float a2[] ) {

        float sum[2];
        for( 0 => int i; i < a1.cap(); i++ ) {
                a1[i] + a2[i] => sum[i];
        }
        
        return sum;
}


fun float[] div_array( float a[], float d ) {

        float quotient[2];
        for( 0 => int i; i < a.cap(); i++ ) {
                a[i] / d => quotient[i];
        }
        
        return quotient;

}

fun float[] mult_array( float a[], float d ) {


        float prod[2];
        for( 0 => int i; i < a.cap(); i++ ) {
                a[i] * d => prod[i];
        }
        
        return prod;

}