LUBackSubst

Solves system of linear equations with LU-factored square matrix.

Syntax

Case 1: Matrix - vector operation

IppStatus ippmLUBackSubst_mv_32f(const Ipp32f* pSrc1, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp32f* pSrc2, int src2Stride2, Ipp32f* pDst, int dstStride2, int widthHeight);

IppStatus ippmLUBackSubst_mv_64f(const Ipp64f* pSrc1, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp64f* pSrc2, int src2Stride2, Ipp64f* pDst, int dstStride2, int widthHeight);

IppStatus ippmLUBackSubst_mv_32f_P(const Ipp32f** ppSrc1, int src1RoiShift, int* pSrcIndex, const Ipp32f** ppSrc2, int src2RoiShift, Ipp32f** ppDst, int dstRoiShift, int widthHeight);

IppStatus ippmLUBackSubst_mv_64f_P(const Ipp64f** ppSrc1, int src1RoiShift, int* pSrcIndex, const Ipp64f** ppSrc2, int src2RoiShift, Ipp64f** ppDst, int dstRoiShift, int widthHeight);

Case 2: Matrix - vector array operation

IppStatus ippmLUBackSubst_mva_32f(const Ipp32f* pSrc1, int src1Stride0, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp32f* pSrc2, int src2Stride0, int src2Stride2, Ipp32f* pDst, int dstStride0, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mva_64f(const Ipp64f* pSrc1, int src1Stride0, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp64f* pSrc2, int src2Stride0, int src2Stride2, Ipp64f* pDst, int dstStride0, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mva_32f_P(const Ipp32f** ppSrc1, int src1RoiShift, int* pSrcIndex, const Ipp32f** ppSrc2, int src2RoiShift, int src2Stride0, Ipp32f** ppDst, int dstRoiShift, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mva_64f_P(const Ipp64f** ppSrc1, int src1RoiShift, int* pSrcIndex, const Ipp64f** ppSrc2, int src2RoiShift, int src2Stride0, Ipp64f** ppDst, int dstRoiShift, int dstStride0, int widthHeight, int count);

IppStatus ippmLUBackSubst_mva_32f_L(const Ipp32f** ppSrc1, int src1RoiShift, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp32f** ppSrc2, int src2RoiShift, int src2Stride2, Ipp32f** ppDst, int dstRoiShift, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mva_64f_L(const Ipp64f** ppSrc1, int src1RoiShift, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp64f** ppSrc2, int src2RoiShift, int src2Stride2, Ipp64f** ppDst, int dstRoiShift, int dstStride2, int widthHeight, int count);

Case 3: Matrix array - vector array operation

IppStatus ippmLUBackSubst_mava_32f(const Ipp32f* pSrc1, int src1Stride0, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp32f* pSrc2, int src2Stride0, int src2Stride2, Ipp32f* pDst, int dstStride0, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mava_64f(const Ipp64f* pSrc1, int src1Stride0, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp64f* pSrc2, int src2Stride0, int src2Stride2, Ipp64f* pDst, int dstStride0, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mava_32f_P(const Ipp32f** ppSrc1, int src1RoiShift, int src1Stride0, int* pSrcIndex, const Ipp32f** ppSrc2, int src2RoiShift, int src2Stride0, Ipp32f** ppDst, int dstRoiShift, int dstStride0, int widthHeight, int count);

IppStatus ippmLUBackSubst_mava_64f_P(const Ipp64f** ppSrc1, int src1RoiShift, int src1Stride0, int* pSrcIndex, const Ipp64f** ppSrc2, int src2RoiShift, int src2Stride0, Ipp64f** ppDst, int dstRoiShift, int dstStride0, int widthHeight, int count);

IppStatus ippmLUBackSubst_mava_32f_L(const Ipp32f** ppSrc1, int src1RoiShift, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp32f** ppSrc2, int src2RoiShift, int src2Stride2, Ipp32f** ppDst, int dstRoiShift, int dstStride2, int widthHeight, int count);

IppStatus ippmLUBackSubst_mava_64f_L(const Ipp64f** ppSrc1, int src1RoiShift, int src1Stride1, int src1Stride2, int* pSrcIndex, const Ipp64f** ppSrc2, int src2RoiShift, int src2Stride2, Ipp64f** ppDst, int dstRoiShift, int dstStride2, int widthHeight, int count);

Parameters

pSrc1, ppSrc1

Pointer to the source matrix or array of matrices. Matrix(ces) must be a result of calling LUDecomp.

src1Stride0

Stride between the matrices in the source matrix array.

src1Stride1

Stride between the rows in the source matrix(ces).

src1Stride2

Stride between the elements in the source matrix(ces).

src1RoiShift

ROI shift in the source matrix(ces).

pSrc2, ppSrc2

Pointer to the source vector or array of vectors.

src2Stride0

Stride between the vectors in the source vector array.

src2Stride2

Stride between the elements in the source vector(s).

src2RoiShift

ROI shift in the source vector(s).

pSrcIndex

Pointer to array of pivot indices. This array must be a result of calling LUDecomp. The array size can be more than or equal to widthHeight. If the operation is performed on an array of matrices, the size of the array of indices must be more than or equal to count *widthHeight.

pDst, ppDst

Pointer to the destination vector or array of vectors.

dstStride0

Stride between the vectors in the destination array.

dstStride2

Stride between the elements in the destination vector(s).

dstRoiShift

ROI shift in the destination vector(s).

widthHeight

Size of the square matrix, the source vector, and the destination vector.

count

Number of matrices and right-hand part vectors in the arrays.

Description

The function ippmLUBackSubst is declared in the ippm.h header file. The function solves for x the following systems of linear equations:

where A is the matrix of linear equations system, stored in pSrc1 or ppSrc1, b is the vector of the right-hand side, stored in pSrc2 or ppSrc2, x is the unknown vector, stored in pDst or ppDst.

You should call the function ippmLUDecomp to compute the LU decomposition of A before calling ippmLUBackSubst.

The following example demonstrates how to use the functions ippmLUDecomp_m_32f and ippmLUBackSubst_mva_32f. For more information, see also examples in Getting Started.

ippmLUFactorization_32f  

IppStatus LUFactorization_32f(void){
    /* Source matrix with widthHeight=3 */
    Ipp32f pSrc[3*3] = { 3, -5, -10,
                        -1,  4,  2 ,
                         1, -2, -3 };
    int srcStride2 = sizeof(Ipp32f);
    int srcStride1 = 3*sizeof(Ipp32f);

    /* Solver right-part is 2 vectors with length=3 */
    Ipp32f pSrc2[3*2] = { 0,  2, 1,
                          1, -1, 2 };
    int src2Stride2 = sizeof(Ipp32f);
    int src2Stride0 = 3*sizeof(Ipp32f);

    Ipp32f pDecomp[3*3]; /* Decomposed matrix location */
    int decompStride2 = sizeof(Ipp32f);
    int decompStride1 = 3*sizeof(Ipp32f);

    Ipp32f pDst[3*2];    /* Solver destination location */
    int dstStride2 = sizeof(Ipp32f);
    int dstStride0 = 3*sizeof(Ipp32f);
    int pIndex[3];       /* Pivoting indeces location */
    int widthHeight = 3;
    int count = 2;

    IppStatus status = ippmLUDecomp_m_32f((const Ipp32f*)pSrc,
        srcStride1, srcStride2, pIndex,
        pDecomp, decompStride1, decompStride2, widthHeight);

    /*
    // It is required for LU decomposition function to check return status
    // for catching wrong result in case of invalid input data 
    */
    if (status == ippStsNoErr) {
        status = ippmLUBackSubst_mva_32f((const Ipp32f*)pDecomp, 
            decompStride1, decompStride2, pIndex, pSrc2, src2Stride0,
            src2Stride2, pDst, dstStride0, dstStride2, widthHeight, count);
 
        printf_m_Ipp32f("LUDecomp result:", pDecomp, 3, 3, status);
        printf_v_int("Pivoting indices:", pIndex, 3);
        printf_va_Ipp32f("2 destination vectors:", pDst, 3, 2, status);
 
    } else {
        printf("Function returns status: %s \n", ippGetStatusString(status));
    }
    return status;
}
 

The program above produces the following output:

LUDecomp result:

3.000000  -5.000000  -10.000000

-0.333333  2.333333  -1.333333

0.333333  -0.142857  0.142857

Pivoting indices:

0  1  2

2 destination vectors:

40.000000  6.000000  9.000000

47.000000  6.000000  11.000000

Return Values

ippStsOk

Returns no error.

ippStsNullPtrErr

Returns an error when at least one input pointer is NULL.

ippStsSizeErr

Returns an error when the size of the source matrix is equal to 0.

ippStsStrideMatrixErr

Returns an error when the stride value is not positive or not divisible by size of data type.

ippStsRoiShiftMatrixErr

Returns an error when the roiShift value is negative or not divisible by size of data type.

ippStsCountMatrixErr

Returns an error when the count value is less or equal to zero.

Submit feedback on this help topic

Copyright © 2000 - 2010, Intel Corporation. All rights reserved.