C++ real_1d_array::setbounds方法代码示例

static void fillsparsede(ap::real_1d_array& d,
     ap::real_1d_array& e,
     int n,
     double sparcity)
{
    int i;
    int j;

    d.setbounds(0, n-1);
    e.setbounds(0, ap::maxint(0, n-2));
    for(i = 0; i <= n-1; i++)
    {
        if( ap::fp_greater_eq(ap::randomreal(),sparcity) )
        {
            d(i) = 2*ap::randomreal()-1;
        }
        else
        {
            d(i) = 0;
        }
    }
    for(i = 0; i <= n-2; i++)
    {
        if( ap::fp_greater_eq(ap::randomreal(),sparcity) )
        {
            e(i) = 2*ap::randomreal()-1;
        }
        else
        {
            e(i) = 0;
        }
    }
}

bool in_out_variable_1D(const ap::boolean_1d_array& in, const ap::real_1d_array& X, ap::real_1d_array& x, ap::real_1d_array& vector,  bool io)
{
// Routine to know the number of variables in/out
	int rows = in.gethighbound(0) + 1;
	int n_invar=0;
	bool flag;
	
	//ap::real_1d_array vector;
	vector.setbounds(0,rows-1);
	unsigned int k=0;
	for (int i=0; i<rows; i++)
	{
		if (in(i)==io) //to know how many variables are in/out
		{
			vector(k) = i;
			k++;
			n_invar++;
		}
	}
	if (n_invar>0)
	{
		// Routine to extract the in/out variables
		x.setbounds(0,n_invar-1);
		for (int i=0; i<n_invar; i++)
		x(i) = X(static_cast<int>(vector(i)));
	
        flag=TRUE;
	}
	else
		flag=FALSE;

	return flag;
}

/*************************************************************************
Serialization of LinearModel strucure

INPUT PARAMETERS:
    LM      -   original

OUTPUT PARAMETERS:
    RA      -   array of real numbers which stores model,
                array[0..RLen-1]
    RLen    -   RA lenght

  -- ALGLIB --
     Copyright 15.03.2009 by Bochkanov Sergey
*************************************************************************/
void lrserialize(const linearmodel& lm, ap::real_1d_array& ra, int& rlen)
{

    rlen = ap::round(lm.w(0))+1;
    ra.setbounds(0, rlen-1);
    ra(0) = lrvnum;
    ap::vmove(&ra(1), &lm.w(0), ap::vlen(1,rlen-1));
}

/*************************************************************************
Solving a system of linear equations with a system matrix given by its
LU decomposition.

The algorithm solves a system of linear equations whose matrix is given by
its LU decomposition. In case of a singular matrix, the algorithm  returns
False.

The algorithm solves systems with a square matrix only.

Input parameters:
    A       -   LU decomposition of a system matrix in compact  form  (the
                result of the RMatrixLU subroutine).
    Pivots  -   row permutation table (the result of a
                RMatrixLU subroutine).
    B       -   right side of a system.
                Array whose index ranges within [0..N-1].
    N       -   size of matrix A.

Output parameters:
    X       -   solution of a system.
                Array whose index ranges within [0..N-1].

Result:
    True, if the matrix is not singular.
    False, if the matrux is singular. In this case, X doesn't contain a
solution.

  -- ALGLIB --
     Copyright 2005-2008 by Bochkanov Sergey
*************************************************************************/
bool rmatrixlusolve(const ap::real_2d_array& a,
     const ap::integer_1d_array& pivots,
     ap::real_1d_array b,
     int n,
     ap::real_1d_array& x)
{
    bool result;
    ap::real_1d_array y;
    int i;
    int j;
    double v;

    y.setbounds(0, n-1);
    x.setbounds(0, n-1);
    result = true;
    for(i = 0; i <= n-1; i++)
    {
        if( a(i,i)==0 )
        {
            result = false;
            return result;
        }
    }
    
    //
    // pivots
    //
    for(i = 0; i <= n-1; i++)
    {
        if( pivots(i)!=i )
        {
            v = b(i);
            b(i) = b(pivots(i));
            b(pivots(i)) = v;
        }
    }
    
    //
    // Ly = b
    //
    y(0) = b(0);
    for(i = 1; i <= n-1; i++)
    {
        v = ap::vdotproduct(&a(i, 0), &y(0), ap::vlen(0,i-1));
        y(i) = b(i)-v;
    }
    
    //
    // Ux = y
    //
    x(n-1) = y(n-1)/a(n-1,n-1);
    for(i = n-2; i >= 0; i--)
    {
        v = ap::vdotproduct(&a(i, i+1), &x(i+1), ap::vlen(i+1,n-1));
        x(i) = (y(i)-v)/a(i,i);
    }
    return result;
}

