COMMENT
ENDCOMMENT
NEURON {
SUFFIX nothing
}
VERBATIM
#ifndef NRN_VERSION_GTEQ_8_2_0
extern double hoc_call_func(Symbol*, int narg);
#endif
ENDVERBATIM
PARAMETER {
MATRIX_INSTALLED=0
}
VERBATIM
static double outprod(void* vv) {
int i, j, nx, ny, nz;
double *x, *y, *z;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
nz = vector_arg_px(2, &z);
if (nx != ny*nz) {
hoc_execerror("Vector size mismatch", 0);
}
for (i=0;i<ny;i++) {
for (j=0;j<nz;j++) {
x[i*nz+j] = y[i]*z[j];
}
}
return nx;
}
ENDVERBATIM
VERBATIM
static double mmult(void* vv) {
int i, j, nx, ny, nz;
double *x, *y, *z;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
nz = vector_arg_px(2, &z);
if (ny != nx*nz) {
hoc_execerror("Vector size mismatch", 0);
}
for (i=0;i<nx;i++) {
x[i] = 0.;
for (j=0;j<nz;j++) {
x[i] += y[i*nz+j]*z[j];
}
}
return nx;
}
ENDVERBATIM
VERBATIM
static double spltp(void* vv) {
int ii, jj, nstpr, nstpo, nw, npr, npo, flag, cnt;
double *stpr, *stpo, *w, *pr, *po;
char func[4] = "ltp";
Symbol* s = hoc_lookup(func);
if (! s) { hoc_execerror("Can't find ltp() func", 0); }
nstpo = vector_instance_px(vv, &stpo);
npr = vector_arg_px(1, &pr);
npo = vector_arg_px(2, &po);
nw = vector_arg_px(3, &w);
nstpr = vector_arg_px(4, &stpr);
for (ii=0,jj=0,cnt=0;ii<nstpo;ii++) {
if (stpo[ii]==1.0) {
for (;po[jj]<ii;jj++) ;
for (;po[jj]==ii;jj++) {
if (stpr[(int)pr[jj]]==1.) {
cnt++; hoc_pushx(1.0);
} else {
cnt--; hoc_pushx(-1.0);
}
hoc_pushx(w[jj]);
w[jj]=hoc_call_func(s, 2);
}
}
}
return cnt;
}
ENDVERBATIM
VERBATIM
static int *pr_int;
static int *po_int;
static int cpfl=0;
ENDVERBATIM
VERBATIM
static double mkspcp(void* vv) {
int j, nw, npr, npo;
double *w, *pr, *po;
if (! ifarg(1)) {
cpfl=0;
if (po_int!=NULL) free(po_int);
if (pr_int!=NULL) free(pr_int);
po_int=(int *)NULL; pr_int=(int *)NULL;
return 0;
}
nw = vector_instance_px(vv, &w);
npr = vector_arg_px(1, &pr);
npo = vector_arg_px(2, &po);
pr_int=(int *)ecalloc(nw, sizeof(int));
po_int=(int *)ecalloc(nw, sizeof(int));
for (j=0;j<nw;j++) {
po_int[j]=(int)po[j];
pr_int[j]=(int)pr[j];
}
cpfl=nw;
return cpfl;
}
ENDVERBATIM
VERBATIM
static double chkspcp(void* vv) {
int j, nw, npr, npo, flag;
double *w, *pr, *po;
nw = vector_instance_px(vv, &w);
npr = vector_arg_px(1, &pr);
npo = vector_arg_px(2, &po);
flag=1;
if (po_int==NULL || pr_int==NULL) { cpfl=0; return 0; }
if (cpfl!=nw) { flag=0;
} else for (j=0;j<nw;j++) {
if (po_int[j]!=(int)po[j] || pr_int[j]!=(int)pr[j]) {flag=0; continue;}
}
if (flag==0) {
cpfl=0; free(po_int); free(pr_int);
po_int=(int *)NULL; pr_int=(int *)NULL;
}
return flag;
}
ENDVERBATIM
VERBATIM
static double spmult(void* vv) {
int i, j, nx, ny, nw, npr, npo, flag;
double *x, *y, *w, *pr, *po, xx;
ny = vector_instance_px(vv, &y);
npr = vector_arg_px(1, &pr);
npo = vector_arg_px(2, &po);
nw = vector_arg_px(3, &w);
nx = vector_arg_px(4, &x);
if (ifarg(5)) {flag=1;} else {flag=0;}
if (nw!=npr || nw!=npo) {
hoc_execerror("Sparse mat must have 3 identical size vecs for pre/post/wt", 0);
}
if (flag==0) for (i=0;i<ny;i++) y[i] = 0.;
if (cpfl==0) {
for (j=0;j<nw;j++) y[(int)po[j]] += (x[(int)pr[j]]*w[j]);
} else if (cpfl!=nw) { hoc_execerror("cpfl!=nw in spmult", 0); } else {
for (j=0;j<nw;j++) if (x[pr_int[j]]!=0) { y[po_int[j]] += ((x[pr_int[j]])*w[j]); }
}
return nx;
}
ENDVERBATIM
VERBATIM
static double spget(void* vv) {
int j, nw, npr, npo;
double *w, *pr, *po, row, col;
nw = vector_instance_px(vv, &w);
npr = vector_arg_px(1, &pr);
npo = vector_arg_px(2, &po);
row = *getarg(3);
col = *getarg(4);
for (j=0;j<nw;j++) if (row==po[j]&&col==pr[j]) break;
if (j==nw) return 0.; else return w[j];
}
ENDVERBATIM
FUNCTION spidget() { VERBATIM
{
int j, npr, npo, nprid, npoid;
double *pr, *po, *prid, *poid, pri, poi, pridi, poidi;
npr = vector_arg_px(1, &pr);
npo = vector_arg_px(2, &po);
nprid = vector_arg_px(3, &prid);
npoid = vector_arg_px(4, &poid);
pri= *getarg(5);
poi= *getarg(6);
pridi= *getarg(7);
poidi= *getarg(8);
for (j=0;j<npr;j++) {
if (poi==po[j]&&pri==pr[j]&&pridi==prid[j]&&poidi==poid[j]) break;
}
if (j==npr) _lspidget=-1.0; else _lspidget=(double)j;
}
ENDVERBATIM
}
VERBATIM
static double transpose(void* vv) {
int i, j, nx, ny, rows, cols;
double *x, *y;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
rows = (int)*getarg(2);
cols = (int)*getarg(3);
if (ny != nx) {
hoc_execerror("Vector size mismatch", 0);
}
for (i=0;i<rows;i++) {
for (j=0;j<cols;j++) {
x[j*rows+i] = y[i*cols+j];
}
}
return nx;
}
ENDVERBATIM
VERBATIM
static double mprintf(void* vv) {
int i, j, nx, rows, cols;
double *x;
nx = vector_instance_px(vv, &x);
rows = (int)*getarg(1);
cols = (int)*getarg(2);
if (nx != rows*cols) {
hoc_execerror("Vector size mismatch", 0);
}
for (i=0;i<rows;i++) {
for (j=0;j<cols;j++) {
printf("%g\t",x[i*cols+j]);
}
printf("\n");
}
return nx;
}
ENDVERBATIM
VERBATIM
static double mget(void* vv) {
int i, j, nx, rows, cols;
double *x;
nx = vector_instance_px(vv, &x);
i = (int)*getarg(1);
j = (int)*getarg(2);
cols = (int)*getarg(3);
if (i*cols+j >= nx) {
hoc_execerror("Indices out of bounds", 0);
}
return x[i*cols+j];
}
ENDVERBATIM
VERBATIM
static double mrow(void* vv) {
int i, j, nx, ny, rows, cols;
double *x, *y;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
i = (int)*getarg(2);
cols = (int)*getarg(3); rows=ny/cols;
if (cols!=nx) {
nx=vector_buffer_size((IvocVect*)vv);
if (cols<=nx) {
vector_resize((IvocVect*)vv, cols); nx=cols;
} else {
printf("%d > %d :: \n",cols,nx);
hoc_execerror("Vector max capacity too small in mrow", 0);
}
}
if (i>=rows) { hoc_execerror("Indices out of bounds", 0); }
for (j=0;j<nx;j++) { x[j] = y[i*cols+j]; }
return nx;
}
ENDVERBATIM
VERBATIM
static double mcol(void* vv) {
int i, j, nx, ny, rows, cols;
double *x, *y;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
j = (int)*getarg(2);
cols = (int)*getarg(3); rows=ny/cols;
if (rows!=nx) {
nx=vector_buffer_size((IvocVect*)vv);
if (rows<=nx) {
vector_resize((IvocVect*)vv, rows); nx=rows;
} else {
printf("%d > %d :: ",rows,nx);
hoc_execerror("Vector max capacity too small in mcol ", 0);
}
}
if (j>=cols) { hoc_execerror("Indices out of bounds", 0); }
for (i=0;i<nx;i++) { x[i] = y[i*cols+j]; }
return nx;
}
ENDVERBATIM
VERBATIM
static double msetrow(void* vv) {
int i, j, nx, ny, rows, cols;
double *x, *y;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
i = (int)*getarg(2);
cols = (int)*getarg(3);
if (cols!=nx || i>=ny/cols) {
hoc_execerror("Indices out of bounds", 0);
}
for (j=0;j<nx;j++) { y[i*cols+j] = x[j]; }
return nx;
}
ENDVERBATIM
VERBATIM
static double msetcol(void* vv) {
int i, j, nx, ny, rows, cols;
double *x, *y;
nx = vector_instance_px(vv, &x);
ny = vector_arg_px(1, &y);
j = (int)*getarg(2);
cols = (int)*getarg(3);
if (cols!=ny/nx || j>=cols) {
hoc_execerror("Indices out of bounds", 0);
}
for (i=0;i<nx;i++) { y[i*cols+j] = x[i]; }
return nx;
}
ENDVERBATIM
VERBATIM
static double mset(void* vv) {
int i, j, nx, rows, cols;
double *x, val;
nx = vector_instance_px(vv, &x);
i = (int)*getarg(1);
j = (int)*getarg(2);
cols = (int)*getarg(3);
val = *getarg(4);
if (i*cols+j >= nx) {
hoc_execerror("Indices out of bounds", 0);
}
return (x[i*cols+j]=val);
}
ENDVERBATIM
PROCEDURE install_matrix() {
MATRIX_INSTALLED=1
VERBATIM
install_vector_method("outprod", outprod);
install_vector_method("mmult", mmult);
install_vector_method("spmult", spmult);
install_vector_method("spget", spget);
install_vector_method("mkspcp", mkspcp);
install_vector_method("chkspcp", chkspcp);
install_vector_method("spltp", spltp);
install_vector_method("transpose", transpose);
install_vector_method("mprintf", mprintf);
install_vector_method("mget", mget);
install_vector_method("mset", mset);
install_vector_method("mrow", mrow);
install_vector_method("mcol", mcol);
install_vector_method("msetrow", msetrow);
install_vector_method("msetcol", msetcol);
ENDVERBATIM
}