testing_zsyrk.c 6.67 KB
Newer Older
1
/**
2 3
 *
 * @file testing_zsyrk.c
4
 *
Mathieu Faverge's avatar
Mathieu Faverge committed
5 6
 * @copyright 2009-2014 The University of Tennessee and The University of
 *                      Tennessee Research Foundation. All rights reserved.
Mathieu Faverge's avatar
Mathieu Faverge committed
7
 * @copyright 2012-2018 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8
 *                      Univ. Bordeaux. All rights reserved.
9
 *
10
 ***
11
 *
12
 * @brief Chameleon zsyrk testing
13
 *
Mathieu Faverge's avatar
Mathieu Faverge committed
14
 * @version 1.0.0
15
 * @comment This file has been automatically generated
16
 *          from Plasma 2.5.0 for CHAMELEON 1.0.0
17 18 19 20 21 22
 * @author Mathieu Faverge
 * @author Emmanuel Agullo
 * @author Cedric Castagnede
 * @date 2010-11-15
 * @precisions normal z -> c d s
 *
23
 */
24 25 26 27 28
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>

29
#include <chameleon.h>
30 31 32
#include <coreblas/cblas.h>
#include <coreblas/lapacke.h>
#include <coreblas.h>
33 34
#include "testing_zauxiliary.h"

35 36
static int check_solution(cham_uplo_t uplo, cham_trans_t trans, int N, int K,
                          CHAMELEON_Complex64_t alpha, CHAMELEON_Complex64_t *A, int LDA,
Mathieu Faverge's avatar
Mathieu Faverge committed
37
                          CHAMELEON_Complex64_t beta,  CHAMELEON_Complex64_t *Cref, CHAMELEON_Complex64_t *Ccham, int LDC);
38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54


int testing_zsyrk(int argc, char **argv)
{
    int hres = 0;
    /* Check for number of arguments*/
    if ( argc != 6){
        USAGE("SYRK", "alpha beta M N LDA LDC",
              "   - alpha : alpha coefficient\n"
              "   - beta : beta coefficient\n"
              "   - N : number of columns and rows of matrix C and number of row of matrix A\n"
              "   - K : number of columns of matrix A\n"
              "   - LDA : leading dimension of matrix A\n"
              "   - LDC : leading dimension of matrix C\n");
        return -1;
    }

55 56
    CHAMELEON_Complex64_t alpha = (CHAMELEON_Complex64_t) atol(argv[0]);
    CHAMELEON_Complex64_t beta  = (CHAMELEON_Complex64_t) atol(argv[1]);
57 58 59 60 61 62 63 64 65 66 67 68
    int N     = atoi(argv[2]);
    int K     = atoi(argv[3]);
    int LDA   = atoi(argv[4]);
    int LDC   = atoi(argv[5]);
    int NKmax = max(N, K);

    double eps;
    int info_solution;
    int u, t;
    size_t LDAxK = LDA*NKmax;
    size_t LDCxN = LDC*N;

69 70 71 72
    CHAMELEON_Complex64_t *A      = (CHAMELEON_Complex64_t *)malloc(LDAxK*sizeof(CHAMELEON_Complex64_t));
    CHAMELEON_Complex64_t *C      = (CHAMELEON_Complex64_t *)malloc(LDCxN*sizeof(CHAMELEON_Complex64_t));
    CHAMELEON_Complex64_t *Cinit  = (CHAMELEON_Complex64_t *)malloc(LDCxN*sizeof(CHAMELEON_Complex64_t));
    CHAMELEON_Complex64_t *Cfinal = (CHAMELEON_Complex64_t *)malloc(LDCxN*sizeof(CHAMELEON_Complex64_t));
73 74

    /* Check if unable to allocate memory */
75 76 77
    if ( (!A) || (!C) || (!Cinit) || (!Cfinal) ){
        free(A); free(C);
        free(Cinit); free(Cfinal);
78 79 80 81 82 83 84
        printf("Out of Memory \n ");
        return -2;
    }

    eps = LAPACKE_dlamch_work('e');

    printf("\n");
85
    printf("------ TESTS FOR CHAMELEON ZSYRK ROUTINE -------  \n");
86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
    printf("            Size of the Matrix A %d by %d\n", N, K);
    printf("\n");
    printf(" The matrix A is randomly generated for each test.\n");
    printf("============\n");
    printf(" The relative machine precision (eps) is to be %e \n",eps);
    printf(" Computational tests pass if scaled residuals are less than 10.\n");

    /*----------------------------------------------------------
    *  TESTING ZSYRK
    */

    /* Initialize A */
    LAPACKE_zlarnv_work(IONE, ISEED, LDAxK, A);

    /* Initialize C */
101
    CHAMELEON_zplgsy( (double)0., ChamUpperLower, N, C, LDC, 51 );
102 103 104

    for (u=0; u<2; u++) {
        for (t=0; t<2; t++) {
105 106
            memcpy(Cinit,  C, LDCxN*sizeof(CHAMELEON_Complex64_t));
            memcpy(Cfinal, C, LDCxN*sizeof(CHAMELEON_Complex64_t));
107

108 109
            /* CHAMELEON ZSYRK */
            CHAMELEON_zsyrk(uplo[u], trans[t], N, K, alpha, A, LDA, beta, Cfinal, LDC);
110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137

            /* Check the solution */
            info_solution = check_solution(uplo[u], trans[t], N, K,
                                           alpha, A, LDA, beta, Cinit, Cfinal, LDC);

            if (info_solution == 0) {
                printf("***************************************************\n");
                printf(" ---- TESTING ZSYRK (%5s, %s) ........... PASSED !\n", uplostr[u], transstr[t]);
                printf("***************************************************\n");
            }
            else {
                printf("************************************************\n");
                printf(" - TESTING ZSYRK (%5s, %s) ... FAILED !\n", uplostr[u], transstr[t]);    hres++;
                printf("************************************************\n");
            }
        }
    }

    free(A); free(C);
    free(Cinit); free(Cfinal);

    return hres;
}

