Free cookie consent management tool by TermsFeed Policy Generator

Changeset 8803


Ignore:
Timestamp:
10/12/12 15:34:05 (11 years ago)
Author:
mkommend
Message:

#1968: Added [ThreadStatic] to the RNG of ALGLIB and removed lock from random forest algorithm.

Location:
trunk/sources
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/sources/HeuristicLab.Algorithms.DataAnalysis/3.4/RandomForest/RandomForestClassification.cs

    r8786 r8803  
    134134      if (m <= 0 || m > 1) throw new ArgumentException("The M parameter in the random forest regression must be between 0 and 1.");
    135135
    136       lock (alglib.math.rndobject) {
    137         alglib.math.rndobject = new System.Random(seed);
    138       }
     136      alglib.math.rndobject = new System.Random(seed);
    139137
    140138      Dataset dataset = problemData.Dataset;
  • trunk/sources/HeuristicLab.Algorithms.DataAnalysis/3.4/RandomForest/RandomForestRegression.cs

    r8786 r8803  
    133133      if (m <= 0 || m > 1) throw new ArgumentException("The M parameter in the random forest regression must be between 0 and 1.");
    134134
    135       lock (alglib.math.rndobject) {
    136         alglib.math.rndobject = new System.Random(seed);
    137       }
     135      alglib.math.rndobject = new System.Random(seed);
    138136
    139137      Dataset dataset = problemData.Dataset;
  • trunk/sources/HeuristicLab.ExtLibs/HeuristicLab.ALGLIB/3.6.0/ALGLIB-3.6.0/ap.cs

    r8418 r8803  
    1919*************************************************************************/
    2020using System;
    21 public partial class alglib
    22 {
    23     /********************************************************************
    24     Callback definitions for optimizers/fitters/solvers.
     21public partial class alglib {
     22  /********************************************************************
     23  Callback definitions for optimizers/fitters/solvers.
    2524   
    26     Callbacks for unparameterized (general) functions:
    27     * ndimensional_func         calculates f(arg), stores result to func
    28     * ndimensional_grad         calculates func = f(arg),
    29                                 grad[i] = df(arg)/d(arg[i])
    30     * ndimensional_hess         calculates func = f(arg),
    31                                 grad[i] = df(arg)/d(arg[i]),
    32                                 hess[i,j] = d2f(arg)/(d(arg[i])*d(arg[j]))
     25  Callbacks for unparameterized (general) functions:
     26  * ndimensional_func         calculates f(arg), stores result to func
     27  * ndimensional_grad         calculates func = f(arg),
     28                              grad[i] = df(arg)/d(arg[i])
     29  * ndimensional_hess         calculates func = f(arg),
     30                              grad[i] = df(arg)/d(arg[i]),
     31                              hess[i,j] = d2f(arg)/(d(arg[i])*d(arg[j]))
    3332   
    34     Callbacks for systems of functions:
    35     * ndimensional_fvec         calculates vector function f(arg),
    36                                 stores result to fi
    37     * ndimensional_jac          calculates f[i] = fi(arg)
    38                                 jac[i,j] = df[i](arg)/d(arg[j])
     33  Callbacks for systems of functions:
     34  * ndimensional_fvec         calculates vector function f(arg),
     35                              stores result to fi
     36  * ndimensional_jac          calculates f[i] = fi(arg)
     37                              jac[i,j] = df[i](arg)/d(arg[j])
    3938                               
    40     Callbacks for  parameterized  functions,  i.e.  for  functions  which
    41     depend on two vectors: P and Q.  Gradient  and Hessian are calculated
    42     with respect to P only.
    43     * ndimensional_pfunc        calculates f(p,q),
    44                                 stores result to func
    45     * ndimensional_pgrad        calculates func = f(p,q),
    46                                 grad[i] = df(p,q)/d(p[i])
    47     * ndimensional_phess        calculates func = f(p,q),
    48                                 grad[i] = df(p,q)/d(p[i]),
    49                                 hess[i,j] = d2f(p,q)/(d(p[i])*d(p[j]))
    50 
    51     Callbacks for progress reports:
    52     * ndimensional_rep          reports current position of optimization algo   
     39  Callbacks for  parameterized  functions,  i.e.  for  functions  which
     40  depend on two vectors: P and Q.  Gradient  and Hessian are calculated
     41  with respect to P only.
     42  * ndimensional_pfunc        calculates f(p,q),
     43                              stores result to func
     44  * ndimensional_pgrad        calculates func = f(p,q),
     45                              grad[i] = df(p,q)/d(p[i])
     46  * ndimensional_phess        calculates func = f(p,q),
     47                              grad[i] = df(p,q)/d(p[i]),
     48                              hess[i,j] = d2f(p,q)/(d(p[i])*d(p[j]))
     49
     50  Callbacks for progress reports:
     51  * ndimensional_rep          reports current position of optimization algo   
    5352   
    54     Callbacks for ODE solvers:
    55     * ndimensional_ode_rp       calculates dy/dx for given y[] and x
     53  Callbacks for ODE solvers:
     54  * ndimensional_ode_rp       calculates dy/dx for given y[] and x
    5655   
    57     Callbacks for integrators:
    58     * integrator1_func          calculates f(x) for given x
    59                                 (additional parameters xminusa and bminusx
    60                                 contain x-a and b-x)
    61     ********************************************************************/
    62     public delegate void ndimensional_func (double[] arg, ref double func, object obj);
    63     public delegate void ndimensional_grad (double[] arg, ref double func, double[] grad, object obj);
    64     public delegate void ndimensional_hess (double[] arg, ref double func, double[] grad, double[,] hess, object obj);
    65    
    66     public delegate void ndimensional_fvec (double[] arg, double[] fi, object obj);
    67     public delegate void ndimensional_jac  (double[] arg, double[] fi, double[,] jac, object obj);
    68    
    69     public delegate void ndimensional_pfunc(double[] p, double[] q, ref double func, object obj);
    70     public delegate void ndimensional_pgrad(double[] p, double[] q, ref double func, double[] grad, object obj);
    71     public delegate void ndimensional_phess(double[] p, double[] q, ref double func, double[] grad, double[,] hess, object obj);
    72    
    73     public delegate void ndimensional_rep(double[] arg, double func, object obj);
    74 
    75     public delegate void ndimensional_ode_rp (double[] y, double x, double[] dy, object obj);
    76 
    77     public delegate void integrator1_func (double x, double xminusa, double bminusx, ref double f, object obj);
    78 
    79     /********************************************************************
    80     Class defining a complex number with double precision.
    81     ********************************************************************/
    82     public struct complex
    83     {
    84         public double x;
    85         public double y;
    86 
    87         public complex(double _x)
    88         {
    89             x = _x;
    90             y = 0;
     56  Callbacks for integrators:
     57  * integrator1_func          calculates f(x) for given x
     58                              (additional parameters xminusa and bminusx
     59                              contain x-a and b-x)
     60  ********************************************************************/
     61  public delegate void ndimensional_func(double[] arg, ref double func, object obj);
     62  public delegate void ndimensional_grad(double[] arg, ref double func, double[] grad, object obj);
     63  public delegate void ndimensional_hess(double[] arg, ref double func, double[] grad, double[,] hess, object obj);
     64
     65  public delegate void ndimensional_fvec(double[] arg, double[] fi, object obj);
     66  public delegate void ndimensional_jac(double[] arg, double[] fi, double[,] jac, object obj);
     67
     68  public delegate void ndimensional_pfunc(double[] p, double[] q, ref double func, object obj);
     69  public delegate void ndimensional_pgrad(double[] p, double[] q, ref double func, double[] grad, object obj);
     70  public delegate void ndimensional_phess(double[] p, double[] q, ref double func, double[] grad, double[,] hess, object obj);
     71
     72  public delegate void ndimensional_rep(double[] arg, double func, object obj);
     73
     74  public delegate void ndimensional_ode_rp(double[] y, double x, double[] dy, object obj);
     75
     76  public delegate void integrator1_func(double x, double xminusa, double bminusx, ref double f, object obj);
     77
     78  /********************************************************************
     79  Class defining a complex number with double precision.
     80  ********************************************************************/
     81  public struct complex {
     82    public double x;
     83    public double y;
     84
     85    public complex(double _x) {
     86      x = _x;
     87      y = 0;
     88    }
     89    public complex(double _x, double _y) {
     90      x = _x;
     91      y = _y;
     92    }
     93    public static implicit operator complex(double _x) {
     94      return new complex(_x);
     95    }
     96    public static bool operator ==(complex lhs, complex rhs) {
     97      return ((double)lhs.x == (double)rhs.x) & ((double)lhs.y == (double)rhs.y);
     98    }
     99    public static bool operator !=(complex lhs, complex rhs) {
     100      return ((double)lhs.x != (double)rhs.x) | ((double)lhs.y != (double)rhs.y);
     101    }
     102    public static complex operator +(complex lhs) {
     103      return lhs;
     104    }
     105    public static complex operator -(complex lhs) {
     106      return new complex(-lhs.x, -lhs.y);
     107    }
     108    public static complex operator +(complex lhs, complex rhs) {
     109      return new complex(lhs.x + rhs.x, lhs.y + rhs.y);
     110    }
     111    public static complex operator -(complex lhs, complex rhs) {
     112      return new complex(lhs.x - rhs.x, lhs.y - rhs.y);
     113    }
     114    public static complex operator *(complex lhs, complex rhs) {
     115      return new complex(lhs.x * rhs.x - lhs.y * rhs.y, lhs.x * rhs.y + lhs.y * rhs.x);
     116    }
     117    public static complex operator /(complex lhs, complex rhs) {
     118      complex result;
     119      double e;
     120      double f;
     121      if (System.Math.Abs(rhs.y) < System.Math.Abs(rhs.x)) {
     122        e = rhs.y / rhs.x;
     123        f = rhs.x + rhs.y * e;
     124        result.x = (lhs.x + lhs.y * e) / f;
     125        result.y = (lhs.y - lhs.x * e) / f;
     126      } else {
     127        e = rhs.x / rhs.y;
     128        f = rhs.y + rhs.x * e;
     129        result.x = (lhs.y + lhs.x * e) / f;
     130        result.y = (-lhs.x + lhs.y * e) / f;
     131      }
     132      return result;
     133    }
     134    public override int GetHashCode() {
     135      return x.GetHashCode() ^ y.GetHashCode();
     136    }
     137    public override bool Equals(object obj) {
     138      if (obj is byte)
     139        return Equals(new complex((byte)obj));
     140      if (obj is sbyte)
     141        return Equals(new complex((sbyte)obj));
     142      if (obj is short)
     143        return Equals(new complex((short)obj));
     144      if (obj is ushort)
     145        return Equals(new complex((ushort)obj));
     146      if (obj is int)
     147        return Equals(new complex((int)obj));
     148      if (obj is uint)
     149        return Equals(new complex((uint)obj));
     150      if (obj is long)
     151        return Equals(new complex((long)obj));
     152      if (obj is ulong)
     153        return Equals(new complex((ulong)obj));
     154      if (obj is float)
     155        return Equals(new complex((float)obj));
     156      if (obj is double)
     157        return Equals(new complex((double)obj));
     158      if (obj is decimal)
     159        return Equals(new complex((double)(decimal)obj));
     160      return base.Equals(obj);
     161    }
     162  }
     163
     164  /********************************************************************
     165  Class defining an ALGLIB exception
     166  ********************************************************************/
     167  public class alglibexception : System.Exception {
     168    public string msg;
     169    public alglibexception(string s) {
     170      msg = s;
     171    }
     172
     173  }
     174
     175  /********************************************************************
     176  reverse communication structure
     177  ********************************************************************/
     178  public class rcommstate {
     179    public rcommstate() {
     180      stage = -1;
     181      ia = new int[0];
     182      ba = new bool[0];
     183      ra = new double[0];
     184      ca = new alglib.complex[0];
     185    }
     186    public int stage;
     187    public int[] ia;
     188    public bool[] ba;
     189    public double[] ra;
     190    public alglib.complex[] ca;
     191  };
     192
     193  /********************************************************************
     194  internal functions
     195  ********************************************************************/
     196  public class ap {
     197    public static int len<T>(T[] a) { return a.Length; }
     198    public static int rows<T>(T[,] a) { return a.GetLength(0); }
     199    public static int cols<T>(T[,] a) { return a.GetLength(1); }
     200    public static void swap<T>(ref T a, ref T b) {
     201      T t = a;
     202      a = b;
     203      b = t;
     204    }
     205
     206    public static void assert(bool cond, string s) {
     207      if (!cond)
     208        throw new alglibexception(s);
     209    }
     210
     211    public static void assert(bool cond) {
     212      assert(cond, "ALGLIB: assertion failed");
     213    }
     214
     215    /****************************************************************
     216    returns dps (digits-of-precision) value corresponding to threshold.
     217    dps(0.9)  = dps(0.5)  = dps(0.1) = 0
     218    dps(0.09) = dps(0.05) = dps(0.01) = 1
     219    and so on
     220    ****************************************************************/
     221    public static int threshold2dps(double threshold) {
     222      int result = 0;
     223      double t;
     224      for (result = 0, t = 1; t / 10 > threshold * (1 + 1E-10); result++, t /= 10) ;
     225      return result;
     226    }
     227
     228    /****************************************************************
     229    prints formatted complex
     230    ****************************************************************/
     231    public static string format(complex a, int _dps) {
     232      int dps = Math.Abs(_dps);
     233      string fmt = _dps >= 0 ? "F" : "E";
     234      string fmtx = String.Format("{{0:" + fmt + "{0}}}", dps);
     235      string fmty = String.Format("{{0:" + fmt + "{0}}}", dps);
     236      string result = String.Format(fmtx, a.x) + (a.y >= 0 ? "+" : "-") + String.Format(fmty, Math.Abs(a.y)) + "i";
     237      result = result.Replace(',', '.');
     238      return result;
     239    }
     240
     241    /****************************************************************
     242    prints formatted array
     243    ****************************************************************/
     244    public static string format(bool[] a) {
     245      string[] result = new string[len(a)];
     246      int i;
     247      for (i = 0; i < len(a); i++)
     248        if (a[i])
     249          result[i] = "true";
     250        else
     251          result[i] = "false";
     252      return "{" + String.Join(",", result) + "}";
     253    }
     254
     255    /****************************************************************
     256    prints formatted array
     257    ****************************************************************/
     258    public static string format(int[] a) {
     259      string[] result = new string[len(a)];
     260      int i;
     261      for (i = 0; i < len(a); i++)
     262        result[i] = a[i].ToString();
     263      return "{" + String.Join(",", result) + "}";
     264    }
     265
     266    /****************************************************************
     267    prints formatted array
     268    ****************************************************************/
     269    public static string format(double[] a, int _dps) {
     270      int dps = Math.Abs(_dps);
     271      string sfmt = _dps >= 0 ? "F" : "E";
     272      string fmt = String.Format("{{0:" + sfmt + "{0}}}", dps);
     273      string[] result = new string[len(a)];
     274      int i;
     275      for (i = 0; i < len(a); i++) {
     276        result[i] = String.Format(fmt, a[i]);
     277        result[i] = result[i].Replace(',', '.');
     278      }
     279      return "{" + String.Join(",", result) + "}";
     280    }
     281
     282    /****************************************************************
     283    prints formatted array
     284    ****************************************************************/
     285    public static string format(complex[] a, int _dps) {
     286      int dps = Math.Abs(_dps);
     287      string fmt = _dps >= 0 ? "F" : "E";
     288      string fmtx = String.Format("{{0:" + fmt + "{0}}}", dps);
     289      string fmty = String.Format("{{0:" + fmt + "{0}}}", dps);
     290      string[] result = new string[len(a)];
     291      int i;
     292      for (i = 0; i < len(a); i++) {
     293        result[i] = String.Format(fmtx, a[i].x) + (a[i].y >= 0 ? "+" : "-") + String.Format(fmty, Math.Abs(a[i].y)) + "i";
     294        result[i] = result[i].Replace(',', '.');
     295      }
     296      return "{" + String.Join(",", result) + "}";
     297    }
     298
     299    /****************************************************************
     300    prints formatted matrix
     301    ****************************************************************/
     302    public static string format(bool[,] a) {
     303      int i, j, m, n;
     304      n = cols(a);
     305      m = rows(a);
     306      bool[] line = new bool[n];
     307      string[] result = new string[m];
     308      for (i = 0; i < m; i++) {
     309        for (j = 0; j < n; j++)
     310          line[j] = a[i, j];
     311        result[i] = format(line);
     312      }
     313      return "{" + String.Join(",", result) + "}";
     314    }
     315
     316    /****************************************************************
     317    prints formatted matrix
     318    ****************************************************************/
     319    public static string format(int[,] a) {
     320      int i, j, m, n;
     321      n = cols(a);
     322      m = rows(a);
     323      int[] line = new int[n];
     324      string[] result = new string[m];
     325      for (i = 0; i < m; i++) {
     326        for (j = 0; j < n; j++)
     327          line[j] = a[i, j];
     328        result[i] = format(line);
     329      }
     330      return "{" + String.Join(",", result) + "}";
     331    }
     332
     333    /****************************************************************
     334    prints formatted matrix
     335    ****************************************************************/
     336    public static string format(double[,] a, int dps) {
     337      int i, j, m, n;
     338      n = cols(a);
     339      m = rows(a);
     340      double[] line = new double[n];
     341      string[] result = new string[m];
     342      for (i = 0; i < m; i++) {
     343        for (j = 0; j < n; j++)
     344          line[j] = a[i, j];
     345        result[i] = format(line, dps);
     346      }
     347      return "{" + String.Join(",", result) + "}";
     348    }
     349
     350    /****************************************************************
     351    prints formatted matrix
     352    ****************************************************************/
     353    public static string format(complex[,] a, int dps) {
     354      int i, j, m, n;
     355      n = cols(a);
     356      m = rows(a);
     357      complex[] line = new complex[n];
     358      string[] result = new string[m];
     359      for (i = 0; i < m; i++) {
     360        for (j = 0; j < n; j++)
     361          line[j] = a[i, j];
     362        result[i] = format(line, dps);
     363      }
     364      return "{" + String.Join(",", result) + "}";
     365    }
     366
     367    /****************************************************************
     368    checks that matrix is symmetric.
     369    max|A-A^T| is calculated; if it is within 1.0E-14 of max|A|,
     370    matrix is considered symmetric
     371    ****************************************************************/
     372    public static bool issymmetric(double[,] a) {
     373      int i, j, n;
     374      double err, mx, v1, v2;
     375      if (rows(a) != cols(a))
     376        return false;
     377      n = rows(a);
     378      if (n == 0)
     379        return true;
     380      mx = 0;
     381      err = 0;
     382      for (i = 0; i < n; i++) {
     383        for (j = i + 1; j < n; j++) {
     384          v1 = a[i, j];
     385          v2 = a[j, i];
     386          if (!math.isfinite(v1))
     387            return false;
     388          if (!math.isfinite(v2))
     389            return false;
     390          err = Math.Max(err, Math.Abs(v1 - v2));
     391          mx = Math.Max(mx, Math.Abs(v1));
     392          mx = Math.Max(mx, Math.Abs(v2));
    91393        }
    92         public complex(double _x, double _y)
    93         {
    94             x = _x;
    95             y = _y;
     394        v1 = a[i, i];
     395        if (!math.isfinite(v1))
     396          return false;
     397        mx = Math.Max(mx, Math.Abs(v1));
     398      }
     399      if (mx == 0)
     400        return true;
     401      return err / mx <= 1.0E-14;
     402    }
     403
     404    /****************************************************************
     405    checks that matrix is Hermitian.
     406    max|A-A^H| is calculated; if it is within 1.0E-14 of max|A|,
     407    matrix is considered Hermitian
     408    ****************************************************************/
     409    public static bool ishermitian(complex[,] a) {
     410      int i, j, n;
     411      double err, mx;
     412      complex v1, v2, vt;
     413      if (rows(a) != cols(a))
     414        return false;
     415      n = rows(a);
     416      if (n == 0)
     417        return true;
     418      mx = 0;
     419      err = 0;
     420      for (i = 0; i < n; i++) {
     421        for (j = i + 1; j < n; j++) {
     422          v1 = a[i, j];
     423          v2 = a[j, i];
     424          if (!math.isfinite(v1.x))
     425            return false;
     426          if (!math.isfinite(v1.y))
     427            return false;
     428          if (!math.isfinite(v2.x))
     429            return false;
     430          if (!math.isfinite(v2.y))
     431            return false;
     432          vt.x = v1.x - v2.x;
     433          vt.y = v1.y + v2.y;
     434          err = Math.Max(err, math.abscomplex(vt));
     435          mx = Math.Max(mx, math.abscomplex(v1));
     436          mx = Math.Max(mx, math.abscomplex(v2));
    96437        }
    97         public static implicit operator complex(double _x)
    98         {
    99             return new complex(_x);
     438        v1 = a[i, i];
     439        if (!math.isfinite(v1.x))
     440          return false;
     441        if (!math.isfinite(v1.y))
     442          return false;
     443        err = Math.Max(err, Math.Abs(v1.y));
     444        mx = Math.Max(mx, math.abscomplex(v1));
     445      }
     446      if (mx == 0)
     447        return true;
     448      return err / mx <= 1.0E-14;
     449    }
     450
     451
     452    /****************************************************************
     453    Forces symmetricity by copying upper half of A to the lower one
     454    ****************************************************************/
     455    public static bool forcesymmetric(double[,] a) {
     456      int i, j, n;
     457      if (rows(a) != cols(a))
     458        return false;
     459      n = rows(a);
     460      if (n == 0)
     461        return true;
     462      for (i = 0; i < n; i++)
     463        for (j = i + 1; j < n; j++)
     464          a[i, j] = a[j, i];
     465      return true;
     466    }
     467
     468    /****************************************************************
     469    Forces Hermiticity by copying upper half of A to the lower one
     470    ****************************************************************/
     471    public static bool forcehermitian(complex[,] a) {
     472      int i, j, n;
     473      complex v;
     474      if (rows(a) != cols(a))
     475        return false;
     476      n = rows(a);
     477      if (n == 0)
     478        return true;
     479      for (i = 0; i < n; i++)
     480        for (j = i + 1; j < n; j++) {
     481          v = a[j, i];
     482          a[i, j].x = v.x;
     483          a[i, j].y = -v.y;
    100484        }
    101         public static bool operator==(complex lhs, complex rhs)
    102         {
    103             return ((double)lhs.x==(double)rhs.x) & ((double)lhs.y==(double)rhs.y);
    104         }
    105         public static bool operator!=(complex lhs, complex rhs)
    106         {
    107             return ((double)lhs.x!=(double)rhs.x) | ((double)lhs.y!=(double)rhs.y);
    108         }
    109         public static complex operator+(complex lhs)
    110         {
    111             return lhs;
    112         }
    113         public static complex operator-(complex lhs)
    114         {
    115             return new complex(-lhs.x,-lhs.y);
    116         }
    117         public static complex operator+(complex lhs, complex rhs)
    118         {
    119             return new complex(lhs.x+rhs.x,lhs.y+rhs.y);
    120         }
    121         public static complex operator-(complex lhs, complex rhs)
    122         {
    123             return new complex(lhs.x-rhs.x,lhs.y-rhs.y);
    124         }
    125         public static complex operator*(complex lhs, complex rhs)
    126         {
    127             return new complex(lhs.x*rhs.x-lhs.y*rhs.y, lhs.x*rhs.y+lhs.y*rhs.x);
    128         }
    129         public static complex operator/(complex lhs, complex rhs)
    130         {
    131             complex result;
    132             double e;
    133             double f;
    134             if( System.Math.Abs(rhs.y)<System.Math.Abs(rhs.x) )
    135             {
    136                 e = rhs.y/rhs.x;
    137                 f = rhs.x+rhs.y*e;
    138                 result.x = (lhs.x+lhs.y*e)/f;
    139                 result.y = (lhs.y-lhs.x*e)/f;
    140             }
    141             else
    142             {
    143                 e = rhs.x/rhs.y;
    144                 f = rhs.y+rhs.x*e;
    145                 result.x = (lhs.y+lhs.x*e)/f;
    146                 result.y = (-lhs.x+lhs.y*e)/f;
    147             }
    148             return result;
    149         }
    150         public override int GetHashCode()
    151         {
    152             return x.GetHashCode() ^ y.GetHashCode();
    153         }
    154         public override bool Equals(object obj)
    155         {
    156             if( obj is byte)
    157                 return Equals(new complex((byte)obj));
    158             if( obj is sbyte)
    159                 return Equals(new complex((sbyte)obj));
    160             if( obj is short)
    161                 return Equals(new complex((short)obj));
    162             if( obj is ushort)
    163                 return Equals(new complex((ushort)obj));
    164             if( obj is int)
    165                 return Equals(new complex((int)obj));
    166             if( obj is uint)
    167                 return Equals(new complex((uint)obj));
    168             if( obj is long)
    169                 return Equals(new complex((long)obj));
    170             if( obj is ulong)
    171                 return Equals(new complex((ulong)obj));
    172             if( obj is float)
    173                 return Equals(new complex((float)obj));
    174             if( obj is double)
    175                 return Equals(new complex((double)obj));
    176             if( obj is decimal)
    177                 return Equals(new complex((double)(decimal)obj));
    178             return base.Equals(obj);
    179         }   
    180     }   
    181    
    182     /********************************************************************
    183     Class defining an ALGLIB exception
    184     ********************************************************************/
    185     public class alglibexception : System.Exception
    186     {
    187         public string msg;
    188         public alglibexception(string s)
    189         {
    190             msg = s;
    191         }
    192        
    193     }
    194    
    195     /********************************************************************
    196     reverse communication structure
    197     ********************************************************************/
    198     public class rcommstate
    199     {
    200         public rcommstate()
    201         {
    202             stage = -1;
    203             ia = new int[0];
    204             ba = new bool[0];
    205             ra = new double[0];
    206             ca = new alglib.complex[0];
    207         }
    208         public int stage;
    209         public int[] ia;
    210         public bool[] ba;
    211         public double[] ra;
    212         public alglib.complex[] ca;
    213     };
    214 
    215     /********************************************************************
    216     internal functions
    217     ********************************************************************/
    218     public class ap
    219     {
    220         public static int len<T>(T[] a)
    221         { return a.Length; }
    222         public static int rows<T>(T[,] a)
    223         { return a.GetLength(0); }
    224         public static int cols<T>(T[,] a)
    225         { return a.GetLength(1); }
    226         public static void swap<T>(ref T a, ref T b)
    227         {
    228             T t = a;
    229             a = b;
    230             b = t;
    231         }
    232        
    233         public static void assert(bool cond, string s)
    234         {
    235             if( !cond )
    236                 throw new alglibexception(s);
    237         }
    238        
    239         public static void assert(bool cond)
    240         {
    241             assert(cond, "ALGLIB: assertion failed");
    242         }
    243        
    244         /****************************************************************
    245         returns dps (digits-of-precision) value corresponding to threshold.
    246         dps(0.9)  = dps(0.5)  = dps(0.1) = 0
    247         dps(0.09) = dps(0.05) = dps(0.01) = 1
    248         and so on
    249         ****************************************************************/
    250         public static int threshold2dps(double threshold)
    251         {
    252             int result = 0;
    253             double t;
    254             for (result = 0, t = 1; t / 10 > threshold*(1+1E-10); result++, t /= 10) ;
    255             return result;
    256         }
    257 
    258         /****************************************************************
    259         prints formatted complex
    260         ****************************************************************/
    261         public static string format(complex a, int _dps)
    262         {
    263             int dps = Math.Abs(_dps);
    264             string fmt = _dps>=0 ? "F" : "E";
    265             string fmtx = String.Format("{{0:"+fmt+"{0}}}", dps);
    266             string fmty = String.Format("{{0:"+fmt+"{0}}}", dps);
    267             string result = String.Format(fmtx, a.x) + (a.y >= 0 ? "+" : "-") + String.Format(fmty, Math.Abs(a.y)) + "i";
    268             result = result.Replace(',', '.');
    269             return result;
    270         }
    271 
    272         /****************************************************************
    273         prints formatted array
    274         ****************************************************************/
    275         public static string format(bool[] a)
    276         {
    277             string[] result = new string[len(a)];
    278             int i;
    279             for(i=0; i<len(a); i++)
    280                 if( a[i] )
    281                     result[i] = "true";
    282                 else
    283                     result[i] = "false";
    284             return "{"+String.Join(",",result)+"}";
    285         }
    286        
    287         /****************************************************************
    288         prints formatted array
    289         ****************************************************************/
    290         public static string format(int[] a)
    291         {
    292             string[] result = new string[len(a)];
    293             int i;
    294             for (i = 0; i < len(a); i++)
    295                 result[i] = a[i].ToString();
    296             return "{" + String.Join(",", result) + "}";
    297         }
    298 
    299         /****************************************************************
    300         prints formatted array
    301         ****************************************************************/
    302         public static string format(double[] a, int _dps)
    303         {
    304             int dps = Math.Abs(_dps);
    305             string sfmt = _dps >= 0 ? "F" : "E";
    306             string fmt = String.Format("{{0:" + sfmt + "{0}}}", dps);
    307             string[] result = new string[len(a)];
    308             int i;
    309             for (i = 0; i < len(a); i++)
    310             {
    311                 result[i] = String.Format(fmt, a[i]);
    312                 result[i] = result[i].Replace(',', '.');
    313             }
    314             return "{" + String.Join(",", result) + "}";
    315         }
    316 
    317         /****************************************************************
    318         prints formatted array
    319         ****************************************************************/
    320         public static string format(complex[] a, int _dps)
    321         {
    322             int dps = Math.Abs(_dps);
    323             string fmt = _dps >= 0 ? "F" : "E";
    324             string fmtx = String.Format("{{0:"+fmt+"{0}}}", dps);
    325             string fmty = String.Format("{{0:"+fmt+"{0}}}", dps);
    326             string[] result = new string[len(a)];
    327             int i;
    328             for (i = 0; i < len(a); i++)
    329             {
    330                 result[i] = String.Format(fmtx, a[i].x) + (a[i].y >= 0 ? "+" : "-") + String.Format(fmty, Math.Abs(a[i].y)) + "i";
    331                 result[i] = result[i].Replace(',', '.');
    332             }
    333             return "{" + String.Join(",", result) + "}";
    334         }
    335 
    336         /****************************************************************
    337         prints formatted matrix
    338         ****************************************************************/
    339         public static string format(bool[,] a)
    340         {
    341             int i, j, m, n;
    342             n = cols(a);
    343             m = rows(a);
    344             bool[] line = new bool[n];
    345             string[] result = new string[m];
    346             for (i = 0; i < m; i++)
    347             {
    348                 for (j = 0; j < n; j++)
    349                     line[j] = a[i, j];
    350                 result[i] = format(line);
    351             }
    352             return "{" + String.Join(",", result) + "}";
    353         }
    354 
    355         /****************************************************************
    356         prints formatted matrix
    357         ****************************************************************/
    358         public static string format(int[,] a)
    359         {
    360             int i, j, m, n;
    361             n = cols(a);
    362             m = rows(a);
    363             int[] line = new int[n];
    364             string[] result = new string[m];
    365             for (i = 0; i < m; i++)
    366             {
    367                 for (j = 0; j < n; j++)
    368                     line[j] = a[i, j];
    369                 result[i] = format(line);
    370             }
    371             return "{" + String.Join(",", result) + "}";
    372         }
    373 
    374         /****************************************************************
    375         prints formatted matrix
    376         ****************************************************************/
    377         public static string format(double[,] a, int dps)
    378         {
    379             int i, j, m, n;
    380             n = cols(a);
    381             m = rows(a);
    382             double[] line = new double[n];
    383             string[] result = new string[m];
    384             for (i = 0; i < m; i++)
    385             {
    386                 for (j = 0; j < n; j++)
    387                     line[j] = a[i, j];
    388                 result[i] = format(line, dps);
    389             }
    390             return "{" + String.Join(",", result) + "}";
    391         }
    392 
    393         /****************************************************************
    394         prints formatted matrix
    395         ****************************************************************/
    396         public static string format(complex[,] a, int dps)
    397         {
    398             int i, j, m, n;
    399             n = cols(a);
    400             m = rows(a);
    401             complex[] line = new complex[n];
    402             string[] result = new string[m];
    403             for (i = 0; i < m; i++)
    404             {
    405                 for (j = 0; j < n; j++)
    406                     line[j] = a[i, j];
    407                 result[i] = format(line, dps);
    408             }
    409             return "{" + String.Join(",", result) + "}";
    410         }
    411 
    412         /****************************************************************
    413         checks that matrix is symmetric.
    414         max|A-A^T| is calculated; if it is within 1.0E-14 of max|A|,
    415         matrix is considered symmetric
    416         ****************************************************************/
    417         public static bool issymmetric(double[,] a)
    418         {
    419             int i, j, n;
    420             double err, mx, v1, v2;
    421             if( rows(a)!=cols(a) )
    422                 return false;
    423             n = rows(a);
    424             if( n==0 )
    425                 return true;
    426             mx = 0;
    427             err = 0;
    428             for( i=0; i<n; i++)
    429             {
    430                 for(j=i+1; j<n; j++)
    431                 {
    432                     v1 = a[i,j];
    433                     v2 = a[j,i];
    434                     if( !math.isfinite(v1) )
    435                         return false;
    436                     if( !math.isfinite(v2) )
    437                         return false;
    438                     err = Math.Max(err, Math.Abs(v1-v2));
    439                     mx  = Math.Max(mx,  Math.Abs(v1));
    440                     mx  = Math.Max(mx,  Math.Abs(v2));
    441                 }
    442                 v1 = a[i,i];
    443                 if( !math.isfinite(v1) )
    444                     return false;
    445                 mx = Math.Max(mx, Math.Abs(v1));
    446             }
    447             if( mx==0 )
    448                 return true;
    449             return err/mx<=1.0E-14;
    450         }
    451        
    452         /****************************************************************
    453         checks that matrix is Hermitian.
    454         max|A-A^H| is calculated; if it is within 1.0E-14 of max|A|,
    455         matrix is considered Hermitian
    456         ****************************************************************/
    457         public static bool ishermitian(complex[,] a)
    458         {
    459             int i, j, n;
    460             double err, mx;
    461             complex v1, v2, vt;
    462             if( rows(a)!=cols(a) )
    463                 return false;
    464             n = rows(a);
    465             if( n==0 )
    466                 return true;
    467             mx = 0;
    468             err = 0;
    469             for( i=0; i<n; i++)
    470             {
    471                 for(j=i+1; j<n; j++)
    472                 {
    473                     v1 = a[i,j];
    474                     v2 = a[j,i];
    475                     if( !math.isfinite(v1.x) )
    476                         return false;
    477                     if( !math.isfinite(v1.y) )
    478                         return false;
    479                     if( !math.isfinite(v2.x) )
    480                         return false;
    481                     if( !math.isfinite(v2.y) )
    482                         return false;
    483                     vt.x = v1.x-v2.x;
    484                     vt.y = v1.y+v2.y;
    485                     err = Math.Max(err, math.abscomplex(vt));
    486                     mx  = Math.Max(mx,  math.abscomplex(v1));
    487                     mx  = Math.Max(mx,  math.abscomplex(v2));
    488                 }
    489                 v1 = a[i,i];
    490                 if( !math.isfinite(v1.x) )
    491                     return false;
    492                 if( !math.isfinite(v1.y) )
    493                     return false;
    494                 err = Math.Max(err, Math.Abs(v1.y));
    495                 mx = Math.Max(mx, math.abscomplex(v1));
    496             }
    497             if( mx==0 )
    498                 return true;
    499             return err/mx<=1.0E-14;
    500         }
    501        
    502        
    503         /****************************************************************
    504         Forces symmetricity by copying upper half of A to the lower one
    505         ****************************************************************/
    506         public static bool forcesymmetric(double[,] a)
    507         {
    508             int i, j, n;
    509             if( rows(a)!=cols(a) )
    510                 return false;
    511             n = rows(a);
    512             if( n==0 )
    513                 return true;
    514             for( i=0; i<n; i++)
    515                 for(j=i+1; j<n; j++)
    516                     a[i,j] = a[j,i];
    517             return true;
    518         }
    519        
    520         /****************************************************************
    521         Forces Hermiticity by copying upper half of A to the lower one
    522         ****************************************************************/
    523         public static bool forcehermitian(complex[,] a)
    524         {
    525             int i, j, n;
    526             complex v;
    527             if( rows(a)!=cols(a) )
    528                 return false;
    529             n = rows(a);
    530             if( n==0 )
    531                 return true;
    532             for( i=0; i<n; i++)
    533                 for(j=i+1; j<n; j++)
    534                 {
    535                     v = a[j,i];
    536                     a[i,j].x = v.x;
    537                     a[i,j].y = -v.y;
    538                 }
    539             return true;
    540         }
    541     };
    542    
    543     /********************************************************************
    544     math functions
    545     ********************************************************************/
    546     public class math
    547     {
    548         //public static System.Random RndObject = new System.Random(System.DateTime.Now.Millisecond);
    549         public static System.Random rndobject = new System.Random(System.DateTime.Now.Millisecond + 1000*System.DateTime.Now.Second + 60*1000*System.DateTime.Now.Minute);
    550 
    551         public const double machineepsilon = 5E-16;
    552         public const double maxrealnumber = 1E300;
    553         public const double minrealnumber = 1E-300;
    554        
    555         public static bool isfinite(double d)
    556         {
    557             return !System.Double.IsNaN(d) && !System.Double.IsInfinity(d);
    558         }
    559        
    560         public static double randomreal()
    561         {
    562             double r = 0;
    563             lock(rndobject){ r = rndobject.NextDouble(); }
    564             return r;
    565         }
    566         public static int randominteger(int N)
    567         {
    568             int r = 0;
    569             lock(rndobject){ r = rndobject.Next(N); }
    570             return r;
    571         }
    572         public static double sqr(double X)
    573         {
    574             return X*X;
    575         }       
    576         public static double abscomplex(complex z)
    577         {
    578             double w;
    579             double xabs;
    580             double yabs;
    581             double v;
    582    
    583             xabs = System.Math.Abs(z.x);
    584             yabs = System.Math.Abs(z.y);
    585             w = xabs>yabs ? xabs : yabs;
    586             v = xabs<yabs ? xabs : yabs;
    587             if( v==0 )
    588                 return w;
    589             else
    590             {
    591                 double t = v/w;
    592                 return w*System.Math.Sqrt(1+t*t);
    593             }
    594         }
    595         public static complex conj(complex z)
    596         {
    597             return new complex(z.x, -z.y);
    598         }   
    599         public static complex csqr(complex z)
    600         {
    601             return new complex(z.x*z.x-z.y*z.y, 2*z.x*z.y);
    602         }
    603 
    604     }
    605    
    606     /********************************************************************
    607     serializer object (should not be used directly)
    608     ********************************************************************/
    609     public class serializer
    610     {
    611         enum SMODE { DEFAULT, ALLOC, TO_STRING, FROM_STRING };
    612         private const int SER_ENTRIES_PER_ROW = 5;
    613         private const int SER_ENTRY_LENGTH    = 11;
    614        
    615         private SMODE mode;
    616         private int entries_needed;
    617         private int entries_saved;
    618         private int bytes_asked;
    619         private int bytes_written;
    620         private int bytes_read;
    621         private char[] out_str;
    622         private char[] in_str;
    623        
    624         public serializer()
    625         {
    626             mode = SMODE.DEFAULT;
    627             entries_needed = 0;
    628             bytes_asked = 0;
    629         }
    630 
    631         public void alloc_start()
    632         {
    633             entries_needed = 0;
    634             bytes_asked = 0;
    635             mode = SMODE.ALLOC;
    636         }
    637 
    638         public void alloc_entry()
    639         {
    640             if( mode!=SMODE.ALLOC )
    641                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    642             entries_needed++;
    643         }
    644 
    645         private int get_alloc_size()
    646         {
    647             int rows, lastrowsize, result;
    648            
    649             // check and change mode
    650             if( mode!=SMODE.ALLOC )
    651                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    652            
    653             // if no entries needes (degenerate case)
    654             if( entries_needed==0 )
    655             {
    656                 bytes_asked = 1;
    657                 return bytes_asked;
    658             }
    659            
    660             // non-degenerate case
    661             rows = entries_needed/SER_ENTRIES_PER_ROW;
    662             lastrowsize = SER_ENTRIES_PER_ROW;
    663             if( entries_needed%SER_ENTRIES_PER_ROW!=0 )
    664             {
    665                 lastrowsize = entries_needed%SER_ENTRIES_PER_ROW;
    666                 rows++;
    667             }
    668            
    669             // calculate result size
    670             result  = ((rows-1)*SER_ENTRIES_PER_ROW+lastrowsize)*SER_ENTRY_LENGTH;
    671             result +=  (rows-1)*(SER_ENTRIES_PER_ROW-1)+(lastrowsize-1);
    672             result += rows*2;
    673             bytes_asked = result;
    674             return result;
    675         }
    676 
    677         public void sstart_str()
    678         {
    679             int allocsize = get_alloc_size();
    680            
    681             // check and change mode
    682             if( mode!=SMODE.ALLOC )
    683                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    684             mode = SMODE.TO_STRING;
    685            
    686             // other preparations
    687             out_str = new char[allocsize];
    688             entries_saved = 0;
    689             bytes_written = 0;
    690         }
    691 
    692         public void ustart_str(string s)
    693         {
    694             // check and change mode
    695             if( mode!=SMODE.DEFAULT )
    696                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    697             mode = SMODE.FROM_STRING;
    698            
    699             in_str = s.ToCharArray();
    700             bytes_read = 0;
    701         }
    702 
    703         public void serialize_bool(bool v)
    704         {
    705             if( mode!=SMODE.TO_STRING )
    706                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    707             bool2str(v, out_str, ref bytes_written);
    708             entries_saved++;
    709             if( entries_saved%SER_ENTRIES_PER_ROW!=0 )
    710             {
    711                 out_str[bytes_written] = ' ';
    712                 bytes_written++;
    713             }
    714             else
    715             {
    716                 out_str[bytes_written+0] = '\r';
    717                 out_str[bytes_written+1] = '\n';
    718                 bytes_written+=2;
    719             }           
    720         }
    721 
    722         public void serialize_int(int v)
    723         {
    724             if( mode!=SMODE.TO_STRING )
    725                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    726             int2str(v, out_str, ref bytes_written);
    727             entries_saved++;
    728             if( entries_saved%SER_ENTRIES_PER_ROW!=0 )
    729             {
    730                 out_str[bytes_written] = ' ';
    731                 bytes_written++;
    732             }
    733             else
    734             {
    735                 out_str[bytes_written+0] = '\r';
    736                 out_str[bytes_written+1] = '\n';
    737                 bytes_written+=2;
    738             }
    739         }
    740 
    741         public void serialize_double(double v)
    742         {
    743             if( mode!=SMODE.TO_STRING )
    744                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    745             double2str(v, out_str, ref bytes_written);
    746             entries_saved++;
    747             if( entries_saved%SER_ENTRIES_PER_ROW!=0 )
    748             {
    749                 out_str[bytes_written] = ' ';
    750                 bytes_written++;
    751             }
    752             else
    753             {
    754                 out_str[bytes_written+0] = '\r';
    755                 out_str[bytes_written+1] = '\n';
    756                 bytes_written+=2;
    757             }
    758         }
    759 
    760         public bool unserialize_bool()
    761         {
    762             if( mode!=SMODE.FROM_STRING )
    763                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    764             return str2bool(in_str, ref bytes_read);
    765         }
    766 
    767         public int unserialize_int()
    768         {
    769             if( mode!=SMODE.FROM_STRING )
    770                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    771             return str2int(in_str, ref bytes_read);
    772         }
    773 
    774         public double unserialize_double()
    775         {
    776             if( mode!=SMODE.FROM_STRING )
    777                 throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
    778             return str2double(in_str, ref bytes_read);
    779         }
    780 
    781         public void stop()
    782         {
    783         }
    784 
    785         public string get_string()
    786         {
    787             return new string(out_str, 0, bytes_written);
    788         }
    789 
    790 
    791         /************************************************************************
    792         This function converts six-bit value (from 0 to 63)  to  character  (only
    793         digits, lowercase and uppercase letters, minus and underscore are used).
    794 
    795         If v is negative or greater than 63, this function returns '?'.
    796         ************************************************************************/
    797         private static char[] _sixbits2char_tbl = new char[64]{
     485      return true;
     486    }
     487  };
     488
     489  /********************************************************************
     490  math functions
     491  ********************************************************************/
     492  public class math {
     493    //public static System.Random RndObject = new System.Random(System.DateTime.Now.Millisecond);
     494    [ThreadStatic] //mkommend: added thread static attribute to RNG to have a separate instance per thread and allow modification of the seed
     495    public static System.Random rndobject = new System.Random(System.DateTime.Now.Millisecond + 1000 * System.DateTime.Now.Second + 60 * 1000 * System.DateTime.Now.Minute);
     496
     497    public const double machineepsilon = 5E-16;
     498    public const double maxrealnumber = 1E300;
     499    public const double minrealnumber = 1E-300;
     500
     501    public static bool isfinite(double d) {
     502      return !System.Double.IsNaN(d) && !System.Double.IsInfinity(d);
     503    }
     504
     505    public static double randomreal() {
     506      double r = 0;
     507      lock (rndobject) { r = rndobject.NextDouble(); }
     508      return r;
     509    }
     510    public static int randominteger(int N) {
     511      int r = 0;
     512      lock (rndobject) { r = rndobject.Next(N); }
     513      return r;
     514    }
     515    public static double sqr(double X) {
     516      return X * X;
     517    }
     518    public static double abscomplex(complex z) {
     519      double w;
     520      double xabs;
     521      double yabs;
     522      double v;
     523
     524      xabs = System.Math.Abs(z.x);
     525      yabs = System.Math.Abs(z.y);
     526      w = xabs > yabs ? xabs : yabs;
     527      v = xabs < yabs ? xabs : yabs;
     528      if (v == 0)
     529        return w;
     530      else {
     531        double t = v / w;
     532        return w * System.Math.Sqrt(1 + t * t);
     533      }
     534    }
     535    public static complex conj(complex z) {
     536      return new complex(z.x, -z.y);
     537    }
     538    public static complex csqr(complex z) {
     539      return new complex(z.x * z.x - z.y * z.y, 2 * z.x * z.y);
     540    }
     541
     542  }
     543
     544  /********************************************************************
     545  serializer object (should not be used directly)
     546  ********************************************************************/
     547  public class serializer {
     548    enum SMODE { DEFAULT, ALLOC, TO_STRING, FROM_STRING };
     549    private const int SER_ENTRIES_PER_ROW = 5;
     550    private const int SER_ENTRY_LENGTH = 11;
     551
     552    private SMODE mode;
     553    private int entries_needed;
     554    private int entries_saved;
     555    private int bytes_asked;
     556    private int bytes_written;
     557    private int bytes_read;
     558    private char[] out_str;
     559    private char[] in_str;
     560
     561    public serializer() {
     562      mode = SMODE.DEFAULT;
     563      entries_needed = 0;
     564      bytes_asked = 0;
     565    }
     566
     567    public void alloc_start() {
     568      entries_needed = 0;
     569      bytes_asked = 0;
     570      mode = SMODE.ALLOC;
     571    }
     572
     573    public void alloc_entry() {
     574      if (mode != SMODE.ALLOC)
     575        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     576      entries_needed++;
     577    }
     578
     579    private int get_alloc_size() {
     580      int rows, lastrowsize, result;
     581
     582      // check and change mode
     583      if (mode != SMODE.ALLOC)
     584        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     585
     586      // if no entries needes (degenerate case)
     587      if (entries_needed == 0) {
     588        bytes_asked = 1;
     589        return bytes_asked;
     590      }
     591
     592      // non-degenerate case
     593      rows = entries_needed / SER_ENTRIES_PER_ROW;
     594      lastrowsize = SER_ENTRIES_PER_ROW;
     595      if (entries_needed % SER_ENTRIES_PER_ROW != 0) {
     596        lastrowsize = entries_needed % SER_ENTRIES_PER_ROW;
     597        rows++;
     598      }
     599
     600      // calculate result size
     601      result = ((rows - 1) * SER_ENTRIES_PER_ROW + lastrowsize) * SER_ENTRY_LENGTH;
     602      result += (rows - 1) * (SER_ENTRIES_PER_ROW - 1) + (lastrowsize - 1);
     603      result += rows * 2;
     604      bytes_asked = result;
     605      return result;
     606    }
     607
     608    public void sstart_str() {
     609      int allocsize = get_alloc_size();
     610
     611      // check and change mode
     612      if (mode != SMODE.ALLOC)
     613        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     614      mode = SMODE.TO_STRING;
     615
     616      // other preparations
     617      out_str = new char[allocsize];
     618      entries_saved = 0;
     619      bytes_written = 0;
     620    }
     621
     622    public void ustart_str(string s) {
     623      // check and change mode
     624      if (mode != SMODE.DEFAULT)
     625        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     626      mode = SMODE.FROM_STRING;
     627
     628      in_str = s.ToCharArray();
     629      bytes_read = 0;
     630    }
     631
     632    public void serialize_bool(bool v) {
     633      if (mode != SMODE.TO_STRING)
     634        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     635      bool2str(v, out_str, ref bytes_written);
     636      entries_saved++;
     637      if (entries_saved % SER_ENTRIES_PER_ROW != 0) {
     638        out_str[bytes_written] = ' ';
     639        bytes_written++;
     640      } else {
     641        out_str[bytes_written + 0] = '\r';
     642        out_str[bytes_written + 1] = '\n';
     643        bytes_written += 2;
     644      }
     645    }
     646
     647    public void serialize_int(int v) {
     648      if (mode != SMODE.TO_STRING)
     649        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     650      int2str(v, out_str, ref bytes_written);
     651      entries_saved++;
     652      if (entries_saved % SER_ENTRIES_PER_ROW != 0) {
     653        out_str[bytes_written] = ' ';
     654        bytes_written++;
     655      } else {
     656        out_str[bytes_written + 0] = '\r';
     657        out_str[bytes_written + 1] = '\n';
     658        bytes_written += 2;
     659      }
     660    }
     661
     662    public void serialize_double(double v) {
     663      if (mode != SMODE.TO_STRING)
     664        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     665      double2str(v, out_str, ref bytes_written);
     666      entries_saved++;
     667      if (entries_saved % SER_ENTRIES_PER_ROW != 0) {
     668        out_str[bytes_written] = ' ';
     669        bytes_written++;
     670      } else {
     671        out_str[bytes_written + 0] = '\r';
     672        out_str[bytes_written + 1] = '\n';
     673        bytes_written += 2;
     674      }
     675    }
     676
     677    public bool unserialize_bool() {
     678      if (mode != SMODE.FROM_STRING)
     679        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     680      return str2bool(in_str, ref bytes_read);
     681    }
     682
     683    public int unserialize_int() {
     684      if (mode != SMODE.FROM_STRING)
     685        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     686      return str2int(in_str, ref bytes_read);
     687    }
     688
     689    public double unserialize_double() {
     690      if (mode != SMODE.FROM_STRING)
     691        throw new alglib.alglibexception("ALGLIB: internal error during (un)serialization");
     692      return str2double(in_str, ref bytes_read);
     693    }
     694
     695    public void stop() {
     696    }
     697
     698    public string get_string() {
     699      return new string(out_str, 0, bytes_written);
     700    }
     701
     702
     703    /************************************************************************
     704    This function converts six-bit value (from 0 to 63)  to  character  (only
     705    digits, lowercase and uppercase letters, minus and underscore are used).
     706
     707    If v is negative or greater than 63, this function returns '?'.
     708    ************************************************************************/
     709    private static char[] _sixbits2char_tbl = new char[64]{
    798710                '0', '1', '2', '3', '4', '5', '6', '7',
    799711                '8', '9', 'A', 'B', 'C', 'D', 'E', 'F',
     
    804716                'm', 'n', 'o', 'p', 'q', 'r', 's', 't',
    805717                'u', 'v', 'w', 'x', 'y', 'z', '-', '_' };
    806         private static char sixbits2char(int v)
    807         {
    808             if( v<0 || v>63 )
    809                 return '?';
    810             return _sixbits2char_tbl[v];
    811         }
    812        
    813         /************************************************************************
    814         This function converts character to six-bit value (from 0 to 63).
    815 
    816         This function is inverse of ae_sixbits2char()
    817         If c is not correct character, this function returns -1.
    818         ************************************************************************/
    819         private static int[] _char2sixbits_tbl = new int[128] {
     718    private static char sixbits2char(int v) {
     719      if (v < 0 || v > 63)
     720        return '?';
     721      return _sixbits2char_tbl[v];
     722    }
     723
     724    /************************************************************************
     725    This function converts character to six-bit value (from 0 to 63).
     726
     727    This function is inverse of ae_sixbits2char()
     728    If c is not correct character, this function returns -1.
     729    ************************************************************************/
     730    private static int[] _char2sixbits_tbl = new int[128] {
    820731            -1, -1, -1, -1, -1, -1, -1, -1,
    821732            -1, -1, -1, -1, -1, -1, -1, -1,
     
    834745            51, 52, 53, 54, 55, 56, 57, 58,
    835746            59, 60, 61, -1, -1, -1, -1, -1 };
    836         private static int char2sixbits(char c)
    837         {
    838             return (c>=0 && c<127) ? _char2sixbits_tbl[c] : -1;
     747    private static int char2sixbits(char c) {
     748      return (c >= 0 && c < 127) ? _char2sixbits_tbl[c] : -1;
     749    }
     750
     751    /************************************************************************
     752    This function converts three bytes (24 bits) to four six-bit values
     753    (24 bits again).
     754
     755    src         array
     756    src_offs    offset of three-bytes chunk
     757    dst         array for ints
     758    dst_offs    offset of four-ints chunk
     759    ************************************************************************/
     760    private static void threebytes2foursixbits(byte[] src, int src_offs, int[] dst, int dst_offs) {
     761      dst[dst_offs + 0] = src[src_offs + 0] & 0x3F;
     762      dst[dst_offs + 1] = (src[src_offs + 0] >> 6) | ((src[src_offs + 1] & 0x0F) << 2);
     763      dst[dst_offs + 2] = (src[src_offs + 1] >> 4) | ((src[src_offs + 2] & 0x03) << 4);
     764      dst[dst_offs + 3] = src[src_offs + 2] >> 2;
     765    }
     766
     767    /************************************************************************
     768    This function converts four six-bit values (24 bits) to three bytes
     769    (24 bits again).
     770
     771    src         pointer to four ints
     772    src_offs    offset of the chunk
     773    dst         pointer to three bytes
     774    dst_offs    offset of the chunk
     775    ************************************************************************/
     776    private static void foursixbits2threebytes(int[] src, int src_offs, byte[] dst, int dst_offs) {
     777      dst[dst_offs + 0] = (byte)(src[src_offs + 0] | ((src[src_offs + 1] & 0x03) << 6));
     778      dst[dst_offs + 1] = (byte)((src[src_offs + 1] >> 2) | ((src[src_offs + 2] & 0x0F) << 4));
     779      dst[dst_offs + 2] = (byte)((src[src_offs + 2] >> 4) | (src[src_offs + 3] << 2));
     780    }
     781
     782    /************************************************************************
     783    This function serializes boolean value into buffer
     784
     785    v           boolean value to be serialized
     786    buf         buffer, at least 11 characters wide
     787    offs        offset in the buffer
     788       
     789    after return from this function, offs points to the char's past the value
     790    being read.
     791    ************************************************************************/
     792    private static void bool2str(bool v, char[] buf, ref int offs) {
     793      char c = v ? '1' : '0';
     794      int i;
     795      for (i = 0; i < SER_ENTRY_LENGTH; i++)
     796        buf[offs + i] = c;
     797      offs += SER_ENTRY_LENGTH;
     798    }
     799
     800    /************************************************************************
     801    This function unserializes boolean value from buffer
     802
     803    buf         buffer which contains value; leading spaces/tabs/newlines are
     804                ignored, traling spaces/tabs/newlines are treated as  end  of
     805                the boolean value.
     806    offs        offset in the buffer
     807       
     808    after return from this function, offs points to the char's past the value
     809    being read.
     810
     811    This function raises an error in case unexpected symbol is found
     812    ************************************************************************/
     813    private static bool str2bool(char[] buf, ref int offs) {
     814      bool was0, was1;
     815      string emsg = "ALGLIB: unable to read boolean value from stream";
     816
     817      was0 = false;
     818      was1 = false;
     819      while (buf[offs] == ' ' || buf[offs] == '\t' || buf[offs] == '\n' || buf[offs] == '\r')
     820        offs++;
     821      while (buf[offs] != ' ' && buf[offs] != '\t' && buf[offs] != '\n' && buf[offs] != '\r' && buf[offs] != 0) {
     822        if (buf[offs] == '0') {
     823          was0 = true;
     824          offs++;
     825          continue;
    839826        }
     827        if (buf[offs] == '1') {
     828          was1 = true;
     829          offs++;
     830          continue;
     831        }
     832        throw new alglib.alglibexception(emsg);
     833      }
     834      if ((!was0) && (!was1))
     835        throw new alglib.alglibexception(emsg);
     836      if (was0 && was1)
     837        throw new alglib.alglibexception(emsg);
     838      return was1 ? true : false;
     839    }
     840
     841    /************************************************************************
     842    This function serializes integer value into buffer
     843
     844    v           integer value to be serialized
     845    buf         buffer, at least 11 characters wide
     846    offs        offset in the buffer
    840847       
    841         /************************************************************************
    842         This function converts three bytes (24 bits) to four six-bit values
    843         (24 bits again).
    844 
    845         src         array
    846         src_offs    offset of three-bytes chunk
    847         dst         array for ints
    848         dst_offs    offset of four-ints chunk
    849         ************************************************************************/
    850         private static void threebytes2foursixbits(byte[] src, int src_offs, int[] dst, int dst_offs)
    851         {
    852             dst[dst_offs+0] =  src[src_offs+0] & 0x3F;
    853             dst[dst_offs+1] = (src[src_offs+0]>>6) | ((src[src_offs+1]&0x0F)<<2);
    854             dst[dst_offs+2] = (src[src_offs+1]>>4) | ((src[src_offs+2]&0x03)<<4);
    855             dst[dst_offs+3] =  src[src_offs+2]>>2;
     848    after return from this function, offs points to the char's past the value
     849    being read.
     850
     851    This function raises an error in case unexpected symbol is found
     852    ************************************************************************/
     853    private static void int2str(int v, char[] buf, ref int offs) {
     854      int i;
     855      byte[] _bytes = System.BitConverter.GetBytes((int)v);
     856      byte[] bytes = new byte[9];
     857      int[] sixbits = new int[12];
     858      byte c;
     859
     860      //
     861      // copy v to array of bytes, sign extending it and
     862      // converting to little endian order. Additionally,
     863      // we set 9th byte to zero in order to simplify
     864      // conversion to six-bit representation
     865      //
     866      if (!System.BitConverter.IsLittleEndian)
     867        System.Array.Reverse(_bytes);
     868      c = v < 0 ? (byte)0xFF : (byte)0x00;
     869      for (i = 0; i < sizeof(int); i++)
     870        bytes[i] = _bytes[i];
     871      for (i = sizeof(int); i < 8; i++)
     872        bytes[i] = c;
     873      bytes[8] = 0;
     874
     875      //
     876      // convert to six-bit representation, output
     877      //
     878      // NOTE: last 12th element of sixbits is always zero, we do not output it
     879      //
     880      threebytes2foursixbits(bytes, 0, sixbits, 0);
     881      threebytes2foursixbits(bytes, 3, sixbits, 4);
     882      threebytes2foursixbits(bytes, 6, sixbits, 8);
     883      for (i = 0; i < SER_ENTRY_LENGTH; i++)
     884        buf[offs + i] = sixbits2char(sixbits[i]);
     885      offs += SER_ENTRY_LENGTH;
     886    }
     887
     888    /************************************************************************
     889    This function unserializes integer value from string
     890
     891    buf         buffer which contains value; leading spaces/tabs/newlines are
     892                ignored, traling spaces/tabs/newlines are treated as  end  of
     893                the integer value.
     894    offs        offset in the buffer
     895       
     896    after return from this function, offs points to the char's past the value
     897    being read.
     898
     899    This function raises an error in case unexpected symbol is found
     900    ************************************************************************/
     901    private static int str2int(char[] buf, ref int offs) {
     902      string emsg = "ALGLIB: unable to read integer value from stream";
     903      string emsg3264 = "ALGLIB: unable to read integer value from stream (value does not fit into 32 bits)";
     904      int[] sixbits = new int[12];
     905      byte[] bytes = new byte[9];
     906      byte[] _bytes = new byte[sizeof(int)];
     907      int sixbitsread, i;
     908      byte c;
     909
     910      //
     911      // 1. skip leading spaces
     912      // 2. read and decode six-bit digits
     913      // 3. set trailing digits to zeros
     914      // 4. convert to little endian 64-bit integer representation
     915      // 5. check that we fit into int
     916      // 6. convert to big endian representation, if needed
     917      //
     918      sixbitsread = 0;
     919      while (buf[offs] == ' ' || buf[offs] == '\t' || buf[offs] == '\n' || buf[offs] == '\r')
     920        offs++;
     921      while (buf[offs] != ' ' && buf[offs] != '\t' && buf[offs] != '\n' && buf[offs] != '\r' && buf[offs] != 0) {
     922        int d;
     923        d = char2sixbits(buf[offs]);
     924        if (d < 0 || sixbitsread >= SER_ENTRY_LENGTH)
     925          throw new alglib.alglibexception(emsg);
     926        sixbits[sixbitsread] = d;
     927        sixbitsread++;
     928        offs++;
     929      }
     930      if (sixbitsread == 0)
     931        throw new alglib.alglibexception(emsg);
     932      for (i = sixbitsread; i < 12; i++)
     933        sixbits[i] = 0;
     934      foursixbits2threebytes(sixbits, 0, bytes, 0);
     935      foursixbits2threebytes(sixbits, 4, bytes, 3);
     936      foursixbits2threebytes(sixbits, 8, bytes, 6);
     937      c = (bytes[sizeof(int) - 1] & 0x80) != 0 ? (byte)0xFF : (byte)0x00;
     938      for (i = sizeof(int); i < 8; i++)
     939        if (bytes[i] != c)
     940          throw new alglib.alglibexception(emsg3264);
     941      for (i = 0; i < sizeof(int); i++)
     942        _bytes[i] = bytes[i];
     943      if (!System.BitConverter.IsLittleEndian)
     944        System.Array.Reverse(_bytes);
     945      return System.BitConverter.ToInt32(_bytes, 0);
     946    }
     947
     948
     949    /************************************************************************
     950    This function serializes double value into buffer
     951
     952    v           double value to be serialized
     953    buf         buffer, at least 11 characters wide
     954    offs        offset in the buffer
     955       
     956    after return from this function, offs points to the char's past the value
     957    being read.
     958    ************************************************************************/
     959    private static void double2str(double v, char[] buf, ref int offs) {
     960      int i;
     961      int[] sixbits = new int[12];
     962      byte[] bytes = new byte[9];
     963
     964      //
     965      // handle special quantities
     966      //
     967      if (System.Double.IsNaN(v)) {
     968        buf[offs + 0] = '.';
     969        buf[offs + 1] = 'n';
     970        buf[offs + 2] = 'a';
     971        buf[offs + 3] = 'n';
     972        buf[offs + 4] = '_';
     973        buf[offs + 5] = '_';
     974        buf[offs + 6] = '_';
     975        buf[offs + 7] = '_';
     976        buf[offs + 8] = '_';
     977        buf[offs + 9] = '_';
     978        buf[offs + 10] = '_';
     979        offs += SER_ENTRY_LENGTH;
     980        return;
     981      }
     982      if (System.Double.IsPositiveInfinity(v)) {
     983        buf[offs + 0] = '.';
     984        buf[offs + 1] = 'p';
     985        buf[offs + 2] = 'o';
     986        buf[offs + 3] = 's';
     987        buf[offs + 4] = 'i';
     988        buf[offs + 5] = 'n';
     989        buf[offs + 6] = 'f';
     990        buf[offs + 7] = '_';
     991        buf[offs + 8] = '_';
     992        buf[offs + 9] = '_';
     993        buf[offs + 10] = '_';
     994        offs += SER_ENTRY_LENGTH;
     995        return;
     996      }
     997      if (System.Double.IsNegativeInfinity(v)) {
     998        buf[offs + 0] = '.';
     999        buf[offs + 1] = 'n';
     1000        buf[offs + 2] = 'e';
     1001        buf[offs + 3] = 'g';
     1002        buf[offs + 4] = 'i';
     1003        buf[offs + 5] = 'n';
     1004        buf[offs + 6] = 'f';
     1005        buf[offs + 7] = '_';
     1006        buf[offs + 8] = '_';
     1007        buf[offs + 9] = '_';
     1008        buf[offs + 10] = '_';
     1009        offs += SER_ENTRY_LENGTH;
     1010        return;
     1011      }
     1012
     1013      //
     1014      // process general case:
     1015      // 1. copy v to array of chars
     1016      // 2. set 9th byte to zero in order to simplify conversion to six-bit representation
     1017      // 3. convert to little endian (if needed)
     1018      // 4. convert to six-bit representation
     1019      //    (last 12th element of sixbits is always zero, we do not output it)
     1020      //
     1021      byte[] _bytes = System.BitConverter.GetBytes((double)v);
     1022      if (!System.BitConverter.IsLittleEndian)
     1023        System.Array.Reverse(_bytes);
     1024      for (i = 0; i < sizeof(double); i++)
     1025        bytes[i] = _bytes[i];
     1026      for (i = sizeof(double); i < 9; i++)
     1027        bytes[i] = 0;
     1028      threebytes2foursixbits(bytes, 0, sixbits, 0);
     1029      threebytes2foursixbits(bytes, 3, sixbits, 4);
     1030      threebytes2foursixbits(bytes, 6, sixbits, 8);
     1031      for (i = 0; i < SER_ENTRY_LENGTH; i++)
     1032        buf[offs + i] = sixbits2char(sixbits[i]);
     1033      offs += SER_ENTRY_LENGTH;
     1034    }
     1035
     1036    /************************************************************************
     1037    This function unserializes double value from string
     1038
     1039    buf         buffer which contains value; leading spaces/tabs/newlines are
     1040                ignored, traling spaces/tabs/newlines are treated as  end  of
     1041                the double value.
     1042    offs        offset in the buffer
     1043       
     1044    after return from this function, offs points to the char's past the value
     1045    being read.
     1046
     1047    This function raises an error in case unexpected symbol is found
     1048    ************************************************************************/
     1049    private static double str2double(char[] buf, ref int offs) {
     1050      string emsg = "ALGLIB: unable to read double value from stream";
     1051      int[] sixbits = new int[12];
     1052      byte[] bytes = new byte[9];
     1053      byte[] _bytes = new byte[sizeof(double)];
     1054      int sixbitsread, i;
     1055
     1056
     1057      //
     1058      // skip leading spaces
     1059      //
     1060      while (buf[offs] == ' ' || buf[offs] == '\t' || buf[offs] == '\n' || buf[offs] == '\r')
     1061        offs++;
     1062
     1063
     1064      //
     1065      // Handle special cases
     1066      //
     1067      if (buf[offs] == '.') {
     1068        string s = new string(buf, offs, SER_ENTRY_LENGTH);
     1069        if (s == ".nan_______") {
     1070          offs += SER_ENTRY_LENGTH;
     1071          return System.Double.NaN;
    8561072        }
    857 
    858         /************************************************************************
    859         This function converts four six-bit values (24 bits) to three bytes
    860         (24 bits again).
    861 
    862         src         pointer to four ints
    863         src_offs    offset of the chunk
    864         dst         pointer to three bytes
    865         dst_offs    offset of the chunk
    866         ************************************************************************/
    867         private static void foursixbits2threebytes(int[] src, int src_offs, byte[] dst, int dst_offs)
    868         {
    869             dst[dst_offs+0] =      (byte)(src[src_offs+0] | ((src[src_offs+1]&0x03)<<6));
    870             dst[dst_offs+1] = (byte)((src[src_offs+1]>>2) | ((src[src_offs+2]&0x0F)<<4));
    871             dst[dst_offs+2] = (byte)((src[src_offs+2]>>4) |  (src[src_offs+3]<<2));
     1073        if (s == ".posinf____") {
     1074          offs += SER_ENTRY_LENGTH;
     1075          return System.Double.PositiveInfinity;
    8721076        }
    873 
    874         /************************************************************************
    875         This function serializes boolean value into buffer
    876 
    877         v           boolean value to be serialized
    878         buf         buffer, at least 11 characters wide
    879         offs        offset in the buffer
    880        
    881         after return from this function, offs points to the char's past the value
    882         being read.
    883         ************************************************************************/
    884         private static void bool2str(bool v, char[] buf, ref int offs)
    885         {
    886             char c = v ? '1' : '0';
    887             int i;
    888             for(i=0; i<SER_ENTRY_LENGTH; i++)
    889                 buf[offs+i] = c;
    890             offs += SER_ENTRY_LENGTH;
     1077        if (s == ".neginf____") {
     1078          offs += SER_ENTRY_LENGTH;
     1079          return System.Double.NegativeInfinity;
    8911080        }
    892 
    893         /************************************************************************
    894         This function unserializes boolean value from buffer
    895 
    896         buf         buffer which contains value; leading spaces/tabs/newlines are
    897                     ignored, traling spaces/tabs/newlines are treated as  end  of
    898                     the boolean value.
    899         offs        offset in the buffer
    900        
    901         after return from this function, offs points to the char's past the value
    902         being read.
    903 
    904         This function raises an error in case unexpected symbol is found
    905         ************************************************************************/
    906         private static bool str2bool(char[] buf, ref int offs)
    907         {
    908             bool was0, was1;
    909             string emsg = "ALGLIB: unable to read boolean value from stream";
    910            
    911             was0 = false;
    912             was1 = false;
    913             while( buf[offs]==' ' || buf[offs]=='\t' || buf[offs]=='\n' || buf[offs]=='\r' )
    914                 offs++;
    915             while( buf[offs]!=' ' && buf[offs]!='\t' && buf[offs]!='\n' && buf[offs]!='\r' && buf[offs]!=0 )
    916             {
    917                 if( buf[offs]=='0' )
    918                 {
    919                     was0 = true;
    920                     offs++;
    921                     continue;
    922                 }
    923                 if( buf[offs]=='1' )
    924                 {
    925                     was1 = true;
    926                     offs++;
    927                     continue;
    928                 }
    929                 throw new alglib.alglibexception(emsg);
    930             }
    931             if( (!was0) && (!was1) )
    932                 throw new alglib.alglibexception(emsg);
    933             if( was0 && was1 )
    934                 throw new alglib.alglibexception(emsg);
    935             return was1 ? true : false;
    936         }
    937 
    938         /************************************************************************
    939         This function serializes integer value into buffer
    940 
    941         v           integer value to be serialized
    942         buf         buffer, at least 11 characters wide
    943         offs        offset in the buffer
    944        
    945         after return from this function, offs points to the char's past the value
    946         being read.
    947 
    948         This function raises an error in case unexpected symbol is found
    949         ************************************************************************/
    950         private static void int2str(int v, char[] buf, ref int offs)
    951         {
    952             int i;
    953             byte[] _bytes = System.BitConverter.GetBytes((int)v);
    954             byte[]  bytes = new byte[9];
    955             int[] sixbits = new int[12];
    956             byte c;
    957            
    958             //
    959             // copy v to array of bytes, sign extending it and
    960             // converting to little endian order. Additionally,
    961             // we set 9th byte to zero in order to simplify
    962             // conversion to six-bit representation
    963             //
    964             if( !System.BitConverter.IsLittleEndian )
    965                 System.Array.Reverse(_bytes);
    966             c = v<0 ? (byte)0xFF : (byte)0x00;
    967             for(i=0; i<sizeof(int); i++)
    968                 bytes[i] = _bytes[i];
    969             for(i=sizeof(int); i<8; i++)
    970                 bytes[i] = c;
    971             bytes[8] = 0;
    972            
    973             //
    974             // convert to six-bit representation, output
    975             //
    976             // NOTE: last 12th element of sixbits is always zero, we do not output it
    977             //
    978             threebytes2foursixbits(bytes, 0, sixbits, 0);
    979             threebytes2foursixbits(bytes, 3, sixbits, 4);
    980             threebytes2foursixbits(bytes, 6, sixbits, 8);       
    981             for(i=0; i<SER_ENTRY_LENGTH; i++)
    982                 buf[offs+i] = sixbits2char(sixbits[i]);
    983             offs += SER_ENTRY_LENGTH;
    984         }
    985 
    986         /************************************************************************
    987         This function unserializes integer value from string
    988 
    989         buf         buffer which contains value; leading spaces/tabs/newlines are
    990                     ignored, traling spaces/tabs/newlines are treated as  end  of
    991                     the integer value.
    992         offs        offset in the buffer
    993        
    994         after return from this function, offs points to the char's past the value
    995         being read.
    996 
    997         This function raises an error in case unexpected symbol is found
    998         ************************************************************************/
    999         private static int str2int(char[] buf, ref int offs)
    1000         {
    1001             string emsg =       "ALGLIB: unable to read integer value from stream";
    1002             string emsg3264 =   "ALGLIB: unable to read integer value from stream (value does not fit into 32 bits)";
    1003             int[] sixbits = new int[12];
    1004             byte[] bytes = new byte[9];
    1005             byte[] _bytes = new byte[sizeof(int)];
    1006             int sixbitsread, i;
    1007             byte c;
    1008            
    1009             //
    1010             // 1. skip leading spaces
    1011             // 2. read and decode six-bit digits
    1012             // 3. set trailing digits to zeros
    1013             // 4. convert to little endian 64-bit integer representation
    1014             // 5. check that we fit into int
    1015             // 6. convert to big endian representation, if needed
    1016             //
    1017             sixbitsread = 0;
    1018             while( buf[offs]==' ' || buf[offs]=='\t' || buf[offs]=='\n' || buf[offs]=='\r' )
    1019                 offs++;
    1020             while( buf[offs]!=' ' && buf[offs]!='\t' && buf[offs]!='\n' && buf[offs]!='\r' && buf[offs]!=0 )
    1021             {
    1022                 int d;
    1023                 d = char2sixbits(buf[offs]);
    1024                 if( d<0 || sixbitsread>=SER_ENTRY_LENGTH )
    1025                     throw new alglib.alglibexception(emsg);
    1026                 sixbits[sixbitsread] = d;
    1027                 sixbitsread++;
    1028                 offs++;
    1029             }
    1030             if( sixbitsread==0 )
    1031                 throw new alglib.alglibexception(emsg);
    1032             for(i=sixbitsread; i<12; i++)
    1033                 sixbits[i] = 0;
    1034             foursixbits2threebytes(sixbits, 0, bytes, 0);
    1035             foursixbits2threebytes(sixbits, 4, bytes, 3);
    1036             foursixbits2threebytes(sixbits, 8, bytes, 6);
    1037             c = (bytes[sizeof(int)-1] & 0x80)!=0 ? (byte)0xFF : (byte)0x00;
    1038             for(i=sizeof(int); i<8; i++)
    1039                 if( bytes[i]!=c )
    1040                     throw new alglib.alglibexception(emsg3264);
    1041             for(i=0; i<sizeof(int); i++)
    1042                 _bytes[i] = bytes[i];       
    1043             if( !System.BitConverter.IsLittleEndian )
    1044                 System.Array.Reverse(_bytes);
    1045             return System.BitConverter.ToInt32(_bytes,0);
    1046         }   
    1047        
    1048        
    1049         /************************************************************************
    1050         This function serializes double value into buffer
    1051 
    1052         v           double value to be serialized
    1053         buf         buffer, at least 11 characters wide
    1054         offs        offset in the buffer
    1055        
    1056         after return from this function, offs points to the char's past the value
    1057         being read.
    1058         ************************************************************************/
    1059         private static void double2str(double v, char[] buf, ref int offs)
    1060         {
    1061             int i;
    1062             int[] sixbits = new int[12];
    1063             byte[] bytes = new byte[9];
    1064 
    1065             //
    1066             // handle special quantities
    1067             //
    1068             if( System.Double.IsNaN(v) )
    1069             {
    1070                 buf[offs+0] = '.';
    1071                 buf[offs+1] = 'n';
    1072                 buf[offs+2] = 'a';
    1073                 buf[offs+3] = 'n';
    1074                 buf[offs+4] = '_';
    1075                 buf[offs+5] = '_';
    1076                 buf[offs+6] = '_';
    1077                 buf[offs+7] = '_';
    1078                 buf[offs+8] = '_';
    1079                 buf[offs+9] = '_';
    1080                 buf[offs+10] = '_';
    1081                 offs += SER_ENTRY_LENGTH;
    1082                 return;
    1083             }
    1084             if( System.Double.IsPositiveInfinity(v) )
    1085             {
    1086                 buf[offs+0] = '.';
    1087                 buf[offs+1] = 'p';
    1088                 buf[offs+2] = 'o';
    1089                 buf[offs+3] = 's';
    1090                 buf[offs+4] = 'i';
    1091                 buf[offs+5] = 'n';
    1092                 buf[offs+6] = 'f';
    1093                 buf[offs+7] = '_';
    1094                 buf[offs+8] = '_';
    1095                 buf[offs+9] = '_';
    1096                 buf[offs+10] = '_';
    1097                 offs += SER_ENTRY_LENGTH;
    1098                 return;
    1099             }
    1100             if( System.Double.IsNegativeInfinity(v) )
    1101             {
    1102                 buf[offs+0] = '.';
    1103                 buf[offs+1] = 'n';
    1104                 buf[offs+2] = 'e';
    1105                 buf[offs+3] = 'g';
    1106                 buf[offs+4] = 'i';
    1107                 buf[offs+5] = 'n';
    1108                 buf[offs+6] = 'f';
    1109                 buf[offs+7] = '_';
    1110                 buf[offs+8] = '_';
    1111                 buf[offs+9] = '_';
    1112                 buf[offs+10] = '_';
    1113                 offs += SER_ENTRY_LENGTH;
    1114                 return;
    1115             }
    1116            
    1117             //
    1118             // process general case:
    1119             // 1. copy v to array of chars
    1120             // 2. set 9th byte to zero in order to simplify conversion to six-bit representation
    1121             // 3. convert to little endian (if needed)
    1122             // 4. convert to six-bit representation
    1123             //    (last 12th element of sixbits is always zero, we do not output it)
    1124             //
    1125             byte[] _bytes = System.BitConverter.GetBytes((double)v);
    1126             if( !System.BitConverter.IsLittleEndian )
    1127                 System.Array.Reverse(_bytes);
    1128             for(i=0; i<sizeof(double); i++)
    1129                 bytes[i] = _bytes[i];
    1130             for(i=sizeof(double); i<9; i++)
    1131                 bytes[i] = 0;
    1132             threebytes2foursixbits(bytes, 0, sixbits, 0);
    1133             threebytes2foursixbits(bytes, 3, sixbits, 4);
    1134             threebytes2foursixbits(bytes, 6, sixbits, 8);
    1135             for(i=0; i<SER_ENTRY_LENGTH; i++)
    1136                 buf[offs+i] = sixbits2char(sixbits[i]);
    1137             offs += SER_ENTRY_LENGTH;
    1138         }
    1139 
    1140         /************************************************************************
    1141         This function unserializes double value from string
    1142 
    1143         buf         buffer which contains value; leading spaces/tabs/newlines are
    1144                     ignored, traling spaces/tabs/newlines are treated as  end  of
    1145                     the double value.
    1146         offs        offset in the buffer
    1147        
    1148         after return from this function, offs points to the char's past the value
    1149         being read.
    1150 
    1151         This function raises an error in case unexpected symbol is found
    1152         ************************************************************************/
    1153         private static double str2double(char[] buf, ref int offs)
    1154         {
    1155             string emsg = "ALGLIB: unable to read double value from stream";
    1156             int[] sixbits = new int[12];
    1157             byte[]  bytes = new byte[9];
    1158             byte[] _bytes = new byte[sizeof(double)];
    1159             int sixbitsread, i;
    1160            
    1161            
    1162             //
    1163             // skip leading spaces
    1164             //
    1165             while( buf[offs]==' ' || buf[offs]=='\t' || buf[offs]=='\n' || buf[offs]=='\r' )
    1166                 offs++;
    1167            
    1168              
    1169             //
    1170             // Handle special cases
    1171             //
    1172             if( buf[offs]=='.' )
    1173             {
    1174                 string s = new string(buf, offs, SER_ENTRY_LENGTH);
    1175                 if( s==".nan_______" )
    1176                 {
    1177                     offs += SER_ENTRY_LENGTH;
    1178                     return System.Double.NaN;
    1179                 }
    1180                 if( s==".posinf____" )
    1181                 {
    1182                     offs += SER_ENTRY_LENGTH;
    1183                     return System.Double.PositiveInfinity;
    1184                 }
    1185                 if( s==".neginf____" )
    1186                 {
    1187                     offs += SER_ENTRY_LENGTH;
    1188                     return System.Double.NegativeInfinity;
    1189                 }
    1190                 throw new alglib.alglibexception(emsg);
    1191             }
    1192            
    1193             //
    1194             // General case:
    1195             // 1. read and decode six-bit digits
    1196             // 2. check that all 11 digits were read
    1197             // 3. set last 12th digit to zero (needed for simplicity of conversion)
    1198             // 4. convert to 8 bytes
    1199             // 5. convert to big endian representation, if needed
    1200             //
    1201             sixbitsread = 0;
    1202             while( buf[offs]!=' ' && buf[offs]!='\t' && buf[offs]!='\n' && buf[offs]!='\r' && buf[offs]!=0 )
    1203             {
    1204                 int d;
    1205                 d = char2sixbits(buf[offs]);
    1206                 if( d<0 || sixbitsread>=SER_ENTRY_LENGTH )
    1207                     throw new alglib.alglibexception(emsg);
    1208                 sixbits[sixbitsread] = d;
    1209                 sixbitsread++;
    1210                 offs++;
    1211             }
    1212             if( sixbitsread!=SER_ENTRY_LENGTH )
    1213                 throw new alglib.alglibexception(emsg);
    1214             sixbits[SER_ENTRY_LENGTH] = 0;
    1215             foursixbits2threebytes(sixbits, 0, bytes, 0);
    1216             foursixbits2threebytes(sixbits, 4, bytes, 3);
    1217             foursixbits2threebytes(sixbits, 8, bytes, 6);
    1218             for(i=0; i<sizeof(double); i++)
    1219                 _bytes[i] = bytes[i];       
    1220             if( !System.BitConverter.IsLittleEndian )
    1221                 System.Array.Reverse(_bytes);       
    1222             return System.BitConverter.ToDouble(_bytes,0);
    1223         }
    1224     }}
     1081        throw new alglib.alglibexception(emsg);
     1082      }
     1083
     1084      //
     1085      // General case:
     1086      // 1. read and decode six-bit digits
     1087      // 2. check that all 11 digits were read
     1088      // 3. set last 12th digit to zero (needed for simplicity of conversion)
     1089      // 4. convert to 8 bytes
     1090      // 5. convert to big endian representation, if needed
     1091      //
     1092      sixbitsread = 0;
     1093      while (buf[offs] != ' ' && buf[offs] != '\t' && buf[offs] != '\n' && buf[offs] != '\r' && buf[offs] != 0) {
     1094        int d;
     1095        d = char2sixbits(buf[offs]);
     1096        if (d < 0 || sixbitsread >= SER_ENTRY_LENGTH)
     1097          throw new alglib.alglibexception(emsg);
     1098        sixbits[sixbitsread] = d;
     1099        sixbitsread++;
     1100        offs++;
     1101      }
     1102      if (sixbitsread != SER_ENTRY_LENGTH)
     1103        throw new alglib.alglibexception(emsg);
     1104      sixbits[SER_ENTRY_LENGTH] = 0;
     1105      foursixbits2threebytes(sixbits, 0, bytes, 0);
     1106      foursixbits2threebytes(sixbits, 4, bytes, 3);
     1107      foursixbits2threebytes(sixbits, 8, bytes, 6);
     1108      for (i = 0; i < sizeof(double); i++)
     1109        _bytes[i] = bytes[i];
     1110      if (!System.BitConverter.IsLittleEndian)
     1111        System.Array.Reverse(_bytes);
     1112      return System.BitConverter.ToDouble(_bytes, 0);
     1113    }
     1114  }
     1115}
Note: See TracChangeset for help on using the changeset viewer.