Max 5 API Reference
00001 /* 00002 Copyright 2001-2005 - Cycling '74 00003 Joshua Kit Clayton jkc@cycling74.com 00004 */ 00005 00006 #include "jit.common.h" 00007 #include "jit.fixmath.h" 00008 00009 #define JIT_LA_UPPERTRI_DEF_THRESH 0.000000001 00010 00011 typedef struct _jit_la_uppertri 00012 { 00013 t_object ob; 00014 double thresh; 00015 long swapcount; 00016 } t_jit_la_uppertri; 00017 00018 void *_jit_la_uppertri_class; 00019 00020 t_jit_la_uppertri *jit_la_uppertri_new(void); 00021 void jit_la_uppertri_free(t_jit_la_uppertri *x); 00022 t_jit_err jit_la_uppertri_matrix_calc(t_jit_la_uppertri *x, void *inputs, void *outputs); 00023 00024 void jit_la_uppertri_calculate_2d_real(t_jit_la_uppertri *x, long *swapcount, t_jit_matrix_info *in_minfo, char *bip); 00025 void jit_la_uppertri_calculate_2d_complex(t_jit_la_uppertri *x, long *swapcount, t_jit_matrix_info *in_minfo, char *bip); 00026 00027 t_jit_err jit_la_uppertri_init(void) 00028 { 00029 long attrflags=0; 00030 t_jit_object *attr,*mop,*o; 00031 t_atom a[2]; 00032 00033 _jit_la_uppertri_class = jit_class_new("jit_la_uppertri",(method)jit_la_uppertri_new,(method)jit_la_uppertri_free, 00034 sizeof(t_jit_la_uppertri),0L); 00035 00036 //add mop 00037 mop = jit_object_new(_jit_sym_jit_mop,1,1); 00038 o = jit_object_method(mop,_jit_sym_getoutput,1); 00039 jit_attr_setlong(o,_jit_sym_maxplanecount,2); 00040 jit_atom_setsym(a,_jit_sym_float32); //default 00041 jit_atom_setsym(a+1,_jit_sym_float64); 00042 jit_object_method(o,_jit_sym_types,2,a); 00043 jit_class_addadornment(_jit_la_uppertri_class,mop); 00044 //add methods 00045 jit_class_addmethod(_jit_la_uppertri_class, (method)jit_la_uppertri_matrix_calc, "matrix_calc", A_CANT, 0L); 00046 //add attributes 00047 attrflags = JIT_ATTR_GET_DEFER_LOW | JIT_ATTR_SET_USURP_LOW; 00048 attr = jit_object_new(_jit_sym_jit_attr_offset,"thresh",_jit_sym_float64,attrflags, 00049 (method)0L,(method)0L,calcoffset(t_jit_la_uppertri,thresh)); 00050 jit_class_addattr(_jit_la_uppertri_class,attr); 00051 attrflags = JIT_ATTR_SET_OPAQUE_USER | JIT_ATTR_GET_DEFER_LOW; 00052 attr = jit_object_new(_jit_sym_jit_attr_offset,"swapcount",_jit_sym_long,attrflags, 00053 (method)0L,(method)0L,calcoffset(t_jit_la_uppertri,swapcount)); 00054 jit_class_addattr(_jit_la_uppertri_class,attr); 00055 00056 jit_class_register(_jit_la_uppertri_class); 00057 00058 return JIT_ERR_NONE; 00059 } 00060 00061 t_jit_err jit_la_uppertri_matrix_calc(t_jit_la_uppertri *x, void *inputs, void *outputs) 00062 { 00063 t_jit_err err=JIT_ERR_NONE; 00064 long in_savelock,out_savelock; 00065 t_jit_matrix_info in_minfo,out_minfo; 00066 char *in_bp,*out_bp; 00067 long i,dimcount,planecount,dim[JIT_MATRIX_MAX_DIMCOUNT]; 00068 long tmpsize; 00069 long *swapcount = &x->swapcount; 00070 void *in_matrix,*out_matrix; 00071 00072 in_matrix = jit_object_method(inputs,_jit_sym_getindex,0); 00073 out_matrix = jit_object_method(outputs,_jit_sym_getindex,0); 00074 00075 if (x&&in_matrix&&out_matrix&&swapcount) { 00076 *swapcount=0; 00077 in_savelock = (long) jit_object_method(in_matrix,_jit_sym_lock,1); 00078 out_savelock = (long) jit_object_method(out_matrix,_jit_sym_lock,1); 00079 jit_object_method(in_matrix,_jit_sym_getinfo,&in_minfo); 00080 jit_object_method(out_matrix,_jit_sym_getinfo,&out_minfo); 00081 jit_object_method(in_matrix,_jit_sym_getdata,&in_bp); 00082 jit_object_method(out_matrix,_jit_sym_getdata,&out_bp); 00083 if (!in_bp) { err=JIT_ERR_INVALID_INPUT; goto out;} 00084 if (!out_bp) { err=JIT_ERR_INVALID_OUTPUT; goto out;} 00085 //compatible types? 00086 if (((in_minfo.type!=_jit_sym_float32)&&(in_minfo.type!=_jit_sym_float64))|| 00087 (in_minfo.type!=out_minfo.type)) { err=JIT_ERR_MISMATCH_TYPE; goto out;} 00088 //compatible planecounts? 00089 //planecount 1=real. planecount 2=complex 00090 if (((in_minfo.planecount!=1)&&(in_minfo.planecount!=2))|| 00091 (in_minfo.planecount!=out_minfo.planecount)) { err=JIT_ERR_MISMATCH_PLANE; goto out;} 00092 //check dimensions 00093 if ((in_minfo.dim[0]!=out_minfo.dim[0])|| 00094 (in_minfo.dim[1]!=out_minfo.dim[1])) { err=JIT_ERR_MISMATCH_DIM; goto out;} 00095 //get dimensions/planecount 00096 dimcount = out_minfo.dimcount; 00097 planecount = out_minfo.planecount; 00098 00099 //only 2d 00100 00101 //first copy matrix, and then calculate upper trianglular matrix inplace on the output buffer 00102 jit_object_method(out_matrix,_jit_sym_frommatrix,in_matrix,NULL); 00103 00104 if (planecount==1) 00105 jit_la_uppertri_calculate_2d_real(x, swapcount, &out_minfo, out_bp); 00106 else 00107 jit_la_uppertri_calculate_2d_complex(x, swapcount, &out_minfo, out_bp); 00108 } else { 00109 return JIT_ERR_INVALID_PTR; 00110 } 00111 00112 out: 00113 jit_object_method(out_matrix,_jit_sym_lock,out_savelock); 00114 jit_object_method(in_matrix,_jit_sym_lock,in_savelock); 00115 return err; 00116 } 00117 00118 void jit_la_uppertri_calculate_2d_real(t_jit_la_uppertri *x, long *swapcount, t_jit_matrix_info *in_minfo, char *bip) 00119 { 00120 long i,j,n,rows; 00121 t_jit_op_info in_opinfo,tmp_opinfo; 00122 double val,val2,val_inv; 00123 char *op; 00124 00125 n = in_minfo->dim[0]; 00126 if (in_minfo->type==_jit_sym_float32) { 00127 tmp_opinfo.stride = in_opinfo.stride = in_minfo->dimstride[0]/4; //row vector 00128 rows = in_minfo->dim[1]; 00129 for (i=0;i<rows;i++) { 00130 in_opinfo.p = bip + i*in_minfo->dimstride[1]; 00131 val = *((float *)(((char *)in_opinfo.p) + i*in_minfo->dimstride[0])); 00132 //perform row swapping if necessary 00133 if (ABS(val)<x->thresh) { 00134 j=i+1; 00135 while ((j<rows)&&(ABS(val)<x->thresh)) { 00136 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00137 val = *((float *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0])); 00138 j++; 00139 } 00140 if (!(ABS(val)<x->thresh)) { 00141 jit_op_vector_swap_float32(n,NULL,&in_opinfo,&tmp_opinfo,NULL); 00142 *swapcount += 1; 00143 } 00144 //else no pivots 00145 } 00146 //eliminate all non zero vectors in column i 00147 if (!(ABS(val)<x->thresh)) { // thresh? 00148 j=i+1; 00149 val_inv=-1./val; 00150 while (j<rows) { 00151 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00152 val2 = *((float *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0])); 00153 if (!(ABS(val2)<x->thresh)) { // thresh? 00154 val2 = val2*val_inv; 00155 jit_op_vector_axpy_float32(n,&val2,&in_opinfo,&tmp_opinfo,&tmp_opinfo); // note in place 00156 } 00157 j++; 00158 } 00159 } 00160 } 00161 } else if (in_minfo->type==_jit_sym_float64) { 00162 tmp_opinfo.stride = in_opinfo.stride = in_minfo->dimstride[0]/8; //row vector 00163 rows = in_minfo->dim[1]; 00164 for (i=0;i<rows;i++) { 00165 in_opinfo.p = bip + i*in_minfo->dimstride[1]; 00166 val = *((double *)(((char *)in_opinfo.p) + i*in_minfo->dimstride[0])); 00167 //perform row swapping if necessary 00168 if ((ABS(val)<x->thresh)) { 00169 j=i+1; 00170 while ((j<rows)&&(ABS(val)<x->thresh)) { 00171 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00172 val = *((double *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0])); 00173 j++; 00174 } 00175 if (!(ABS(val)<x->thresh)) { 00176 jit_op_vector_swap_float64(n,NULL,&in_opinfo,&tmp_opinfo,NULL); 00177 *swapcount += 1; 00178 } 00179 //else no pivots 00180 } 00181 //eliminate all non zero vectors in column i 00182 if (!(ABS(val)<x->thresh)) { // thresh? 00183 j=i+1; 00184 val_inv=-1./val; 00185 while (j<rows) { 00186 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00187 val2 = *((double *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0])); 00188 if (!(ABS(val2)<x->thresh)) { // thresh? 00189 val2 = val2*val_inv; 00190 jit_op_vector_axpy_float64(n,&val2,&in_opinfo,&tmp_opinfo,&tmp_opinfo); // note in place 00191 } 00192 j++; 00193 } 00194 } 00195 } 00196 } 00197 } 00198 00199 #define COMPLEX_IS_ZERO(a,b) (((a)==0)&&((b)==0)) 00200 #define COMPLEX_IS_ZERO_THRESH(a,b,thresh) ((ABS(a)<(thresh))&&(ABS(b)<(thresh))) 00201 #define COMPLEX_RECIP_REAL(a,b) ((a)/(((a)*(a))+((b)*(b)))) 00202 #define COMPLEX_RECIP_IMAG(a,b) (-(b)/(((a)*(a))+((b)*(b)))) 00203 #define COMPLEX_MULT_REAL(a,b,c,d) (((a)*(c))-((b)*(d))) 00204 #define COMPLEX_MULT_IMAG(a,b,c,d) (((a)*(d))+((b)*(c))) 00205 00206 void jit_la_uppertri_calculate_2d_complex(t_jit_la_uppertri *x, long *swapcount, t_jit_matrix_info *in_minfo, char *bip) 00207 { 00208 long i,j,n,rows; 00209 t_jit_op_info in_opinfo,tmp_opinfo; 00210 double val[2],val2[2],val3[2],val_inv[2]; 00211 char *op; 00212 00213 n = in_minfo->dim[0]; 00214 if (in_minfo->type==_jit_sym_float32) { 00215 tmp_opinfo.stride = in_opinfo.stride = in_minfo->dimstride[0]/4; //row vector 00216 rows = in_minfo->dim[1]; 00217 for (i=0;i<rows;i++) { 00218 in_opinfo.p = bip + i*in_minfo->dimstride[1]; 00219 val[0] = ((float *)(((char *)in_opinfo.p) + i*in_minfo->dimstride[0]))[0]; 00220 val[1] = ((float *)(((char *)in_opinfo.p) + i*in_minfo->dimstride[0]))[1]; 00221 00222 //perform row swapping if necessary 00223 if (COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { 00224 j=i+1; 00225 while ((j<rows)&&COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { 00226 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00227 val[0] = ((float *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[0]; 00228 val[1] = ((float *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[1]; 00229 j++; 00230 } 00231 if (!COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { 00232 jit_op_vector_swap_float32_complex(n,NULL,&in_opinfo,&tmp_opinfo,NULL); 00233 *swapcount += 1; 00234 } 00235 //else no pivots 00236 } 00237 //eliminate all non zero vectors in column i 00238 if (!COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { // thresh? 00239 j=i+1; 00240 val_inv[0] = -COMPLEX_RECIP_REAL(val[0],val[1]); 00241 val_inv[1] = -COMPLEX_RECIP_IMAG(val[0],val[1]); 00242 while (j<rows) { 00243 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00244 val2[0] = ((float *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[0]; 00245 val2[1] = ((float *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[1]; 00246 if (!COMPLEX_IS_ZERO_THRESH(val2[0],val2[1],x->thresh)) { // thresh? 00247 val3[0] = COMPLEX_MULT_REAL(val2[0],val2[1],val_inv[0],val_inv[1]); 00248 val3[1] = COMPLEX_MULT_IMAG(val2[0],val2[1],val_inv[0],val_inv[1]); 00249 jit_op_vector_axpy_float32_complex(n,val3,&in_opinfo,&tmp_opinfo,&tmp_opinfo); // note in place 00250 } 00251 j++; 00252 } 00253 } 00254 } 00255 } else if (in_minfo->type==_jit_sym_float64) { 00256 tmp_opinfo.stride = in_opinfo.stride = in_minfo->dimstride[0]/8; //row vector 00257 rows = in_minfo->dim[1]; 00258 for (i=0;i<rows;i++) { 00259 in_opinfo.p = bip + i*in_minfo->dimstride[1]; 00260 val[0] = ((double *)(((char *)in_opinfo.p) + i*in_minfo->dimstride[0]))[0]; 00261 val[1] = ((double *)(((char *)in_opinfo.p) + i*in_minfo->dimstride[0]))[1]; 00262 00263 //perform row swapping if necessary 00264 if (COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { 00265 j=i+1; 00266 while ((j<rows)&&COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { 00267 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00268 val[0] = ((double *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[0]; 00269 val[1] = ((double *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[1]; 00270 j++; 00271 } 00272 if (!COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { 00273 jit_op_vector_swap_float64_complex(n,NULL,&in_opinfo,&tmp_opinfo,NULL); 00274 *swapcount += 1; 00275 } 00276 //else no pivots 00277 } 00278 //eliminate all non zero vectors in column i 00279 if (!COMPLEX_IS_ZERO_THRESH(val[0],val[1],x->thresh)) { // thresh? 00280 j=i+1; 00281 val_inv[0] = -COMPLEX_RECIP_REAL(val[0],val[1]); 00282 val_inv[1] = -COMPLEX_RECIP_IMAG(val[0],val[1]); 00283 while (j<rows) { 00284 tmp_opinfo.p = bip + j*in_minfo->dimstride[1]; //row j 00285 val2[0] = ((double *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[0]; 00286 val2[1] = ((double *)(((char *)tmp_opinfo.p) + i*in_minfo->dimstride[0]))[1]; 00287 if (!COMPLEX_IS_ZERO_THRESH(val2[0],val2[1],x->thresh)) { // thresh? 00288 val3[0] = COMPLEX_MULT_REAL(val2[0],val2[1],val_inv[0],val_inv[1]); 00289 val3[1] = COMPLEX_MULT_IMAG(val2[0],val2[1],val_inv[0],val_inv[1]); 00290 jit_op_vector_axpy_float64_complex(n,val3,&in_opinfo,&tmp_opinfo,&tmp_opinfo); // note in place 00291 } 00292 j++; 00293 } 00294 } 00295 } 00296 } 00297 } 00298 00299 00300 t_jit_la_uppertri *jit_la_uppertri_new(void) 00301 { 00302 t_jit_la_uppertri *x; 00303 00304 if (x=(t_jit_la_uppertri *)jit_object_alloc(_jit_la_uppertri_class)) { 00305 x->thresh = JIT_LA_UPPERTRI_DEF_THRESH; 00306 } else { 00307 x = NULL; 00308 } 00309 return x; 00310 } 00311 00312 void jit_la_uppertri_free(t_jit_la_uppertri *x) 00313 { 00314 }
Copyright © 2008, Cycling '74