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

SHOGUN Machine Learning Toolbox - Documentation