/*************************************************************************
Unpacks coefficients of linear model.

INPUT PARAMETERS:
    LM          -   linear model in ALGLIB format

OUTPUT PARAMETERS:
    V           -   coefficients, array[0..NVars]
    NVars       -   number of independent variables (one less than number
                    of coefficients)

  -- ALGLIB --
     Copyright 30.08.2008 by Bochkanov Sergey
*************************************************************************/
void lrunpack(const linearmodel& lm, ap::real_1d_array& v, int& nvars)
{
    int offs;

    ap::ap_error::make_assertion(ap::round(lm.w(1))==lrvnum, "LINREG: Incorrect LINREG version!");
    nvars = ap::round(lm.w(2));
    offs = ap::round(lm.w(3));
    v.setbounds(0, nvars);
    ap::vmove(&v(0), &lm.w(offs), ap::vlen(0,nvars));
}

/*************************************************************************
Conversion of a series of Chebyshev polynomials to a power series.

Represents A[0]*T0(x) + A[1]*T1(x) + ... + A[N]*Tn(x) as
B[0] + B[1]*X + ... + B[N]*X^N.

Input parameters:
    A   -   Chebyshev series coefficients
    N   -   degree, N>=0
    
Output parameters
    B   -   power series coefficients
*************************************************************************/
void fromchebyshev(const ap::real_1d_array& a,
     const int& n,
     ap::real_1d_array& b)
{
    int i;
    int k;
    double e;
    double d;

    b.setbounds(0, n);
    for(i = 0; i <= n; i++)
    {
        b(i) = 0;
    }
    d = 0;
    i = 0;
    do
    {
        k = i;
        do
        {
            e = b(k);
            b(k) = 0;
            if( i<=1&&k==i )
            {
                b(k) = 1;
            }
            else
            {
                if( i!=0 )
                {
                    b(k) = 2*d;
                }
                if( k>i+1 )
                {
                    b(k) = b(k)-b(k-2);
                }
            }
            d = e;
            k = k+1;
        }
        while(k<=n);
        d = b(i);
        e = 0;
        k = i;
        while(k<=n)
        {
            e = e+b(k)*a(k);
            k = k+2;
        }
        b(i) = e;
        i = i+1;
    }
    while(i<=n);
}

/*************************************************************************
Representation of Ln as C[0] + C[1]*X + ... + C[N]*X^N

Input parameters:
    N   -   polynomial degree, n>=0

Output parameters:
    C   -   coefficients
*************************************************************************/
void laguerrecoefficients(const int& n, ap::real_1d_array& c)
{
    int i;

    c.setbounds(0, n);
    c(0) = 1;
    for(i = 0; i <= n-1; i++)
    {
        c(i+1) = -c(i)*(n-i)/(i+1)/(i+1);
    }
}

/*************************************************************************
L-BFGS algorithm results

Called after MinLBFGSIteration() returned False.

INPUT PARAMETERS:
    State   -   algorithm state (used by MinLBFGSIteration).

OUTPUT PARAMETERS:
    X       -   array[0..N-1], solution
    Rep     -   optimization report:
                * Rep.TerminationType completetion code:
                    * -2    rounding errors prevent further improvement.
                            X contains best point found.
                    * -1    incorrect parameters were specified
                    *  1    relative function improvement is no more than
                            EpsF.
                    *  2    relative step is no more than EpsX.
                    *  4    gradient norm is no more than EpsG
                    *  5    MaxIts steps was taken
                    *  7    stopping conditions are too stringent,
                            further improvement is impossible
                * Rep.IterationsCount contains iterations count
                * NFEV countains number of function calculations

  -- ALGLIB --
     Copyright 02.04.2010 by Bochkanov Sergey
*************************************************************************/
void minlbfgsresults(const minlbfgsstate& state,
     ap::real_1d_array& x,
     minlbfgsreport& rep)
{

    x.setbounds(0, state.n-1);
    ap::vmove(&x(0), 1, &state.x(0), 1, ap::vlen(0,state.n-1));
    rep.iterationscount = state.repiterationscount;
    rep.nfev = state.repnfev;
    rep.terminationtype = state.repterminationtype;
}

