SHOGUN  v3.2.0
Math.h
Go to the documentation of this file.
1 /*
2  * This program is free software; you can redistribute it and/or modify
3  * it under the terms of the GNU General Public License as published by
4  * the Free Software Foundation; either version 3 of the License, or
5  * (at your option) any later version.
6  *
7  * Written (W) 2013 Thoralf Klein
8  * Written (W) 2013 Soumyajit De
9  * Written (W) 2012 Fernando Jose Iglesias Garcia
10  * Written (W) 2011 Siddharth Kherada
11  * Written (W) 2011 Justin Patera
12  * Written (W) 2011 Alesis Novik
13  * Written (W) 2011-2013 Heiko Strathmann
14  * Written (W) 1999-2009 Soeren Sonnenburg
15  * Written (W) 1999-2008 Gunnar Raetsch
16  * Written (W) 2007 Konrad Rieck
17  * Copyright (C) 1999-2009 Fraunhofer Institute FIRST and Max-Planck-Society
18  */
19 
20 #ifndef __MATHEMATICS_H_
21 #define __MATHEMATICS_H_
22 
23 #include <shogun/base/SGObject.h>
24 #include <shogun/lib/common.h>
25 #include <shogun/io/SGIO.h>
26 #include <shogun/base/Parallel.h>
28 
29 #ifndef _USE_MATH_DEFINES
30 #define _USE_MATH_DEFINES
31 #endif
32 
33 #include <math.h>
34 #include <stdio.h>
35 #include <float.h>
36 #include <sys/types.h>
37 #ifndef _WIN32
38 #include <unistd.h>
39 #endif
40 
41 #ifdef HAVE_PTHREAD
42 #include <pthread.h>
43 #endif
44 
45 #ifdef SUNOS
46 #include <ieeefp.h>
47 #endif
48 
50 #ifdef log2
51 #define cygwin_log2 log2
52 #undef log2
53 #endif
54 
55 #ifndef M_PI
56 #define M_PI 3.14159265358979323846
57 #endif
58 
59 #ifdef _WIN32
60 #ifndef isnan
61 #define isnan _isnan
62 #endif
63 
64 #ifndef isfinite
65 #define isfinite _isfinite
66 #endif
67 #endif //_WIN32
68 
69 /* Size of RNG seed */
70 #define RNG_SEED_SIZE 256
71 
72 /* Maximum stack size */
73 #define RADIX_STACK_SIZE 512
74 
75 /* Stack macros */
76 #define radix_push(a, n, i) sp->sa = a, sp->sn = n, (sp++)->si = i
77 #define radix_pop(a, n, i) a = (--sp)->sa, n = sp->sn, i = sp->si
78 
79 #ifndef DOXYGEN_SHOULD_SKIP_THIS
80 
81 template <class T> struct radix_stack_t
82 {
84  T *sa;
86  size_t sn;
88  uint16_t si;
89 };
90 
92 
94 template <class T1, class T2> struct thread_qsort
95 {
97  T1* output;
99  T2* index;
101  uint32_t size;
102 
104  int32_t* qsort_threads;
106  int32_t sort_limit;
108  int32_t num_threads;
109 };
110 #endif // DOXYGEN_SHOULD_SKIP_THIS
111 
112 #define COMPLEX128_ERROR_ONEARG(function) \
113 static inline complex128_t function(complex128_t a) \
114 { \
115  SG_SERROR("CMath::%s():: Not supported for complex128_t\n",\
116  #function);\
117  return complex128_t(0.0, 0.0); \
118 }
119 
120 #define COMPLEX128_STDMATH(function) \
121 static inline complex128_t function(complex128_t a) \
122 { \
123  return std::function(a); \
124 }
125 
126 namespace shogun
127 {
129  extern CRandom* sg_rand;
130  class CSGObject;
133 class CMath : public CSGObject
134 {
135  public:
139  CMath();
141 
143  virtual ~CMath();
145 
149 
151  //
152  template <class T>
153  static inline T min(T a, T b)
154  {
155  return (a<=b) ? a : b;
156  }
157 
159  template <class T>
160  static inline T max(T a, T b)
161  {
162  return (a>=b) ? a : b;
163  }
164 
166  template <class T>
167  static inline T clamp(T value, T lb, T ub)
168  {
169  if (value<=lb)
170  return lb;
171  else if (value>=ub)
172  return ub;
173  else
174  return value;
175  }
176 
178  template <class T>
179  static inline T abs(T a)
180  {
181  // can't be a>=0?(a):(-a), because compiler complains about
182  // 'comparison always true' when T is unsigned
183  if (a==0)
184  return 0;
185  else if (a>0)
186  return a;
187  else
188  return -a;
189  }
190 
192  static inline float64_t abs(complex128_t a)
193  {
194  float64_t a_real=a.real();
195  float64_t a_imag=a.imag();
196  return (CMath::sqrt(a_real*a_real+a_imag*a_imag));
197  }
199 
202 
203  static inline float64_t round(float64_t d)
204  {
205  return ::floor(d+0.5);
206  }
207 
208  static inline float64_t floor(float64_t d)
209  {
210  return ::floor(d);
211  }
212 
213  static inline float64_t ceil(float64_t d)
214  {
215  return ::ceil(d);
216  }
217 
219  template <class T>
220  static inline T sign(T a)
221  {
222  if (a==0)
223  return 0;
224  else return (a<0) ? (-1) : (+1);
225  }
226 
228  template <class T>
229  static inline void swap(T &a,T &b)
230  {
231  /* register for fast swaps */
232  register T c=a;
233  a=b;
234  b=c;
235  }
236 
238  template <class T>
239  static inline T sq(T x)
240  {
241  return x*x;
242  }
243 
245  static inline float32_t sqrt(float32_t x)
246  {
247  return ::sqrtf(x);
248  }
249 
251  static inline float64_t sqrt(float64_t x)
252  {
253  return ::sqrt(x);
254  }
255 
257  static inline floatmax_t sqrt(floatmax_t x)
258  {
259  //fall back to double precision sqrt if sqrtl is not
260  //available
261 #ifdef HAVE_SQRTL
262  return ::sqrtl(x);
263 #else
264  return ::sqrt(x);
265 #endif
266  }
267 
269  COMPLEX128_STDMATH(sqrt)
270 
271 
272  static inline float32_t invsqrt(float32_t x)
273  {
274  union float_to_int
275  {
276  float32_t f;
277  int32_t i;
278  };
279 
280  float_to_int tmp;
281  tmp.f=x;
282 
283  float32_t xhalf = 0.5f * x;
284  // store floating-point bits in integer tmp.i
285  // and do initial guess for Newton's method
286  tmp.i = 0x5f3759d5 - (tmp.i >> 1);
287  x = tmp.f; // convert new bits into float
288  x = x*(1.5f - xhalf*x*x); // One round of Newton's method
289  return x;
290  }
291 
293  static inline floatmax_t powl(floatmax_t x, floatmax_t n)
294  {
295  //fall back to double precision pow if powl is not
296  //available
297 #ifdef HAVE_POWL
298  return ::powl((long double) x, (long double) n);
299 #else
300  return ::pow((double) x, (double) n);
301 #endif
302  }
303 
304  static inline int32_t pow(bool x, int32_t n)
305  {
306  return 0;
307  }
308 
309  static inline int32_t pow(int32_t x, int32_t n)
310  {
311  ASSERT(n>=0)
312  int32_t result=1;
313  while (n--)
314  result*=x;
315 
316  return result;
317  }
318 
319  static inline float64_t pow(float64_t x, int32_t n)
320  {
321  if (n>=0)
322  {
323  float64_t result=1;
324  while (n--)
325  result*=x;
326 
327  return result;
328  }
329  else
330  return ::pow((double)x, (double)n);
331  }
332 
333  static inline float64_t pow(float64_t x, float64_t n)
334  {
335  return ::pow((double) x, (double) n);
336  }
337 
339  static inline complex128_t pow(complex128_t x, int32_t n)
340  {
341  return std::pow(x, n);
342  }
343 
345  {
346  return std::pow(x, n);
347  }
348 
350  {
351  return std::pow(x, n);
352  }
353 
355  {
356  return std::pow(x, n);
357  }
358 
359  static inline float64_t exp(float64_t x)
360  {
361  return ::exp((double) x);
362  }
363 
365  COMPLEX128_STDMATH(exp)
366 
367 
368  static inline float64_t tan(float64_t x)
369  {
370  return ::tan((double) x);
371  }
372 
374  COMPLEX128_STDMATH(tan)
375 
376 
377  static inline float64_t atan(float64_t x)
378  {
379  return ::atan((double) x);
380  }
381 
384 
385 
386  static inline float64_t atan2(float64_t x, float64_t y)
387  {
388  return ::atan2((double) x, (double) y);
389  }
390 
393 
394 
395  static inline float64_t tanh(float64_t x)
396  {
397  return ::tanh((double) x);
398  }
399 
401  COMPLEX128_STDMATH(tanh)
402 
403  static inline float64_t log10(float64_t v)
404  {
405  return ::log(v)/::log(10.0);
406  }
407 
409  COMPLEX128_STDMATH(log10)
410 
411  static inline float64_t log2(float64_t v)
412  {
413 #ifdef HAVE_LOG2
414  return ::log2(v);
415 #else
416  return ::log(v)/::log(2.0);
417 #endif //HAVE_LOG2
418  }
419 
420  static inline float64_t log(float64_t v)
421  {
422  return ::log(v);
423  }
424 
426  COMPLEX128_STDMATH(log)
427 
428  static inline index_t floor_log(index_t n)
429  {
430  index_t i;
431  for (i = 0; n != 0; i++)
432  n >>= 1;
433 
434  return i;
435  }
436 
437  static inline float64_t sin(float64_t x)
438  {
439  return ::sin(x);
440  }
441 
443  COMPLEX128_STDMATH(sin)
444 
445  static inline float64_t asin(float64_t x)
446  {
447  return ::asin(x);
448  }
449 
452 
453  static inline float64_t sinh(float64_t x)
454  {
455  return ::sinh(x);
456  }
457 
459  COMPLEX128_STDMATH(sinh)
460 
461  static inline float64_t cos(float64_t x)
462  {
463  return ::cos(x);
464  }
465 
467  COMPLEX128_STDMATH(cos)
468 
469  static inline float64_t acos(float64_t x)
470  {
471  return ::acos(x);
472  }
473 
476 
477  static inline float64_t cosh(float64_t x)
478  {
479  return ::cosh(x);
480  }
481 
483  COMPLEX128_STDMATH(cosh)
484 
485  static float64_t area_under_curve(float64_t* xy, int32_t len, bool reversed)
486  {
487  ASSERT(len>0 && xy)
488 
489  float64_t area = 0.0;
490 
491  if (!reversed)
492  {
493  for (int i=1; i<len; i++)
494  area += 0.5*(xy[2*i]-xy[2*(i-1)])*(xy[2*i+1]+xy[2*(i-1)+1]);
495  }
496  else
497  {
498  for (int i=1; i<len; i++)
499  area += 0.5*(xy[2*i+1]-xy[2*(i-1)+1])*(xy[2*i]+xy[2*(i-1)]);
500  }
501 
502  return area;
503  }
504 
505  static bool strtof(const char* str, float32_t* float_result);
506  static bool strtod(const char* str, float64_t* double_result);
507  static bool strtold(const char* str, floatmax_t* long_double_result);
508 
509  static inline int64_t factorial(int32_t n)
510  {
511  int64_t res=1;
512  for (int i=2; i<=n; i++)
513  res*=i ;
514  return res ;
515  }
516 
517  static void init_random(uint32_t initseed=0)
518  {
519  if (initseed==0)
520  seed = CRandom::generate_seed();
521  else
522  seed=initseed;
523 
524  sg_rand->set_seed(seed);
525  }
526 
527  static inline uint64_t random()
528  {
529  return sg_rand->random_64();
530  }
531 
532  static inline uint64_t random(uint64_t min_value, uint64_t max_value)
533  {
534  return sg_rand->random(min_value, max_value);
535  }
536 
537  static inline int64_t random(int64_t min_value, int64_t max_value)
538  {
539  return sg_rand->random(min_value, max_value);
540  }
541 
542  static inline uint32_t random(uint32_t min_value, uint32_t max_value)
543  {
544  return sg_rand->random(min_value, max_value);
545  }
546 
547  static inline int32_t random(int32_t min_value, int32_t max_value)
548  {
549  return sg_rand->random(min_value, max_value);
550  }
551 
552  static inline float32_t random(float32_t min_value, float32_t max_value)
553  {
554  return sg_rand->random(min_value, max_value);
555  }
556 
557  static inline float64_t random(float64_t min_value, float64_t max_value)
558  {
559  return sg_rand->random(min_value, max_value);
560  }
561 
562  static inline floatmax_t random(floatmax_t min_value, floatmax_t max_value)
563  {
564  return sg_rand->random(min_value, max_value);
565  }
566 
570  static inline float32_t normal_random(float32_t mean, float32_t std_dev)
571  {
572  // sets up variables & makes sure rand_s.range == (0,1)
573  float32_t ret;
574  float32_t rand_u;
575  float32_t rand_v;
576  float32_t rand_s;
577  do
578  {
579  rand_u = CMath::random(-1.0, 1.0);
580  rand_v = CMath::random(-1.0, 1.0);
581  rand_s = rand_u*rand_u + rand_v*rand_v;
582  } while ((rand_s == 0) || (rand_s >= 1));
583 
584  // the meat & potatos, and then the mean & standard deviation shifting...
585  ret = rand_u*sqrt(-2.0*log(rand_s)/rand_s);
586  ret = std_dev*ret + mean;
587  return ret;
588  }
589 
593  static inline float64_t normal_random(float64_t mean, float64_t std_dev)
594  {
595  return sg_rand->normal_distrib(mean, std_dev);
596  }
597 
600  static inline float32_t randn_float()
601  {
602  return normal_random(0.0, 1.0);
603  }
604 
607  static inline float64_t randn_double()
608  {
609  return sg_rand->std_normal_distrib();
610  }
611 
612  template <class T>
613  static int32_t get_num_nonzero(T* vec, int32_t len)
614  {
615  int32_t nnz = 0;
616  for (index_t i=0; i<len; ++i)
617  nnz += vec[i] != 0;
618 
619  return nnz;
620  }
621 
622  static int32_t get_num_nonzero(complex128_t* vec, int32_t len)
623  {
624  int32_t nnz=0;
625  for (index_t i=0; i<len; ++i)
626  nnz+=vec[i]!=0.0;
627  return nnz;
628  }
629 
630  static inline int64_t nchoosek(int32_t n, int32_t k)
631  {
632  int64_t res=1;
633 
634  for (int32_t i=n-k+1; i<=n; i++)
635  res*=i;
636 
637  return res/factorial(k);
638  }
639 
647  static void linspace(float64_t* output, float64_t start, float64_t end, int32_t n = 100);
648 
655  template <class T>
656  static T log_sum_exp(SGVector<T> values)
657  {
658  REQUIRE(values.vector, "Values are empty");
659 
660  /* find minimum element index */
661  index_t min_index=0;
662  T X0=values[0];
663  if (values.vlen>1)
664  {
665  for (index_t i=1; i<values.vlen; ++i)
666  {
667  if (values[i]<X0)
668  {
669  X0=values[i];
670  min_index=i;
671  }
672  }
673  }
674 
675  /* remove element from vector copy and compute log sum exp */
676  SGVector<T> values_without_X0(values.vlen-1);
677  index_t from_idx=0;
678  index_t to_idx=0;
679  for (from_idx=0; from_idx<values.vlen; ++from_idx)
680  {
681  if (from_idx!=min_index)
682  {
683  values_without_X0[to_idx]=exp(values[from_idx]-X0);
684  to_idx++;
685  }
686  }
687 
688  return X0+log(SGVector<T>::sum(values_without_X0)+1);
689  }
690 
697  template <class T>
698  static T log_mean_exp(SGVector<T> values)
699  {
700  return log_sum_exp(values) - log(values.vlen);
701  }
702 
706  static void sort(int32_t *a, int32_t cols, int32_t sort_col=0);
707  static void sort(float64_t *a, int32_t*idx, int32_t N);
708 
711  template <class T>
712  static void qsort(T* output, int32_t size)
713  {
714  if (size<=1)
715  return;
716 
717  if (size==2)
718  {
719  if (output[0] > output [1])
720  CMath::swap(output[0],output[1]);
721  return;
722  }
723  //T split=output[random(0,size-1)];
724  T split=output[size/2];
725 
726  int32_t left=0;
727  int32_t right=size-1;
728 
729  while (left<=right)
730  {
731  while (output[left] < split)
732  left++;
733  while (output[right] > split)
734  right--;
735 
736  if (left<=right)
737  {
738  CMath::swap(output[left],output[right]);
739  left++;
740  right--;
741  }
742  }
743 
744  if (right+1> 1)
745  qsort(output,right+1);
746 
747  if (size-left> 1)
748  qsort(&output[left],size-left);
749  }
750 
753  template <class T>
754  static void insertion_sort(T* output, int32_t size)
755  {
756  for (int32_t i=0; i<size-1; i++)
757  {
758  int32_t j=i-1;
759  T value=output[i];
760  while (j >= 0 && output[j] > value)
761  {
762  output[j+1] = output[j];
763  j--;
764  }
765  output[j+1]=value;
766  }
767  }
768 
770  template <class T>
771  inline static void radix_sort(T* array, int32_t size)
772  {
773  radix_sort_helper(array,size,0);
774  }
775 
776  /*
777  * Inline function to extract the byte at position p (from left)
778  * of an 64 bit integer. The function is somewhat identical to
779  * accessing an array of characters via [].
780  */
781  template <class T>
782  static inline uint8_t byte(T word, uint16_t p)
783  {
784  return (word >> (sizeof(T)-p-1) * 8) & 0xff;
785  }
786 
788  static inline uint8_t byte(complex128_t word, uint16_t p)
789  {
790  SG_SERROR("CMath::byte():: Not supported for complex128_t\n");
791  return uint8_t(0);
792  }
793 
794  template <class T>
795  static void radix_sort_helper(T* array, int32_t size, uint16_t i)
796  {
797  static size_t count[256], nc, cmin;
798  T *ak;
799  uint8_t c=0;
800  radix_stack_t<T> s[RADIX_STACK_SIZE], *sp, *olds, *bigs;
801  T *an, *aj, *pile[256];
802  size_t *cp, cmax;
803 
804  /* Push initial array to stack */
805  sp = s;
806  radix_push(array, size, i);
807 
808  /* Loop until all digits have been sorted */
809  while (sp>s) {
810  radix_pop(array, size, i);
811  an = array + size;
812 
813  /* Make character histogram */
814  if (nc == 0) {
815  cmin = 0xff;
816  for (ak = array; ak < an; ak++) {
817  c = byte(*ak, i);
818  count[c]++;
819  if (count[c] == 1) {
820  /* Determine smallest character */
821  if (c < cmin)
822  cmin = c;
823  nc++;
824  }
825  }
826 
827  /* Sort recursively if stack size too small */
828  if (sp + nc > s + RADIX_STACK_SIZE) {
829  radix_sort_helper(array, size, i);
830  continue;
831  }
832  }
833 
834  /*
835  * Set pile[]; push incompletely sorted bins onto stack.
836  * pile[] = pointers to last out-of-place element in bins.
837  * Before permuting: pile[c-1] + count[c] = pile[c];
838  * during deal: pile[c] counts down to pile[c-1].
839  */
840  olds = bigs = sp;
841  cmax = 2;
842  ak = array;
843  pile[0xff] = an;
844  for (cp = count + cmin; nc > 0; cp++) {
845  /* Find next non-empty pile */
846  while (*cp == 0)
847  cp++;
848  /* Pile with several entries */
849  if (*cp > 1) {
850  /* Determine biggest pile */
851  if (*cp > cmax) {
852  cmax = *cp;
853  bigs = sp;
854  }
855 
856  if (i < sizeof(T)-1)
857  radix_push(ak, *cp, (uint16_t) (i + 1));
858  }
859  pile[cp - count] = ak += *cp;
860  nc--;
861  }
862 
863  /* Play it safe -- biggest bin last. */
864  swap(*olds, *bigs);
865 
866  /*
867  * Permute misplacements home. Already home: everything
868  * before aj, and in pile[c], items from pile[c] on.
869  * Inner loop:
870  * r = next element to put in place;
871  * ak = pile[r[i]] = location to put the next element.
872  * aj = bottom of 1st disordered bin.
873  * Outer loop:
874  * Once the 1st disordered bin is done, ie. aj >= ak,
875  * aj<-aj + count[c] connects the bins in array linked list;
876  * reset count[c].
877  */
878  aj = array;
879  while (aj<an)
880  {
881  T r;
882 
883  for (r = *aj; aj < (ak = --pile[c = byte(r, i)]);)
884  swap(*ak, r);
885 
886  *aj = r;
887  aj += count[c];
888  count[c] = 0;
889  }
890  }
891  }
892 
894  static void radix_sort_helper(complex128_t* array, int32_t size, uint16_t i)
895  {
896  SG_SERROR("CMath::radix_sort_helper():: Not supported for complex128_t\n");
897  }
898 
908  template <class T>
909  static void qsort(T** vector, index_t length)
910  {
911  if (length<=1)
912  return;
913 
914  if (length==2)
915  {
916  if (*vector[0]>*vector[1])
917  swap(vector[0],vector[1]);
918  return;
919  }
920  T* split=vector[length/2];
921 
922  int32_t left=0;
923  int32_t right=length-1;
924 
925  while (left<=right)
926  {
927  while (*vector[left]<*split)
928  ++left;
929  while (*vector[right]>*split)
930  --right;
931 
932  if (left<=right)
933  {
934  swap(vector[left],vector[right]);
935  ++left;
936  --right;
937  }
938  }
939 
940  if (right+1>1)
941  qsort(vector,right+1);
942 
943  if (length-left>1)
944  qsort(&vector[left],length-left);
945  }
946 
948  static void qsort(complex128_t** vector, index_t length)
949  {
950  SG_SERROR("CMath::qsort():: Not supported for complex128_t\n");
951  }
952 
954  template <class T> static void display_bits(T word, int32_t width=8*sizeof(T))
955  {
956  ASSERT(width>=0)
957  for (int i=0; i<width; i++)
958  {
959  T mask = ((T) 1)<<(sizeof(T)*8-1);
960  while (mask)
961  {
962  if (mask & word)
963  SG_SPRINT("1")
964  else
965  SG_SPRINT("0")
966 
967  mask>>=1;
968  }
969  }
970  }
971 
973  static void display_bits(complex128_t word,
974  int32_t width=8*sizeof(complex128_t))
975  {
976  SG_SERROR("CMath::display_bits():: Not supported for complex128_t\n");
977  }
978 
984  template <class T1,class T2>
985  static void qsort_index(T1* output, T2* index, uint32_t size);
986 
988  template <class T>
989  static void qsort_index(complex128_t* output, T* index, uint32_t size)
990  {
991  SG_SERROR("CMath::qsort_index():: Not supported for complex128_t\n");
992  }
993 
999  template <class T1,class T2>
1000  static void qsort_backward_index(
1001  T1* output, T2* index, int32_t size);
1002 
1004  template <class T>
1006  complex128_t* output, T* index, uint32_t size)
1007  {
1008  SG_SERROR("CMath::qsort_backword_index():: \
1009  Not supported for complex128_t\n");
1010  }
1011 
1019  template <class T1,class T2>
1020  inline static void parallel_qsort_index(T1* output, T2* index, uint32_t size, int32_t n_threads, int32_t limit=262144)
1021  {
1022  int32_t n=0;
1023  thread_qsort<T1,T2> t;
1024  t.output=output;
1025  t.index=index;
1026  t.size=size;
1027  t.qsort_threads=&n;
1028  t.sort_limit=limit;
1029  t.num_threads=n_threads;
1030  parallel_qsort_index<T1,T2>(&t);
1031  }
1032 
1034  template <class T>
1035  inline static void parallel_qsort_index(complex128_t* output, T* index,
1036  uint32_t size, int32_t n_threads, int32_t limit=0)
1037  {
1038  SG_SERROR("CMath::parallel_qsort_index():: Not supported for complex128_t\n");
1039  }
1040 
1041 
1042  template <class T1,class T2>
1043  static void* parallel_qsort_index(void* p);
1044 
1045 
1046  /* finds the smallest element in output and puts that element as the
1047  first element */
1048  template <class T>
1049  static void min(float64_t* output, T* index, int32_t size);
1050 
1052  static void min(float64_t* output, complex128_t* index, int32_t size)
1053  {
1054  SG_SERROR("CMath::min():: Not supported for complex128_t\n");
1055  }
1056 
1057  /* finds the n smallest elements in output and puts these elements as the
1058  first n elements */
1059  template <class T>
1060  static void nmin(
1061  float64_t* output, T* index, int32_t size, int32_t n);
1062 
1064  static void nmin(float64_t* output, complex128_t* index,
1065  int32_t size, int32_t n)
1066  {
1067  SG_SERROR("CMath::nmin():: Not supported for complex128_t\n");
1068  }
1069 
1070 
1071 
1072  /* finds an element in a sorted array via binary search
1073  * returns -1 if not found */
1074  template <class T>
1075  static int32_t binary_search_helper(T* output, int32_t size, T elem)
1076  {
1077  int32_t start=0;
1078  int32_t end=size-1;
1079 
1080  if (size<1)
1081  return 0;
1082 
1083  while (start<end)
1084  {
1085  int32_t middle=(start+end)/2;
1086 
1087  if (output[middle]>elem)
1088  end=middle-1;
1089  else if (output[middle]<elem)
1090  start=middle+1;
1091  else
1092  return middle;
1093  }
1094 
1095  return start;
1096  }
1097 
1099  static int32_t binary_search_helper(complex128_t* output, int32_t size, complex128_t elem)
1100  {
1101  SG_SERROR("CMath::binary_search_helper():: Not supported for complex128_t\n");
1102  return int32_t(0);
1103  }
1104 
1105  /* finds an element in a sorted array via binary search
1106  * * returns -1 if not found */
1107  template <class T>
1108  static inline int32_t binary_search(T* output, int32_t size, T elem)
1109  {
1110  int32_t ind = binary_search_helper(output, size, elem);
1111  if (ind >= 0 && output[ind] == elem)
1112  return ind;
1113  return -1;
1114  }
1115 
1117  static inline int32_t binary_search(complex128_t* output, int32_t size, complex128_t elem)
1118  {
1119  SG_SERROR("CMath::binary_search():: Not supported for complex128_t\n");
1120  return int32_t(-1);
1121  }
1122 
1123 
1124  /* Finds an element in a sorted array of pointers via binary search
1125  * Every element is dereferenced once before being compared
1126  *
1127  * @param array array of pointers to search in (assumed being sorted)
1128  * @param length length of array
1129  * @param elem pointer to element to search for
1130  * @return index of elem, -1 if not found */
1131  template<class T>
1132  static inline int32_t binary_search(T** vector, index_t length,
1133  T* elem)
1134  {
1135  int32_t start=0;
1136  int32_t end=length-1;
1137 
1138  if (length<1)
1139  return -1;
1140 
1141  while (start<end)
1142  {
1143  int32_t middle=(start+end)/2;
1144 
1145  if (*vector[middle]>*elem)
1146  end=middle-1;
1147  else if (*vector[middle]<*elem)
1148  start=middle+1;
1149  else
1150  {
1151  start=middle;
1152  break;
1153  }
1154  }
1155 
1156  if (start>=0&&*vector[start]==*elem)
1157  return start;
1158 
1159  return -1;
1160  }
1161 
1163  static inline int32_t binary_search(complex128_t** vector, index_t length, complex128_t* elem)
1164  {
1165  SG_SERROR("CMath::binary_search():: Not supported for complex128_t\n");
1166  return int32_t(-1);
1167  }
1168 
1169  template <class T>
1171  T* output, int32_t size, T elem)
1172  {
1173  int32_t ind = binary_search_helper(output, size, elem);
1174 
1175  if (output[ind]<=elem)
1176  return ind;
1177  if (ind>0 && output[ind-1] <= elem)
1178  return ind-1;
1179  return -1;
1180  }
1181 
1184  int32_t size, complex128_t elem)
1185  {
1186  SG_SERROR("CMath::binary_search_max_lower_equal():: \
1187  Not supported for complex128_t\n");
1188  return int32_t(-1);
1189  }
1190 
1193  static float64_t Align(
1194  char * seq1, char* seq2, int32_t l1, int32_t l2, float64_t gapCost);
1195 
1196 
1198 
1201  {
1202  return c.real();
1203  }
1204 
1207  {
1208  return c.imag();
1209  }
1210 
1212  inline static uint32_t get_seed()
1213  {
1214  return CMath::seed;
1215  }
1216 
1218  inline static uint32_t get_log_range()
1219  {
1220  return CMath::LOGRANGE;
1221  }
1222 
1223 #ifdef USE_LOGCACHE
1224  inline static uint32_t get_log_accuracy()
1226  {
1227  return CMath::LOGACCURACY;
1228  }
1229 #endif
1230 
1232  static int is_finite(double f);
1233 
1235  static int is_infinity(double f);
1236 
1238  static int is_nan(double f);
1239 
1250 #ifdef USE_LOGCACHE
1251  static inline float64_t logarithmic_sum(float64_t p, float64_t q)
1252  {
1253  float64_t diff;
1254 
1255  if (!CMath::is_finite(p))
1256  return q;
1257 
1258  if (!CMath::is_finite(q))
1259  {
1260  SG_SWARNING("INVALID second operand to logsum(%f,%f) expect undefined results\n", p, q)
1261  return NOT_A_NUMBER;
1262  }
1263  diff = p - q;
1264  if (diff > 0)
1265  return diff > LOGRANGE? p : p + logtable[(int)(diff * LOGACCURACY)];
1266  return -diff > LOGRANGE? q : q + logtable[(int)(-diff * LOGACCURACY)];
1267  }
1268 
1270  static void init_log_table();
1271 
1273  static int32_t determine_logrange();
1274 
1276  static int32_t determine_logaccuracy(int32_t range);
1277 #else
1279  float64_t p, float64_t q)
1280  {
1281  float64_t diff;
1282 
1283  if (!CMath::is_finite(p))
1284  return q;
1285  if (!CMath::is_finite(q))
1286  return p;
1287  diff = p - q;
1288  if (diff > 0)
1289  return diff > LOGRANGE? p : p + log(1 + exp(-diff));
1290  return -diff > LOGRANGE? q : q + log(1 + exp(diff));
1291  }
1292 #endif
1293 #ifdef USE_LOGSUMARRAY
1294 
1299  static inline float64_t logarithmic_sum_array(
1300  float64_t *p, int32_t len)
1301  {
1302  if (len<=2)
1303  {
1304  if (len==2)
1305  return logarithmic_sum(p[0],p[1]) ;
1306  if (len==1)
1307  return p[0];
1308  return -INFTY ;
1309  }
1310  else
1311  {
1312  register float64_t *pp=p ;
1313  if (len%2==1) pp++ ;
1314  for (register int32_t j=0; j < len>>1; j++)
1315  pp[j]=logarithmic_sum(pp[j<<1], pp[1+(j<<1)]) ;
1316  }
1317  return logarithmic_sum_array(p,len%2 + (len>>1)) ;
1318  }
1319 #endif //USE_LOGSUMARRAY
1320 
1321 
1323  virtual const char* get_name() const { return "Math"; }
1324  public:
1327  static const float64_t NOT_A_NUMBER;
1330  static const float64_t INFTY;
1331  static const float64_t ALMOST_INFTY;
1332 
1335 
1337  static const float64_t PI;
1338 
1341 
1342  /* largest and smallest possible float64_t */
1345 
1346  protected:
1348  static int32_t LOGRANGE;
1349 
1351  static uint32_t seed;
1352 
1353 #ifdef USE_LOGCACHE
1354 
1356  static int32_t LOGACCURACY;
1358  static float64_t* logtable;
1360 #endif
1361 };
1362 
1363 //implementations of template functions
1364 template <class T1,class T2>
1365 void* CMath::parallel_qsort_index(void* p)
1366  {
1367  struct thread_qsort<T1,T2>* ps=(thread_qsort<T1,T2>*) p;
1368  T1* output=ps->output;
1369  T2* index=ps->index;
1370  int32_t size=ps->size;
1371  int32_t* qsort_threads=ps->qsort_threads;
1372  int32_t sort_limit=ps->sort_limit;
1373  int32_t num_threads=ps->num_threads;
1374 
1375  if (size<2)
1376  {
1377  return NULL;
1378  }
1379 
1380  if (size==2)
1381  {
1382  if (output[0] > output [1])
1383  {
1384  swap(output[0], output[1]);
1385  swap(index[0], index[1]);
1386  }
1387  return NULL;
1388  }
1389  //T1 split=output[random(0,size-1)];
1390  T1 split=output[size/2];
1391 
1392  int32_t left=0;
1393  int32_t right=size-1;
1394 
1395  while (left<=right)
1396  {
1397  while (output[left] < split)
1398  left++;
1399  while (output[right] > split)
1400  right--;
1401 
1402  if (left<=right)
1403  {
1404  swap(output[left], output[right]);
1405  swap(index[left], index[right]);
1406  left++;
1407  right--;
1408  }
1409  }
1410  bool lthread_start=false;
1411  bool rthread_start=false;
1412  pthread_t lthread;
1413  pthread_t rthread;
1414  struct thread_qsort<T1,T2> t1;
1415  struct thread_qsort<T1,T2> t2;
1416 
1417  if (right+1> 1 && (right+1< sort_limit || *qsort_threads >= num_threads-1))
1418  qsort_index(output,index,right+1);
1419  else if (right+1> 1)
1420  {
1421  (*qsort_threads)++;
1422  lthread_start=true;
1423  t1.output=output;
1424  t1.index=index;
1425  t1.size=right+1;
1426  t1.qsort_threads=qsort_threads;
1427  t1.sort_limit=sort_limit;
1428  t1.num_threads=num_threads;
1429  if (pthread_create(&lthread, NULL, parallel_qsort_index<T1,T2>, &t1) != 0)
1430  {
1431  lthread_start=false;
1432  (*qsort_threads)--;
1433  qsort_index(output,index,right+1);
1434  }
1435  }
1436 
1437 
1438  if (size-left> 1 && (size-left< sort_limit || *qsort_threads >= num_threads-1))
1439  qsort_index(&output[left],&index[left], size-left);
1440  else if (size-left> 1)
1441  {
1442  (*qsort_threads)++;
1443  rthread_start=true;
1444  t2.output=&output[left];
1445  t2.index=&index[left];
1446  t2.size=size-left;
1447  t2.qsort_threads=qsort_threads;
1448  t2.sort_limit=sort_limit;
1449  t2.num_threads=num_threads;
1450  if (pthread_create(&rthread, NULL, parallel_qsort_index<T1,T2>, &t2) != 0)
1451  {
1452  rthread_start=false;
1453  (*qsort_threads)--;
1454  qsort_index(&output[left],&index[left], size-left);
1455  }
1456  }
1457 
1458  if (lthread_start)
1459  {
1460  pthread_join(lthread, NULL);
1461  (*qsort_threads)--;
1462  }
1463 
1464  if (rthread_start)
1465  {
1466  pthread_join(rthread, NULL);
1467  (*qsort_threads)--;
1468  }
1469 
1470  return NULL;
1471  }
1472 
1473  template <class T1,class T2>
1474 void CMath::qsort_index(T1* output, T2* index, uint32_t size)
1475 {
1476  if (size<=1)
1477  return;
1478 
1479  if (size==2)
1480  {
1481  if (output[0] > output [1])
1482  {
1483  swap(output[0],output[1]);
1484  swap(index[0],index[1]);
1485  }
1486  return;
1487  }
1488  //T1 split=output[random(0,size-1)];
1489  T1 split=output[size/2];
1490 
1491  int32_t left=0;
1492  int32_t right=size-1;
1493 
1494  while (left<=right)
1495  {
1496  while (output[left] < split)
1497  left++;
1498  while (output[right] > split)
1499  right--;
1500 
1501  if (left<=right)
1502  {
1503  swap(output[left],output[right]);
1504  swap(index[left],index[right]);
1505  left++;
1506  right--;
1507  }
1508  }
1509 
1510  if (right+1> 1)
1511  qsort_index(output,index,right+1);
1512 
1513  if (size-left> 1)
1514  qsort_index(&output[left],&index[left], size-left);
1515 }
1516 
1517  template <class T1,class T2>
1518 void CMath::qsort_backward_index(T1* output, T2* index, int32_t size)
1519 {
1520  if (size<=1)
1521  return;
1522 
1523  if (size==2)
1524  {
1525  if (output[0] < output [1])
1526  {
1527  swap(output[0],output[1]);
1528  swap(index[0],index[1]);
1529  }
1530  return;
1531  }
1532 
1533  //T1 split=output[random(0,size-1)];
1534  T1 split=output[size/2];
1535 
1536  int32_t left=0;
1537  int32_t right=size-1;
1538 
1539  while (left<=right)
1540  {
1541  while (output[left] > split)
1542  left++;
1543  while (output[right] < split)
1544  right--;
1545 
1546  if (left<=right)
1547  {
1548  swap(output[left],output[right]);
1549  swap(index[left],index[right]);
1550  left++;
1551  right--;
1552  }
1553  }
1554 
1555  if (right+1> 1)
1556  qsort_backward_index(output,index,right+1);
1557 
1558  if (size-left> 1)
1559  qsort_backward_index(&output[left],&index[left], size-left);
1560 }
1561 
1562  template <class T>
1563 void CMath::nmin(float64_t* output, T* index, int32_t size, int32_t n)
1564 {
1565  if (6*n*size<13*size*CMath::log(size))
1566  for (int32_t i=0; i<n; i++)
1567  min(&output[i], &index[i], size-i) ;
1568  else
1569  qsort_index(output, index, size) ;
1570 }
1571 
1572 /* move the smallest entry in the array to the beginning */
1573  template <class T>
1574 void CMath::min(float64_t* output, T* index, int32_t size)
1575 {
1576  if (size<=1)
1577  return;
1578  float64_t min_elem=output[0];
1579  int32_t min_index=0;
1580  for (int32_t i=1; i<size; i++)
1581  {
1582  if (output[i]<min_elem)
1583  {
1584  min_index=i;
1585  min_elem=output[i];
1586  }
1587  }
1588  swap(output[0], output[min_index]);
1589  swap(index[0], index[min_index]);
1590 }
1591 
1592 #define COMPLEX128_ERROR_ONEARG_T(function) \
1593 template <> \
1594 inline complex128_t CMath::function<complex128_t>(complex128_t a) \
1595 { \
1596  SG_SERROR("CMath::%s():: Not supported for complex128_t\n",\
1597  #function);\
1598  return complex128_t(0.0, 0.0); \
1599 }
1600 
1601 #define COMPLEX128_ERROR_TWOARGS_T(function) \
1602 template <> \
1603 inline complex128_t CMath::function<complex128_t>(complex128_t a, complex128_t b) \
1604 { \
1605  SG_SERROR("CMath::%s():: Not supported for complex128_t\n",\
1606  #function);\
1607  return complex128_t(0.0, 0.0); \
1608 }
1609 
1610 #define COMPLEX128_ERROR_THREEARGS_T(function) \
1611 template <> \
1612 inline complex128_t CMath::function<complex128_t>(complex128_t a, complex128_t b, complex128_t c) \
1613 { \
1614  SG_SERROR("CMath::%s():: Not supported for complex128_t\n",\
1615  #function);\
1616  return complex128_t(0.0, 0.0); \
1617 }
1618 
1619 #define COMPLEX128_ERROR_SORT_T(function) \
1620 template <> \
1621 inline void CMath::function<complex128_t>(complex128_t* output, int32_t b) \
1622 { \
1623  SG_SERROR("CMath::%s():: Not supported for complex128_t\n",\
1624  #function);\
1625 }
1626 
1629 
1632 
1635 
1637 // COMPLEX128_ERROR_ONEARG_T(sign)
1638 
1641 
1643 COMPLEX128_ERROR_SORT_T(insertion_sort)
1644 
1646 COMPLEX128_ERROR_SORT_T(radix_sort)
1647 
1648 }
1649 #undef COMPLEX128_ERROR_ONEARG
1650 #undef COMPLEX128_ERROR_ONEARG_T
1651 #undef COMPLEX128_ERROR_TWOARGS_T
1652 #undef COMPLEX128_ERROR_THREEARGS_T
1653 #undef COMPLEX128_STDMATH
1654 #undef COMPLEX128_ERROR_SORT_T
1655 #endif
static float64_t sin(float64_t x)
Definition: Math.h:437
static float64_t normal_random(float64_t mean, float64_t std_dev)
Definition: Math.h:593
static const float64_t MACHINE_EPSILON
Definition: Math.h:1340
static int32_t binary_search(complex128_t *output, int32_t size, complex128_t elem)
binary_search not implemented for complex128_t
Definition: Math.h:1117
static void parallel_qsort_index(T1 *output, T2 *index, uint32_t size, int32_t n_threads, int32_t limit=262144)
Definition: Math.h:1020
static uint32_t seed
random generator seed
Definition: Math.h:1351
std::complex< float64_t > complex128_t
Definition: common.h:65
uint64_t random(uint64_t min_value, uint64_t max_value)
Definition: Random.h:100
float64_t std_normal_distrib() const
Definition: Random.cpp:238
static float64_t sqrt(float64_t x)
x^0.5
Definition: Math.h:251
static floatmax_t random(floatmax_t min_value, floatmax_t max_value)
Definition: Math.h:562
static floatmax_t sqrt(floatmax_t x)
x^0.5
Definition: Math.h:257
int32_t index_t
Definition: common.h:60
static float64_t ceil(float64_t d)
Definition: Math.h:213
static complex128_t pow(complex128_t x, int32_t n)
x^n, x or n being a complex128_t
Definition: Math.h:339
static const float64_t INFTY
infinity
Definition: Math.h:1330
static int32_t binary_search_helper(T *output, int32_t size, T elem)
Definition: Math.h:1075
#define COMPLEX128_ERROR_THREEARGS_T(function)
Definition: Math.h:1610
static void qsort_index(complex128_t *output, T *index, uint32_t size)
qsort_index not implemented for complex128_t
Definition: Math.h:989
#define COMPLEX128_ERROR_TWOARGS_T(function)
Definition: Math.h:1601
#define SG_SWARNING(...)
Definition: SGIO.h:180
static float32_t normal_random(float32_t mean, float32_t std_dev)
Definition: Math.h:570
uint64_t random_64() const
Definition: Random.cpp:135
static float64_t random(float64_t min_value, float64_t max_value)
Definition: Math.h:557
static int32_t binary_search(T *output, int32_t size, T elem)
Definition: Math.h:1108
static T sq(T x)
x^2
Definition: Math.h:239
static uint32_t get_seed()
returns number generator seed
Definition: Math.h:1212
static float32_t randn_float()
Definition: Math.h:600
static const float64_t MIN_REAL_NUMBER
Definition: Math.h:1344
static float64_t randn_double()
Definition: Math.h:607
static int32_t binary_search_max_lower_equal(T *output, int32_t size, T elem)
Definition: Math.h:1170
#define REQUIRE(x,...)
Definition: SGIO.h:208
static uint32_t get_log_range()
returns range of logtable
Definition: Math.h:1218
void split(v_array< ds_node< P > > &point_set, v_array< ds_node< P > > &far_set, int max_scale)
Definition: JLCoverTree.h:147
static float32_t random(float32_t min_value, float32_t max_value)
Definition: Math.h:552
CRandom * sg_rand
Definition: init.cpp:29
#define RADIX_STACK_SIZE
Definition: Math.h:73
static float64_t imag(complex128_t c)
returns imag part of a complex128_t number
Definition: Math.h:1206
static float64_t floor(float64_t d)
Definition: Math.h:208
static uint64_t random()
Definition: Math.h:527
static int32_t get_num_nonzero(complex128_t *vec, int32_t len)
Definition: Math.h:622
static float64_t real(complex128_t c)
returns real part of a complex128_t number
Definition: Math.h:1200
static uint8_t byte(complex128_t word, uint16_t p)
byte not implemented for complex128_t
Definition: Math.h:788
static void qsort(T *output, int32_t size)
Definition: Math.h:712
static int32_t LOGRANGE
range for logtable: log(1+exp(x)) -LOGRANGE <= x <= 0
Definition: Math.h:1348
static const float64_t ALMOST_NEG_INFTY
almost neg (log) infinity
Definition: Math.h:1334
static complex128_t pow(complex128_t x, complex128_t n)
Definition: Math.h:344
static uint8_t byte(T word, uint16_t p)
Definition: Math.h:782
float64_t normal_distrib(float64_t mu, float64_t sigma) const
Definition: Random.cpp:233
static complex128_t pow(float64_t x, complex128_t n)
Definition: Math.h:354
#define SG_SPRINT(...)
Definition: SGIO.h:182
static float64_t pow(float64_t x, int32_t n)
Definition: Math.h:319
#define ASSERT(x)
Definition: SGIO.h:203
Class SGObject is the base class of all shogun objects.
Definition: SGObject.h:102
static int32_t random(int32_t min_value, int32_t max_value)
Definition: Math.h:547
static float64_t pow(float64_t x, float64_t n)
Definition: Math.h:333
static void qsort(T **vector, index_t length)
Definition: Math.h:909
static void init_random(uint32_t initseed=0)
Definition: Math.h:517
static int32_t pow(int32_t x, int32_t n)
Definition: Math.h:309
#define radix_pop(a, n, i)
Definition: Math.h:77
double float64_t
Definition: common.h:48
static void parallel_qsort_index(complex128_t *output, T *index, uint32_t size, int32_t n_threads, int32_t limit=0)
parallel_qsort_index not implemented for complex128_t
Definition: Math.h:1035
static void min(float64_t *output, complex128_t *index, int32_t size)
complex128_t cannot be used as index
Definition: Math.h:1052
long double floatmax_t
Definition: common.h:49
#define COMPLEX128_ERROR_ONEARG(function)
Definition: Math.h:112
static void qsort_backword_index(complex128_t *output, T *index, uint32_t size)
qsort_backword_index not implemented for complex128_t
Definition: Math.h:1005
static int32_t get_num_nonzero(T *vec, int32_t len)
Definition: Math.h:613
static T log_sum_exp(SGVector< T > values)
Definition: Math.h:656
static T max(T a, T b)
return the maximum of two integers
Definition: Math.h:160
static void display_bits(T word, int32_t width=8 *sizeof(T))
display bits (useful for debugging)
Definition: Math.h:954
static uint64_t random(uint64_t min_value, uint64_t max_value)
Definition: Math.h:532
static int64_t factorial(int32_t n)
Definition: Math.h:509
static void qsort(complex128_t **vector, index_t length)
qsort not implemented for complex128_t
Definition: Math.h:948
static int32_t binary_search_max_lower_equal(complex128_t *output, int32_t size, complex128_t elem)
binary_search_max_lower_equal not implemented for complex128_t
Definition: Math.h:1183
float float32_t
Definition: common.h:47
static int64_t nchoosek(int32_t n, int32_t k)
Definition: Math.h:630
all of classes and functions are contained in the shogun namespace
Definition: class_list.h:16
static int32_t binary_search(T **vector, index_t length, T *elem)
Definition: Math.h:1132
void set_seed(uint32_t seed)
Definition: Random.cpp:57
static float64_t abs(complex128_t a)
return the absolute value of a complex number
Definition: Math.h:192
static void display_bits(complex128_t word, int32_t width=8 *sizeof(complex128_t))
disply_bits not implemented for complex128_t
Definition: Math.h:973
#define COMPLEX128_STDMATH(function)
Definition: Math.h:120
static void radix_sort_helper(T *array, int32_t size, uint16_t i)
Definition: Math.h:795
static T sign(T a)
signum of type T variable a
Definition: Math.h:220
virtual const char * get_name() const
Definition: Math.h:1323
static T min(T a, T b)
return the minimum of two integers
Definition: Math.h:153
#define SG_SERROR(...)
Definition: SGIO.h:181
static float64_t exp(float64_t x)
Definition: Math.h:359
static float64_t log(float64_t v)
Definition: Math.h:420
static void radix_sort_helper(complex128_t *array, int32_t size, uint16_t i)
radix_sort_helper not implemented for complex128_t
Definition: Math.h:894
static const float64_t ALMOST_INFTY
Definition: Math.h:1331
Class which collects generic mathematical functions.
Definition: Math.h:133
static T log_mean_exp(SGVector< T > values)
Definition: Math.h:698
static void insertion_sort(T *output, int32_t size)
Definition: Math.h:754
static void swap(T &a, T &b)
swap e.g. floats a and b
Definition: Math.h:229
static complex128_t pow(complex128_t x, float64_t n)
Definition: Math.h:349
static int32_t binary_search_helper(complex128_t *output, int32_t size, complex128_t elem)
binary_search_helper not implemented for complex128_t
Definition: Math.h:1099
static float64_t round(float64_t d)
Definition: Math.h:203
static float32_t sqrt(float32_t x)
x^0.5
Definition: Math.h:245
static void radix_sort(T *array, int32_t size)
Definition: Math.h:771
static uint32_t random(uint32_t min_value, uint32_t max_value)
Definition: Math.h:542
static floatmax_t powl(floatmax_t x, floatmax_t n)
x^n
Definition: Math.h:293
static float64_t logarithmic_sum(float64_t p, float64_t q)
Definition: Math.h:1278
#define radix_push(a, n, i)
Definition: Math.h:76
static T clamp(T value, T lb, T ub)
return the value clamped to interval [lb,ub]
Definition: Math.h:167
static int32_t binary_search(complex128_t **vector, index_t length, complex128_t *elem)
binary_search not implemented for complex128_t
Definition: Math.h:1163
#define COMPLEX128_ERROR_SORT_T(function)
Definition: Math.h:1619
static int32_t pow(bool x, int32_t n)
Definition: Math.h:304
static const float64_t MAX_REAL_NUMBER
Definition: Math.h:1343
index_t vlen
Definition: SGVector.h:706
static T abs(T a)
return the absolute value of a number
Definition: Math.h:179
static void nmin(float64_t *output, complex128_t *index, int32_t size, int32_t n)
complex128_t cannot be used as index
Definition: Math.h:1064
static int64_t random(int64_t min_value, int64_t max_value)
Definition: Math.h:537
static const float64_t PI
Definition: Math.h:1337

SHOGUN Machine Learning Toolbox - Documentation