/*************************************************************************
        This function calculates values of the RBF model at the given point.

        This is general function which can be used for arbitrary NX (dimension  of 
        the space of arguments) and NY (dimension of the function itself). However
        when  you  have  NY=1  you  may  find more convenient to use RBFCalc2() or 
        RBFCalc3().

        This function returns 0.0 when model is not initialized.

        INPUT PARAMETERS:
            S       -   RBF model
            X       -   coordinates, array[NX].
                        X may have more than NX elements, in this case only 
                        leading NX will be used.

        OUTPUT PARAMETERS:
            Y       -   function value, array[NY]. Y is out-parameter and 
                        reallocated after call to this function. In case you  want
                        to reuse previously allocated Y, you may use RBFCalcBuf(),
                        which reallocates Y only when it is too small.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfcalc(rbfmodel s,
            double[] x,
            ref double[] y)
        {
            y = new double[0];

            alglib.ap.assert(alglib.ap.len(x)>=s.nx, "RBFCalc: Length(X)<NX");
            alglib.ap.assert(apserv.isfinitevector(x, s.nx), "RBFCalc: X contains infinite or NaN values");
            rbfcalcbuf(s, x, ref y);
        }

        /*************************************************************************
        This function calculates values of the RBF model in the given point.

        This function should be used when we have NY=1 (scalar function) and  NX=2
        (2-dimensional space). If you have 3-dimensional space, use RBFCalc3(). If
        you have general situation (NX-dimensional space, NY-dimensional function)
        you should use general, less efficient implementation RBFCalc().

        If  you  want  to  calculate  function  values  many times, consider using 
        RBFGridCalc2(), which is far more efficient than many subsequent calls  to
        RBFCalc2().

        This function returns 0.0 when:
        * model is not initialized
        * NX<>2
         *NY<>1

        INPUT PARAMETERS:
            S       -   RBF model
            X0      -   first coordinate, finite number
            X1      -   second coordinate, finite number

        RESULT:
            value of the model or 0.0 (as defined above)

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static double rbfcalc2(rbfmodel s,
            double x0,
            double x1)
        {
            double result = 0;
            int i = 0;
            int j = 0;
            int lx = 0;
            int tg = 0;
            double d2 = 0;
            double t = 0;
            double bfcur = 0;
            double rcur = 0;

            alglib.ap.assert(math.isfinite(x0), "RBFCalc2: invalid value for X0 (X0 is Inf)!");
            alglib.ap.assert(math.isfinite(x1), "RBFCalc2: invalid value for X1 (X1 is Inf)!");
            if( s.ny!=1 || s.nx!=2 )
            {
                result = 0;
                return result;
            }
            result = s.v[0,0]*x0+s.v[0,1]*x1+s.v[0,mxnx];
            if( s.nc==0 )
            {
                return result;
            }
            apserv.rvectorsetlengthatleast(ref s.calcbufxcx, mxnx);
            for(i=0; i<=mxnx-1; i++)
            {
                s.calcbufxcx[i] = 0.0;
            }
            s.calcbufxcx[0] = x0;
            s.calcbufxcx[1] = x1;
            lx = nearestneighbor.kdtreequeryrnn(s.tree, s.calcbufxcx, s.rmax*rbffarradius, true);
            nearestneighbor.kdtreequeryresultsx(s.tree, ref s.calcbufx);
            nearestneighbor.kdtreequeryresultstags(s.tree, ref s.calcbuftags);
            for(i=0; i<=lx-1; i++)
            {
                tg = s.calcbuftags[i];
                d2 = math.sqr(x0-s.calcbufx[i,0])+math.sqr(x1-s.calcbufx[i,1]);
                rcur = s.wr[tg,0];
                bfcur = Math.Exp(-(d2/(rcur*rcur)));
                for(j=0; j<=s.nl-1; j++)
                {
                    result = result+bfcur*s.wr[tg,1+j];
                    rcur = 0.5*rcur;
                    t = bfcur*bfcur;
                    bfcur = t*t;
                }
            }
            return result;
        }

        /*************************************************************************
        This function calculates values of the RBF model in the given point.

        This function should be used when we have NY=1 (scalar function) and  NX=3
        (3-dimensional space). If you have 2-dimensional space, use RBFCalc2(). If
        you have general situation (NX-dimensional space, NY-dimensional function)
        you should use general, less efficient implementation RBFCalc().

        This function returns 0.0 when:
        * model is not initialized
        * NX<>3
         *NY<>1

        INPUT PARAMETERS:
            S       -   RBF model
            X0      -   first coordinate, finite number
            X1      -   second coordinate, finite number
            X2      -   third coordinate, finite number

        RESULT:
            value of the model or 0.0 (as defined above)

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static double rbfcalc3(rbfmodel s,
            double x0,
            double x1,
            double x2)
        {
            double result = 0;
            int i = 0;
            int j = 0;
            int lx = 0;
            int tg = 0;
            double t = 0;
            double rcur = 0;
            double bf = 0;

            alglib.ap.assert(math.isfinite(x0), "RBFCalc3: invalid value for X0 (X0 is Inf or NaN)!");
            alglib.ap.assert(math.isfinite(x1), "RBFCalc3: invalid value for X1 (X1 is Inf or NaN)!");
            alglib.ap.assert(math.isfinite(x2), "RBFCalc3: invalid value for X2 (X2 is Inf or NaN)!");
            if( s.ny!=1 || s.nx!=3 )
            {
                result = 0;
                return result;
            }
            result = s.v[0,0]*x0+s.v[0,1]*x1+s.v[0,2]*x2+s.v[0,mxnx];
            if( s.nc==0 )
            {
                return result;
            }
            
            //
            // calculating value for F(X)
            //
            apserv.rvectorsetlengthatleast(ref s.calcbufxcx, mxnx);
            for(i=0; i<=mxnx-1; i++)
            {
                s.calcbufxcx[i] = 0.0;
            }
            s.calcbufxcx[0] = x0;
            s.calcbufxcx[1] = x1;
            s.calcbufxcx[2] = x2;
            lx = nearestneighbor.kdtreequeryrnn(s.tree, s.calcbufxcx, s.rmax*rbffarradius, true);
            nearestneighbor.kdtreequeryresultsx(s.tree, ref s.calcbufx);
            nearestneighbor.kdtreequeryresultstags(s.tree, ref s.calcbuftags);
            for(i=0; i<=lx-1; i++)
            {
                tg = s.calcbuftags[i];
                rcur = s.wr[tg,0];
                bf = Math.Exp(-((math.sqr(x0-s.calcbufx[i,0])+math.sqr(x1-s.calcbufx[i,1])+math.sqr(x2-s.calcbufx[i,2]))/math.sqr(rcur)));
                for(j=0; j<=s.nl-1; j++)
                {
                    result = result+bf*s.wr[tg,1+j];
                    t = bf*bf;
                    bf = t*t;
                }
            }
            return result;
        }

        /*************************************************************************
        This function changes centers allocation algorithm to one which  allocates
        centers exactly at the dataset points (one input point = one center). This
        function won't have effect until next call to RBFBuildModel().

        INPUT PARAMETERS:
            S       -   RBF model, initialized by RBFCreate() call

        NOTE: this   function  has   some   serialization-related  subtleties.  We
              recommend you to study serialization examples from ALGLIB  Reference
              Manual if you want to perform serialization of your models.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        private static void rbfgridpoints(rbfmodel s)
        {
            s.gridtype = 2;
        }

        /*************************************************************************
        This function sets centers, defined by user.

        This function overrides results of the previous calls, i.e. multiple calls
        of this function will result in only the last set being added.

        INPUT PARAMETERS:
            S      -   RBF model, initialized by RBFCreate() call.
          CXY      -   centers, array[N,NX]. Centers must be distinct, 
                       non-distinct centers are not supported.
           NC      -   number of centers

        NOTE: this   function  has   some   serialization-related  subtleties.  We
              recommend you to study serialization examples from ALGLIB  Reference
              Manual if you want to perform serialization of your models.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        private static void rbfgridspecial(rbfmodel s,
            double[,] cxy,
            int nc)
        {
            int i = 0;
            int j = 0;

            alglib.ap.assert(nc>0, "RBFGridSpecial: N<0");
            alglib.ap.assert(nc<=alglib.ap.rows(cxy), "RBFGridSpecial: Length(CXY)<NC");
            s.gridtype = 3;
            s.nc = nc;
            s.xc = new double[s.nc, mxnx];
            for(i=0; i<=s.nc-1; i++)
            {
                for(j=0; j<=mxnx-1; j++)
                {
                    s.xc[i,j] = 0;
                }
                for(j=0; j<=s.nx-1; j++)
                {
                    s.xc[i,j] = cxy[i,j];
                }
            }
        }

        /*************************************************************************
        This function "unpacks" RBF model by extracting its coefficients.

        INPUT PARAMETERS:
            S       -   RBF model

        OUTPUT PARAMETERS:
            NX      -   dimensionality of argument
            NY      -   dimensionality of the target function
            XWR     -   model information, array[NC,NX+NY+1].
                        One row of the array corresponds to one basis function:
                        * first NX columns  - coordinates of the center 
                        * next NY columns   - weights, one per dimension of the 
                                              function being modelled
                        * last column       - radius, same for all dimensions of
                                              the function being modelled
            NC      -   number of the centers
            V       -   polynomial  term , array[NY,NX+1]. One row per one 
                        dimension of the function being modelled. First NX 
                        elements are linear coefficients, V[NX] is equal to the 
                        constant part.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfunpack(rbfmodel s,
            ref int nx,
            ref int ny,
            ref double[,] xwr,
            ref int nc,
            ref double[,] v)
        {
            int i = 0;
            int j = 0;
            double rcur = 0;
            int i_ = 0;
            int i1_ = 0;

            nx = 0;
            ny = 0;
            xwr = new double[0,0];
            nc = 0;
            v = new double[0,0];

            nx = s.nx;
            ny = s.ny;
            nc = s.nc;
            
            //
            // Fill V
            //
            v = new double[s.ny, s.nx+1];
            for(i=0; i<=s.ny-1; i++)
            {
                for(i_=0; i_<=s.nx-1;i_++)
                {
                    v[i,i_] = s.v[i,i_];
                }
                v[i,s.nx] = s.v[i,mxnx];
            }
            
            //
            // Fill XWR and V
            //
            if( nc*s.nl>0 )
            {
                xwr = new double[s.nc*s.nl, s.nx+s.ny+1];
                for(i=0; i<=s.nc-1; i++)
                {
                    rcur = s.wr[i,0];
                    for(j=0; j<=s.nl-1; j++)
                    {
                        for(i_=0; i_<=s.nx-1;i_++)
                        {
                            xwr[i*s.nl+j,i_] = s.xc[i,i_];
                        }
                        i1_ = (1+j*s.ny) - (s.nx);
                        for(i_=s.nx; i_<=s.nx+s.ny-1;i_++)
                        {
                            xwr[i*s.nl+j,i_] = s.wr[i,i_+i1_];
                        }
                        xwr[i*s.nl+j,s.nx+s.ny] = rcur;
                        rcur = 0.5*rcur;
                    }
                }
            }
        }

        /*************************************************************************
        Serializer: serialization

          -- ALGLIB --
             Copyright 02.02.2012 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfserialize(alglib.serializer s,
            rbfmodel model)
        {
            
            //
            // Header
            //
            s.serialize_int(scodes.getrbfserializationcode());
            s.serialize_int(rbffirstversion);
            
            //
            // Data
            //
            s.serialize_int(model.nx);
            s.serialize_int(model.ny);
            s.serialize_int(model.nc);
            s.serialize_int(model.nl);
            nearestneighbor.kdtreeserialize(s, model.tree);
            apserv.serializerealmatrix(s, model.xc, -1, -1);
            apserv.serializerealmatrix(s, model.wr, -1, -1);
            s.serialize_double(model.rmax);
            apserv.serializerealmatrix(s, model.v, -1, -1);
        }

        /*************************************************************************
        This  function  sets  RBF interpolation algorithm. ALGLIB supports several
        RBF algorithms with different properties.

        This  algorithm is called RBF-ML. It builds  multilayer  RBF  model,  i.e.
        model with subsequently decreasing  radii,  which  allows  us  to  combine
        smoothness (due to  large radii of  the first layers) with  exactness (due
        to small radii of the last layers) and fast convergence.

        Internally RBF-ML uses many different  means  of acceleration, from sparse
        matrices  to  KD-trees,  which  results in algorithm whose working time is
        roughly proportional to N*log(N)*Density*RBase^2*NLayers,  where  N  is  a
        number of points, Density is an average density if points per unit of  the
        interpolation space, RBase is an initial radius, NLayers is  a  number  of
        layers.

        RBF-ML is good for following kinds of interpolation problems:
        1. "exact" problems (perfect fit) with well separated points
        2. least squares problems with arbitrary distribution of points (algorithm
           gives  perfect  fit  where it is possible, and resorts to least squares
           fit in the hard areas).
        3. noisy problems where  we  want  to  apply  some  controlled  amount  of
           smoothing.

        INPUT PARAMETERS:
            S       -   RBF model, initialized by RBFCreate() call
            RBase   -   RBase parameter, RBase>0
            NLayers -   NLayers parameter, NLayers>0, recommended value  to  start
                        with - about 5.
            LambdaV -   regularization value, can be useful when  solving  problem
                        in the least squares sense.  Optimal  lambda  is  problem-
                        dependent and require trial and error. In our  experience,
                        good lambda can be as large as 0.1, and you can use  0.001
                        as initial guess.
                        Default  value  - 0.01, which is used when LambdaV is  not
                        given.  You  can  specify  zero  value,  but  it  is   not
                        recommended to do so.

        TUNING ALGORITHM

        In order to use this algorithm you have to choose three parameters:
        * initial radius RBase
        * number of layers in the model NLayers
        * regularization coefficient LambdaV

        Initial radius is easy to choose - you can pick any number  several  times
        larger  than  the  average  distance between points. Algorithm won't break
        down if you choose radius which is too large (model construction time will
        increase, but model will be built correctly).

        Choose such number of layers that RLast=RBase/2^(NLayers-1)  (radius  used
        by  the  last  layer)  will  be  smaller than the typical distance between
        points.  In  case  model  error  is  too large, you can increase number of
        layers.  Having  more  layers  will make model construction and evaluation
        proportionally slower, but it will allow you to have model which precisely
        fits your data. From the other side, if you want to  suppress  noise,  you
        can DECREASE number of layers to make your model less flexible.

        Regularization coefficient LambdaV controls smoothness of  the  individual
        models built for each layer. We recommend you to use default value in case
        you don't want to tune this parameter,  because  having  non-zero  LambdaV
        accelerates and stabilizes internal iterative algorithm. In case you  want
        to suppress noise you can use  LambdaV  as  additional  parameter  (larger
        value = more smoothness) to tune.

        TYPICAL ERRORS

        1. Using  initial  radius  which is too large. Memory requirements  of the
           RBF-ML are roughly proportional to N*Density*RBase^2 (where Density  is
           an average density of points per unit of the interpolation  space).  In
           the extreme case of the very large RBase we will need O(N^2)  units  of
           memory - and many layers in order to decrease radius to some reasonably
           small value.

        2. Using too small number of layers - RBF models with large radius are not
           flexible enough to reproduce small variations in the  target  function.
           You  need  many  layers  with  different radii, from large to small, in
           order to have good model.

        3. Using  initial  radius  which  is  too  small.  You will get model with
           "holes" in the areas which are too far away from interpolation centers.
           However, algorithm will work correctly (and quickly) in this case.

        4. Using too many layers - you will get too large and too slow model. This
           model  will  perfectly  reproduce  your function, but maybe you will be
           able to achieve similar results with less layers (and less memory).
           
          -- ALGLIB --
             Copyright 02.03.2012 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfsetalgomultilayer(rbfmodel s,
            double rbase,
            int nlayers,
            double lambdav)
        {
            alglib.ap.assert(math.isfinite(rbase), "RBFSetAlgoMultiLayer: RBase is infinite or NaN");
            alglib.ap.assert((double)(rbase)>(double)(0), "RBFSetAlgoMultiLayer: RBase<=0");
            alglib.ap.assert(nlayers>=0, "RBFSetAlgoMultiLayer: NLayers<0");
            alglib.ap.assert(math.isfinite(lambdav), "RBFSetAlgoMultiLayer: LambdaV is infinite or NAN");
            alglib.ap.assert((double)(lambdav)>=(double)(0), "RBFSetAlgoMultiLayer: LambdaV<0");
//.........这里部分代码省略.........

        /*************************************************************************
        This function sets linear term (model is a sum of radial  basis  functions
        plus linear polynomial). This function won't have effect until  next  call 
        to RBFBuildModel().

        INPUT PARAMETERS:
            S       -   RBF model, initialized by RBFCreate() call

        NOTE: this   function  has   some   serialization-related  subtleties.  We
              recommend you to study serialization examples from ALGLIB  Reference
              Manual if you want to perform serialization of your models.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfsetlinterm(rbfmodel s)
        {
            s.aterm = 1;
        }

        /*************************************************************************
        This function adds dataset.

        This function overrides results of the previous calls, i.e. multiple calls
        of this function will result in only the last set being added.

        INPUT PARAMETERS:
            S       -   RBF model, initialized by RBFCreate() call.
            XY      -   points, array[N,NX+NY]. One row corresponds to  one  point
                        in the dataset. First NX elements  are  coordinates,  next
                        NY elements are function values. Array may  be larger than 
                        specific,  in  this  case  only leading [N,NX+NY] elements 
                        will be used.
            N       -   number of points in the dataset

        After you've added dataset and (optionally) tuned algorithm  settings  you
        should call RBFBuildModel() in order to build a model for you.

        NOTE: this   function  has   some   serialization-related  subtleties.  We
              recommend you to study serialization examples from ALGLIB  Reference
              Manual if you want to perform serialization of your models.
              

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfsetpoints(rbfmodel s,
            double[,] xy,
            int n)
        {
            int i = 0;
            int j = 0;

            alglib.ap.assert(n>0, "RBFSetPoints: N<0");
            alglib.ap.assert(alglib.ap.rows(xy)>=n, "RBFSetPoints: Rows(XY)<N");
            alglib.ap.assert(alglib.ap.cols(xy)>=s.nx+s.ny, "RBFSetPoints: Cols(XY)<NX+NY");
            s.n = n;
            s.x = new double[s.n, mxnx];
            s.y = new double[s.n, s.ny];
            for(i=0; i<=s.n-1; i++)
            {
                for(j=0; j<=mxnx-1; j++)
                {
                    s.x[i,j] = 0;
                }
                for(j=0; j<=s.nx-1; j++)
                {
                    s.x[i,j] = xy[i,j];
                }
                for(j=0; j<=s.ny-1; j++)
                {
                    s.y[i,j] = xy[i,j+s.nx];
                }
            }
        }

        /*************************************************************************
        This  function  sets  RBF interpolation algorithm. ALGLIB supports several
        RBF algorithms with different properties.

        This algorithm is called RBF-QNN and  it  is  good  for  point  sets  with
        following properties:
        a) all points are distinct
        b) all points are well separated.
        c) points  distribution  is  approximately  uniform.  There is no "contour
           lines", clusters of points, or other small-scale structures.

        Algorithm description:
        1) interpolation centers are allocated to data points
        2) interpolation radii are calculated as distances to the  nearest centers
           times Q coefficient (where Q is a value from [0.75,1.50]).
        3) after  performing (2) radii are transformed in order to avoid situation
           when single outlier has very large radius and  influences  many  points
           across all dataset. Transformation has following form:
               new_r[i] = min(r[i],Z*median(r[]))
           where r[i] is I-th radius, median()  is a median  radius across  entire
           dataset, Z is user-specified value which controls amount  of  deviation
           from median radius.

        When (a) is violated,  we  will  be unable to build RBF model. When (b) or
        (c) are violated, model will be built, but interpolation quality  will  be
        low. See http://www.alglib.net/interpolation/ for more information on this
        subject.

        This algorithm is used by default.

        Additional Q parameter controls smoothness properties of the RBF basis:
        * Q<0.75 will give perfectly conditioned basis,  but  terrible  smoothness
          properties (RBF interpolant will have sharp peaks around function values)
        * Q around 1.0 gives good balance between smoothness and condition number
        * Q>1.5 will lead to badly conditioned systems and slow convergence of the
          underlying linear solver (although smoothness will be very good)
        * Q>2.0 will effectively make optimizer useless because it won't  converge
          within reasonable amount of iterations. It is possible to set such large
          Q, but it is advised not to do so.

        INPUT PARAMETERS:
            S       -   RBF model, initialized by RBFCreate() call
            Q       -   Q parameter, Q>0, recommended value - 1.0
            Z       -   Z parameter, Z>0, recommended value - 5.0

        NOTE: this   function  has   some   serialization-related  subtleties.  We
              recommend you to study serialization examples from ALGLIB  Reference
              Manual if you want to perform serialization of your models.


          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfsetalgoqnn(rbfmodel s,
            double q,
            double z)
        {
            alglib.ap.assert(math.isfinite(q), "RBFSetAlgoQNN: Q is infinite or NAN");
            alglib.ap.assert((double)(q)>(double)(0), "RBFSetAlgoQNN: Q<=0");
            rbfgridpoints(s);
            rbfradnn(s, q, z);
            s.algorithmtype = 1;
        }

        /*************************************************************************
        This function creates RBF  model  for  a  scalar (NY=1)  or  vector (NY>1)
        function in a NX-dimensional space (NX=2 or NX=3).

        Newly created model is empty. It can be used for interpolation right after
        creation, but it just returns zeros. You have to add points to the  model,
        tune interpolation settings, and then  call  model  construction  function
        RBFBuildModel() which will update model according to your specification.

        USAGE:
        1. User creates model with RBFCreate()
        2. User adds dataset with RBFSetPoints() (points do NOT have to  be  on  a
           regular grid)
        3. (OPTIONAL) User chooses polynomial term by calling:
           * RBFLinTerm() to set linear term
           * RBFConstTerm() to set constant term
           * RBFZeroTerm() to set zero term
           By default, linear term is used.
        4. User chooses specific RBF algorithm to use: either QNN (RBFSetAlgoQNN)
           or ML (RBFSetAlgoMultiLayer).
        5. User calls RBFBuildModel() function which rebuilds model  according  to
           the specification
        6. User may call RBFCalc() to calculate model value at the specified point,
           RBFGridCalc() to  calculate   model  values at the points of the regular
           grid. User may extract model coefficients with RBFUnpack() call.
           
        INPUT PARAMETERS:
            NX      -   dimension of the space, NX=2 or NX=3
            NY      -   function dimension, NY>=1

        OUTPUT PARAMETERS:
            S       -   RBF model (initially equals to zero)

        NOTE 1: memory requirements. RBF models require amount of memory  which is
                proportional  to  the  number  of data points. Memory is allocated 
                during model construction, but most of this memory is freed  after
                model coefficients are calculated.
                
                Some approximate estimates for N centers with default settings are
                given below:
                * about 250*N*(sizeof(double)+2*sizeof(int)) bytes  of  memory  is
                  needed during model construction stage.
                * about 15*N*sizeof(double) bytes is needed after model is built.
                For example, for N=100000 we may need 0.6 GB of memory  to  build
                model, but just about 0.012 GB to store it.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfcreate(int nx,
            int ny,
            rbfmodel s)
        {
            int i = 0;
            int j = 0;

            alglib.ap.assert(nx==2 || nx==3, "RBFCreate: NX<>2 and NX<>3");
            alglib.ap.assert(ny>=1, "RBFCreate: NY<1");
            s.nx = nx;
            s.ny = ny;
            s.nl = 0;
            s.nc = 0;
            s.v = new double[ny, mxnx+1];
            for(i=0; i<=ny-1; i++)
            {
                for(j=0; j<=mxnx; j++)
                {
                    s.v[i,j] = 0;
                }
            }
            s.n = 0;
            s.rmax = 0;
            s.gridtype = 2;
            s.fixrad = false;
            s.radvalue = 1;
            s.radzvalue = 5;
            s.aterm = 1;
            s.algorithmtype = 1;
            
            //
            // stopping criteria
            //
            s.epsort = eps;
            s.epserr = eps;
            s.maxits = 0;
        }

 public override alglib.apobject make_copy()
 {
     rbfmodel _result = new rbfmodel();
     _result.ny = ny;
     _result.nx = nx;
     _result.nc = nc;
     _result.nl = nl;
     _result.tree = (nearestneighbor.kdtree)tree.make_copy();
     _result.xc = (double[,])xc.Clone();
     _result.wr = (double[,])wr.Clone();
     _result.rmax = rmax;
     _result.v = (double[,])v.Clone();
     _result.gridtype = gridtype;
     _result.fixrad = fixrad;
     _result.lambdav = lambdav;
     _result.radvalue = radvalue;
     _result.radzvalue = radzvalue;
     _result.nlayers = nlayers;
     _result.aterm = aterm;
     _result.algorithmtype = algorithmtype;
     _result.epsort = epsort;
     _result.epserr = epserr;
     _result.maxits = maxits;
     _result.h = h;
     _result.n = n;
     _result.x = (double[,])x.Clone();
     _result.y = (double[,])y.Clone();
     _result.calcbufxcx = (double[])calcbufxcx.Clone();
     _result.calcbufx = (double[,])calcbufx.Clone();
     _result.calcbuftags = (int[])calcbuftags.Clone();
     return _result;
 }

        /*************************************************************************
        This function changes radii calculation algorithm to one which  makes  all
        radii equal to the same fixed value R. This  function  won't  have  effect
        until next call to RBFBuildModel().

        IMPORTANT: you should use this function with caution because too  large  R
        will make model fitting algorithm unstable, while too small  R  will  make
        perfect, but useless fit (it will be non-smooth, with sharp  peaks  around
        dataset points).

        INPUT PARAMETERS:
            S       -   RBF model, initialized by RBFCreate() call
            R       -   radius, R>0

        NOTE: this   function  has   some   serialization-related  subtleties.  We
              recommend you to study serialization examples from ALGLIB  Reference
              Manual if you want to perform serialization of your models.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        private static void rbfradfixed(rbfmodel s,
            double r)
        {
            alglib.ap.assert(math.isfinite(r) && (double)(r)>(double)(0), "RBFRadFixed: R<=0, infinite or NAN");
            s.fixrad = true;
            s.radvalue = r;
        }

        /*************************************************************************
        This function calculates values of the RBF model at the given point.

        Same as RBFCalc(), but does not reallocate Y when in is large enough to 
        store function values.

        INPUT PARAMETERS:
            S       -   RBF model
            X       -   coordinates, array[NX].
                        X may have more than NX elements, in this case only 
                        leading NX will be used.
            Y       -   possibly preallocated array

        OUTPUT PARAMETERS:
            Y       -   function value, array[NY]. Y is not reallocated when it
                        is larger than NY.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfcalcbuf(rbfmodel s,
            double[] x,
            ref double[] y)
        {
            int i = 0;
            int j = 0;
            int k = 0;
            int lx = 0;
            int tg = 0;
            double t = 0;
            double rcur = 0;
            double bf = 0;

            alglib.ap.assert(alglib.ap.len(x)>=s.nx, "RBFCalcBuf: Length(X)<NX");
            alglib.ap.assert(apserv.isfinitevector(x, s.nx), "RBFCalcBuf: X contains infinite or NaN values");
            if( alglib.ap.len(y)<s.ny )
            {
                y = new double[s.ny];
            }
            for(i=0; i<=s.ny-1; i++)
            {
                y[i] = s.v[i,mxnx];
                for(j=0; j<=s.nx-1; j++)
                {
                    y[i] = y[i]+s.v[i,j]*x[j];
                }
            }
            if( s.nc==0 )
            {
                return;
            }
            apserv.rvectorsetlengthatleast(ref s.calcbufxcx, mxnx);
            for(i=0; i<=mxnx-1; i++)
            {
                s.calcbufxcx[i] = 0.0;
            }
            for(i=0; i<=s.nx-1; i++)
            {
                s.calcbufxcx[i] = x[i];
            }
            lx = nearestneighbor.kdtreequeryrnn(s.tree, s.calcbufxcx, s.rmax*rbffarradius, true);
            nearestneighbor.kdtreequeryresultsx(s.tree, ref s.calcbufx);
            nearestneighbor.kdtreequeryresultstags(s.tree, ref s.calcbuftags);
            for(i=0; i<=s.ny-1; i++)
            {
                for(j=0; j<=lx-1; j++)
                {
                    tg = s.calcbuftags[j];
                    rcur = s.wr[tg,0];
                    bf = Math.Exp(-((math.sqr(s.calcbufxcx[0]-s.calcbufx[j,0])+math.sqr(s.calcbufxcx[1]-s.calcbufx[j,1])+math.sqr(s.calcbufxcx[2]-s.calcbufx[j,2]))/math.sqr(rcur)));
                    for(k=0; k<=s.nl-1; k++)
                    {
                        y[i] = y[i]+bf*s.wr[tg,1+k*s.ny+i];
                        t = bf*bf;
                        bf = t*t;
                    }
                }
            }
        }

        /*************************************************************************
        This function calculates values of the RBF model at the regular grid.

        Grid have N0*N1 points, with Point[I,J] = (X0[I], X1[J])

        This function returns 0.0 when:
        * model is not initialized
        * NX<>2
         *NY<>1

        INPUT PARAMETERS:
            S       -   RBF model
            X0      -   array of grid nodes, first coordinates, array[N0]
            N0      -   grid size (number of nodes) in the first dimension
            X1      -   array of grid nodes, second coordinates, array[N1]
            N1      -   grid size (number of nodes) in the second dimension

        OUTPUT PARAMETERS:
            Y       -   function values, array[N0,N1]. Y is out-variable and 
                        is reallocated by this function.

          -- ALGLIB --
             Copyright 13.12.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfgridcalc2(rbfmodel s,
            double[] x0,
            int n0,
            double[] x1,
            int n1,
            ref double[,] y)
        {
            double[] cpx0 = new double[0];
            double[] cpx1 = new double[0];
            int[] p01 = new int[0];
            int[] p11 = new int[0];
            int[] p2 = new int[0];
            double rlimit = 0;
            double xcnorm2 = 0;
            int hp01 = 0;
            double hcpx0 = 0;
            double xc0 = 0;
            double xc1 = 0;
            double omega = 0;
            double radius = 0;
            int i = 0;
            int j = 0;
            int k = 0;
            int d = 0;
            int i00 = 0;
            int i01 = 0;
            int i10 = 0;
            int i11 = 0;

            y = new double[0,0];

            alglib.ap.assert(n0>0, "RBFGridCalc2: invalid value for N0 (N0<=0)!");
            alglib.ap.assert(n1>0, "RBFGridCalc2: invalid value for N1 (N1<=0)!");
            alglib.ap.assert(alglib.ap.len(x0)>=n0, "RBFGridCalc2: Length(X0)<N0");
            alglib.ap.assert(alglib.ap.len(x1)>=n1, "RBFGridCalc2: Length(X1)<N1");
            alglib.ap.assert(apserv.isfinitevector(x0, n0), "RBFGridCalc2: X0 contains infinite or NaN values!");
            alglib.ap.assert(apserv.isfinitevector(x1, n1), "RBFGridCalc2: X1 contains infinite or NaN values!");
            y = new double[n0, n1];
            for(i=0; i<=n0-1; i++)
            {
                for(j=0; j<=n1-1; j++)
                {
                    y[i,j] = 0;
                }
            }
            if( (s.ny!=1 || s.nx!=2) || s.nc==0 )
            {
                return;
            }
            
            //
            //create and sort arrays
            //
            cpx0 = new double[n0];
            for(i=0; i<=n0-1; i++)
            {
                cpx0[i] = x0[i];
            }
            tsort.tagsort(ref cpx0, n0, ref p01, ref p2);
            cpx1 = new double[n1];
            for(i=0; i<=n1-1; i++)
            {
                cpx1[i] = x1[i];
            }
            tsort.tagsort(ref cpx1, n1, ref p11, ref p2);
            
            //
            //calculate function's value
            //
            for(i=0; i<=s.nc-1; i++)
            {
                radius = s.wr[i,0];
                for(d=0; d<=s.nl-1; d++)
                {
                    omega = s.wr[i,1+d];
                    rlimit = radius*rbffarradius;
//.........这里部分代码省略.........

        /*************************************************************************
        Serializer: allocation

          -- ALGLIB --
             Copyright 02.02.2012 by Bochkanov Sergey
        *************************************************************************/
        public static void rbfalloc(alglib.serializer s,
            rbfmodel model)
        {
            
            //
            // Header
            //
            s.alloc_entry();
            s.alloc_entry();
            
            //
            // Data
            //
            s.alloc_entry();
            s.alloc_entry();
            s.alloc_entry();
            s.alloc_entry();
            nearestneighbor.kdtreealloc(s, model.tree);
            apserv.allocrealmatrix(s, model.xc, -1, -1);
            apserv.allocrealmatrix(s, model.wr, -1, -1);
            s.alloc_entry();
            apserv.allocrealmatrix(s, model.v, -1, -1