/*************************************************************************
Serialization of DecisionForest strucure

INPUT PARAMETERS:
    DF      -   original

OUTPUT PARAMETERS:
    RA      -   array of real numbers which stores decision forest,
                array[0..RLen-1]
    RLen    -   RA lenght

  -- ALGLIB --
     Copyright 13.02.2009 by Bochkanov Sergey
*************************************************************************/
void dfserialize(const decisionforest& df, ap::real_1d_array& ra, int& rlen)
{

    ra.setbounds(0, df.bufsize+5-1);
    ra(0) = dfvnum;
    ra(1) = df.nvars;
    ra(2) = df.nclasses;
    ra(3) = df.ntrees;
    ra(4) = df.bufsize;
    ap::vmove(&ra(5), 1, &df.trees(0), 1, ap::vlen(5,5+df.bufsize-1));
    rlen = 5+df.bufsize;
}

/*************************************************************************
Unpacking of the main and secondary diagonals of bidiagonal decomposition
of matrix A.

Input parameters:
    B   -   output of RMatrixBD subroutine.
    M   -   number of rows in matrix B.
    N   -   number of columns in matrix B.

Output parameters:
    IsUpper -   True, if the matrix is upper bidiagonal.
                otherwise IsUpper is False.
    D       -   the main diagonal.
                Array whose index ranges within [0..Min(M,N)-1].
    E       -   the secondary diagonal (upper or lower, depending on
                the value of IsUpper).
                Array index ranges within [0..Min(M,N)-1], the last
                element is not used.

  -- ALGLIB --
     Copyright 2005-2007 by Bochkanov Sergey
*************************************************************************/
void rmatrixbdunpackdiagonals(const ap::real_2d_array& b,
     int m,
     int n,
     bool& isupper,
     ap::real_1d_array& d,
     ap::real_1d_array& e)
{
    int i;

    isupper = m>=n;
    if( m<=0||n<=0 )
    {
        return;
    }
    if( isupper )
    {
        d.setbounds(0, n-1);
        e.setbounds(0, n-1);
        for(i = 0; i <= n-2; i++)
        {
            d(i) = b(i,i);
            e(i) = b(i,i+1);
        }
        d(n-1) = b(n-1,n-1);
    }
    else
    {
        d.setbounds(0, m-1);
        e.setbounds(0, m-1);
        for(i = 0; i <= m-2; i++)
        {
            d(i) = b(i,i);
            e(i) = b(i+1,i);
        }
        d(m-1) = b(m-1,m-1);
    }
}

/*************************************************************************
Obsolete 1-based subroutine.
See RMatrixBDUnpackDiagonals for 0-based replacement.
*************************************************************************/
void unpackdiagonalsfrombidiagonal(const ap::real_2d_array& b,
     int m,
     int n,
     bool& isupper,
     ap::real_1d_array& d,
     ap::real_1d_array& e)
{
    int i;

    isupper = m>=n;
    if( m==0||n==0 )
    {
        return;
    }
    if( isupper )
    {
        d.setbounds(1, n);
        e.setbounds(1, n);
        for(i = 1; i <= n-1; i++)
        {
            d(i) = b(i,i);
            e(i) = b(i,i+1);
        }
        d(n) = b(n,n);
    }
    else
    {
        d.setbounds(1, m);
        e.setbounds(1, m);
        for(i = 1; i <= m-1; i++)
        {
            d(i) = b(i,i);
            e(i) = b(i+1,i);
        }
        d(m) = b(m,m);
    }
}

/*************************************************************************
Representation of Hn as C[0] + C[1]*X + ... + C[N]*X^N

Input parameters:
    N   -   polynomial degree, n>=0

Output parameters:
    C   -   coefficients
*************************************************************************/
void hermitecoefficients(const int& n, ap::real_1d_array& c)
{
    int i;

    c.setbounds(0, n);
    for(i = 0; i <= n; i++)
    {
        c(i) = 0;
    }
    c(n) = exp(n*log(double(2)));
    for(i = 0; i <= n/2-1; i++)
    {
        c(n-2*(i+1)) = -c(n-2*i)*(n-2*i)*(n-2*i-1)/4/(i+1);
    }
}

/*************************************************************************
Multiclass Fisher LDA

Subroutine finds coefficients of linear combination which optimally separates
training set on classes.

INPUT PARAMETERS:
    XY          -   training set, array[0..NPoints-1,0..NVars].
                    First NVars columns store values of independent
                    variables, next column stores number of class (from 0
                    to NClasses-1) which dataset element belongs to. Fractional
                    values are rounded to nearest integer.
    NPoints     -   training set size, NPoints>=0
    NVars       -   number of independent variables, NVars>=1
    NClasses    -   number of classes, NClasses>=2


OUTPUT PARAMETERS:
    Info        -   return code:
                    * -4, if internal EVD subroutine hasn't converged
                    * -2, if there is a point with class number
                          outside of [0..NClasses-1].
                    * -1, if incorrect parameters was passed (NPoints<0,
                          NVars<1, NClasses<2)
                    *  1, if task has been solved
                    *  2, if there was a multicollinearity in training set,
                          but task has been solved.
    W           -   linear combination coefficients, array[0..NVars-1]

  -- ALGLIB --
     Copyright 31.05.2008 by Bochkanov Sergey
*************************************************************************/
void fisherlda(const ap::real_2d_array& xy,
     int npoints,
     int nvars,
     int nclasses,
     int& info,
     ap::real_1d_array& w)
{
    ap::real_2d_array w2;

    fisherldan(xy, npoints, nvars, nclasses, info, w2);
    if( info>0 )
    {
        w.setbounds(0, nvars-1);
        ap::vmove(&w(0), 1, &w2(0, 0), w2.getstride(), ap::vlen(0,nvars-1));
    }
}

