pdtrsmrv.c

Go to the documentation of this file.

#include "common.h"

#define A(m,n) BLKADDR(A, double, m, n)
#define W(m)   BLKADDR(W, double, m, 0)

/***************************************************************************/
void plasma_pdtrsmrv_quark(PLASMA_enum side, PLASMA_enum uplo, PLASMA_enum trans, PLASMA_enum diag,
                           double alpha, PLASMA_desc A, PLASMA_desc W, 
                           PLASMA_sequence *sequence, PLASMA_request *request)
{
    plasma_context_t *plasma;
    Quark_Task_Flags task_flags = Quark_Task_Flags_Initializer;

    int k, m;
    int ldak, ldam, ldwk, ldwm;
    int tempm, tempkm, tempkn, tempmm;

    plasma = plasma_context_self();
    if (sequence->status != PLASMA_SUCCESS)
        return;
    QUARK_Task_Flag_Set(&task_flags, TASK_SEQUENCE, (intptr_t)sequence->quark_sequence);

    for (k = A.mt-1; k >= 0; k--) {
        tempm = (k+1)*A.mb;
        tempkm = k == A.mt-1 ? A.m-k*A.mb : A.mb;
        tempkn = k == A.nt-1 ? A.n-k*A.nb : A.nb;
        ldak = BLKLDD(A, k);
        ldwk = BLKLDD(W, k);
        
        QUARK_CORE_dlacpy(
            plasma->quark, &task_flags,
            PlasmaLower, tempkm, tempkn, A.mb,
            A(k, k), ldak, W(k), ldwk );

        QUARK_CORE_dlaset2(
            plasma->quark, &task_flags,
            PlasmaLower, tempkn, tempkn, 
            0.0, A(k, k), ldak );

        for(m=k+1; m<A.mt; m++) {
            tempmm = m == A.mt-1 ? A.m-m*A.mb : A.mb;
            ldam = BLKLDD(A, m);
            ldwm = BLKLDD(W, m);

            QUARK_CORE_dlacpy(
                plasma->quark, &task_flags,
                PlasmaUpperLower, tempmm, tempkn, A.mb,
                A(m, k), ldam, W(m), ldwm );
            
            QUARK_CORE_dlaset(
                plasma->quark, &task_flags,
                PlasmaUpperLower, tempmm, tempkn,
                0.0, 0.0, A(m, k), ldam );
        }            

        if (k*A.mb+tempkn < A.m) {
           plasma_pdgemm_quark(
               PlasmaNoTrans, PlasmaNoTrans, -alpha,
               plasma_desc_submatrix(A, 0,     tempm,  tempm,     A.n-tempm),
               plasma_desc_submatrix(W, tempm, 0,      W.m-tempm, tempkn   ),
               alpha,
               plasma_desc_submatrix(A, 0,     k*A.mb, tempm,     tempkn   ),
               sequence, request);

           plasma_pdgemm_quark(
               PlasmaNoTrans, PlasmaNoTrans, -alpha,
               plasma_desc_submatrix(A, tempm, tempm,  A.m-tempm, A.n-tempm),
               plasma_desc_submatrix(W, tempm, 0,      W.m-tempm, tempkn   ),
               (double)0.,
               plasma_desc_submatrix(A, tempm, k*A.mb, A.m-tempm, tempkn   ),
               sequence, request);
        }

        plasma_pdtrsm_quark(
            PlasmaRight, PlasmaLower, 
            PlasmaNoTrans, PlasmaUnit, 
            alpha,
            plasma_desc_submatrix(W, k*W.mb, 0,      tempkm, tempkn),
            plasma_desc_submatrix(A, 0,      k*A.mb, A.m,    tempkn),
            sequence, request);
    }
}