C++ real_2d_array::setlength方法代码示例

/*************************************************************************
Generate matrix with given condition number C (2-norm)
*************************************************************************/
static void rmatrixgenzero(ap::real_2d_array& a0, int n)
{
    int i;
    int j;

    a0.setlength(n, n);
    for(i = 0; i <= n-1; i++)
    {
        for(j = 0; j <= n-1; j++)
        {
            a0(i,j) = 0;
        }
    }
}

static void mheapresize(ap::real_2d_array& heap,
     int& heapsize,
     int newheapsize,
     int heapwidth)
{
    ap::real_2d_array tmp;
    int i;

    tmp.setlength(heapsize, heapwidth);
    for(i = 0; i <= heapsize-1; i++)
    {
        ap::vmove(&tmp(i, 0), &heap(i, 0), ap::vlen(0,heapwidth-1));
    }
    heap.setlength(newheapsize, heapwidth);
    for(i = 0; i <= heapsize-1; i++)
    {
        ap::vmove(&heap(i, 0), &tmp(i, 0), ap::vlen(0,heapwidth-1));
    }
    heapsize = newheapsize;
}

/*************************************************************************
Dense solver.

This  subroutine  solves  a  system  A*X=B,  where A is NxN non-denegerate
real matrix, X and B are NxM real matrices.

Additional features include:
* automatic detection of degenerate cases
* iterative improvement

INPUT PARAMETERS
    A       -   array[0..N-1,0..N-1], system matrix
    N       -   size of A
    B       -   array[0..N-1,0..M-1], right part
    M       -   size of right part
    
OUTPUT PARAMETERS
    Info    -   return code:
                * -3    if A is singular, or VERY close to singular.
                        X is filled by zeros in such cases.
                * -1    if N<=0 or M<=0 was passed
                *  1    if task is solved (matrix A may be near  singular,
                        check R1/RInf parameters for condition numbers).
    Rep     -   solver report, see below for more info
    X       -   array[0..N-1,0..M-1], it contains:
                * solution of A*X=B if A is non-singular (well-conditioned
                  or ill-conditioned, but not very close to singular)
                * zeros,  if  A  is  singular  or  VERY  close to singular
                  (in this case Info=-3).

SOLVER REPORT

Subroutine sets following fields of the Rep structure:
* R1        reciprocal of condition number: 1/cond(A), 1-norm.
* RInf      reciprocal of condition number: 1/cond(A), inf-norm.

SEE ALSO:
    DenseSolverR() - solves A*x = b, where x and b are Nx1 matrices.

  -- ALGLIB --
     Copyright 24.08.2009 by Bochkanov Sergey
*************************************************************************/
void rmatrixsolvem(const ap::real_2d_array& a,
     int n,
     const ap::real_2d_array& b,
     int m,
     int& info,
     densesolverreport& rep,
     ap::real_2d_array& x)
{
    int i;
    int j;
    int k;
    int rfs;
    int nrfs;
    ap::integer_1d_array p;
    ap::real_1d_array xc;
    ap::real_1d_array y;
    ap::real_1d_array bc;
    ap::real_1d_array xa;
    ap::real_1d_array xb;
    ap::real_1d_array tx;
    ap::real_2d_array da;
    double v;
    double verr;
    bool smallerr;
    bool terminatenexttime;

    
    //
    // prepare: check inputs, allocate space...
    //
    if( n<=0||m<=0 )
    {
        info = -1;
        return;
    }
    da.setlength(n, n);
    x.setlength(n, m);
    y.setlength(n);
    xc.setlength(n);
    bc.setlength(n);
    tx.setlength(n+1);
    xa.setlength(n+1);
    xb.setlength(n+1);
    
    //
    // factorize matrix, test for exact/near singularity
    //
    for(i = 0; i <= n-1; i++)
    {
        ap::vmove(&da(i, 0), &a(i, 0), ap::vlen(0,n-1));
    }
    rmatrixlu(da, n, n, p);
    rep.r1 = rmatrixlurcond1(da, n);
    rep.rinf = rmatrixlurcondinf(da, n);
    if( ap::fp_less(rep.r1,10*ap::machineepsilon)||ap::fp_less(rep.rinf,10*ap::machineepsilon) )
    {
        for(i = 0; i <= n-1; i++)
        {
//.........这里部分代码省略.........