/*--------------------------------------------------------------
 * Check the solution
 */

138 139
static int check_solution(cham_uplo_t uplo, cham_trans_t trans, int N, int K,
                          CHAMELEON_Complex64_t alpha, CHAMELEON_Complex64_t *A, int LDA,
Mathieu Faverge's avatar
Mathieu Faverge committed
140
                          CHAMELEON_Complex64_t beta,  CHAMELEON_Complex64_t *Cref, CHAMELEON_Complex64_t *Ccham, int LDC)
141 142
{
    int info_solution;
Mathieu Faverge's avatar
Mathieu Faverge committed
143
    double Anorm, Cinitnorm, Cchamnorm, Clapacknorm, Rnorm;
144
    double eps;
145
    CHAMELEON_Complex64_t beta_const;
146 147 148 149
    double result;
    double *work = (double *)malloc(max(N, K)* sizeof(double));

    beta_const  = -1.0;
150
    Anorm       = LAPACKE_zlange_work(LAPACK_COL_MAJOR, 'I',
151 152
                                (trans == ChamNoTrans) ? N : K,
                                (trans == ChamNoTrans) ? K : N, A, LDA, work);
153
    Cinitnorm   = LAPACKE_zlange_work(LAPACK_COL_MAJOR, 'I', N, N, Cref,    LDC, work);
Mathieu Faverge's avatar
Mathieu Faverge committed
154
    Cchamnorm = LAPACKE_zlange_work(LAPACK_COL_MAJOR, 'I', N, N, Ccham, LDC, work);
155 156 157 158

    cblas_zsyrk(CblasColMajor, (CBLAS_UPLO)uplo, (CBLAS_TRANSPOSE)trans,
                N, K, CBLAS_SADDR(alpha), A, LDA, CBLAS_SADDR(beta), Cref, LDC);

159
    Clapacknorm = LAPACKE_zlange_work(LAPACK_COL_MAJOR, 'I', N, N, Cref, LDC, work);
160

Mathieu Faverge's avatar
Mathieu Faverge committed
161
    cblas_zaxpy(LDC*N, CBLAS_SADDR(beta_const), Ccham, 1, Cref, 1);
162

163
    Rnorm = LAPACKE_zlange_work(LAPACK_COL_MAJOR, 'I', N, N, Cref, LDC, work);
164 165 166

    eps = LAPACKE_dlamch_work('e');

Mathieu Faverge's avatar
Mathieu Faverge committed
167 168
    printf("Rnorm %e, Anorm %e, Cinitnorm %e, Cchamnorm %e, Clapacknorm %e\n",
           Rnorm, Anorm, Cinitnorm, Cchamnorm, Clapacknorm);
169 170 171 172 173

    result = Rnorm / ((Anorm + Cinitnorm) * N * eps);

    printf("============\n");
    printf("Checking the norm of the difference against reference ZSYRK \n");
Mathieu Faverge's avatar
Mathieu Faverge committed
174
    printf("-- ||Ccham - Clapack||_oo/((||A||_oo+||C||_oo).N.eps) = %e \n", result);
175

Mathieu Faverge's avatar
Mathieu Faverge committed
176
    if ( isinf(Clapacknorm) || isinf(Cchamnorm) || isnan(result) || isinf(result) || (result > 10.0) ) {
177 178 179 180 181 182 183 184 185 186 187 188
        printf("-- The solution is suspicious ! \n");
        info_solution = 1;
    }
    else {
        printf("-- The solution is CORRECT ! \n");
        info_solution= 0 ;
    }

    free(work);

    return info_solution;
}