void qrdecomposition(ap::real_2d_array& a,
     int m,
     int n,
     ap::real_1d_array& tau)
{
    ap::real_1d_array work;
    ap::real_1d_array t;
    int i;
    int k;
    int mmip1;
    int minmn;
    double tmp;

    minmn = ap::minint(m, n);
    work.setbounds(1, n);
    t.setbounds(1, m);
    tau.setbounds(1, minmn);
    
    //
    // Test the input arguments
    //
    k = ap::minint(m, n);
    for(i = 1; i <= k; i++)
    {
        
        //
        // Generate elementary reflector H(i) to annihilate A(i+1:m,i)
        //
        mmip1 = m-i+1;
        ap::vmove(t.getvector(1, mmip1), a.getcolumn(i, i, m));
        generatereflection(t, mmip1, tmp);
        tau(i) = tmp;
        ap::vmove(a.getcolumn(i, i, m), t.getvector(1, mmip1));
        t(1) = 1;
        if( i<n )
        {
            
            //
            // Apply H(i) to A(i:m,i+1:n) from the left
            //
            applyreflectionfromtheleft(a, tau(i), t, i, m, i+1, n, work);
        }
    }
}

/*************************************************************************
QR decomposition of a rectangular matrix of size MxN

Input parameters:
    A   -   matrix A whose indexes range within [0..M-1, 0..N-1].
    M   -   number of rows in matrix A.
    N   -   number of columns in matrix A.

Output parameters:
    A   -   matrices Q and R in compact form (see below).
    Tau -   array of scalar factors which are used to form
            matrix Q. Array whose index ranges within [0.. Min(M-1,N-1)].

Matrix A is represented as A = QR, where Q is an orthogonal matrix of size
MxM, R - upper triangular (or upper trapezoid) matrix of size M x N.

The elements of matrix R are located on and above the main diagonal of
matrix A. The elements which are located in Tau array and below the main
diagonal of matrix A are used to form matrix Q as follows:

Matrix Q is represented as a product of elementary reflections

Q = H(0)*H(2)*...*H(k-1),

where k = min(m,n), and each H(i) is in the form

H(i) = 1 - tau * v * (v^T)

where tau is a scalar stored in Tau[I]; v - real vector,
so that v(0:i-1) = 0, v(i) = 1, v(i+1:m-1) stored in A(i+1:m-1,i).

  -- LAPACK routine (version 3.0) --
     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
     Courant Institute, Argonne National Lab, and Rice University
     February 29, 1992.
     Translation from FORTRAN to pseudocode (AlgoPascal)
     by Sergey Bochkanov, ALGLIB project, 2005-2007.
*************************************************************************/
void rmatrixqr(ap::real_2d_array& a, int m, int n, ap::real_1d_array& tau)
{
    ap::real_1d_array work;
    ap::real_1d_array t;
    int i;
    int k;
    int minmn;
    double tmp;

    if( m<=0||n<=0 )
    {
        return;
    }
    minmn = ap::minint(m, n);
    work.setbounds(0, n-1);
    t.setbounds(1, m);
    tau.setbounds(0, minmn-1);
    
    //
    // Test the input arguments
    //
    k = minmn;
    for(i = 0; i <= k-1; i++)
    {
        
        //
        // Generate elementary reflector H(i) to annihilate A(i+1:m,i)
        //
        ap::vmove(t.getvector(1, m-i), a.getcolumn(i, i, m-1));
        generatereflection(t, m-i, tmp);
        tau(i) = tmp;
        ap::vmove(a.getcolumn(i, i, m-1), t.getvector(1, m-i));
        t(1) = 1;
        if( i<n )
        {
            
            //
            // Apply H(i) to A(i:m-1,i+1:n-1) from the left
            //
            applyreflectionfromtheleft(a, tau(i), t, i, m-1, i+1, n-1, work);
        }
    }
}