MKL

#include <stdio.h>  
#include <stdlib.h>  
#include "mkl.h"  
  
int ()  
{  
    double *A, *B, *C;  
    int m, n, p, i, j;  
    double alpha, beta;  
  
    printf ("n This example computes real matrix C=alpha*A*B+beta*C using n"  
            " Intel(R) MKL function dgemm, where A, B, and  C are matrices and n"  
            " alpha and beta are double precision scalarsnn");  
  
    m = 2000, p = 200, n = 1000;  
    printf (" Initializing data for matrix multiplication C=A*B for matrix n"  
            " A(%ix%i) and matrix B(%ix%i)nn", m, p, p, n);  
    alpha = 1.0; beta = 0.0;  
    printf (" Allocating memory for matrices aligned on 64-byte boundary for better n"  
            " performance nn");  
    A = (double *)mkl_malloc( m*p*sizeof( double ), 64 );  
    B = (double *)mkl_malloc( p*n*sizeof( double ), 64 );  
    C = (double *)mkl_malloc( m*n*sizeof( double ), 64 );  
    if (A == NULL || B == NULL || C == NULL) {  
        printf( "n ERROR: Can't allocate memory for matrices. Aborting... nn");  
        mkl_free(A);  
        mkl_free(B);  
        mkl_free(C);  
        return 1;  
    }  
  
    printf (" Intializing matrix data nn");  
    for (i = 0; i < (m*p); i++) {  
        A[i] = (double)(i+1);  
    }  
  
    for (i = 0; i < (p*n); i++) {  
        B[i] = (double)(-i-1);  
    }  
  
    for (i = 0; i < (m*n); i++) {  
        C[i] = 0.0;  
    }  
  
    printf (" Computing matrix product using Intel(R) MKL dgemm function via CBLAS interface nn");  
    cblas_dgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,   
                m, n, p, alpha, A, p, B, n, beta, C, n);  
    printf ("n Computations completed.nn");  
  
    printf (" Top left corner of matrix A: n");  
    for (i=0; i<min(m,6); i++) {  
        for (j=0; j<min(p,6); j++) {  
            printf ("%12.0f", A[j+i*p]);  
        }  
        printf ("n");  
    }  
  
    printf ("n Top left corner of matrix B: n");  
    for (i=0; i<min(p,6); i++) {  
        for (j=0; j<min(n,6); j++) {  
            printf ("%12.0f", B[j+i*n]);  
        }  
        printf ("n");  
    }  
      
    printf ("n Top left corner of matrix C: n");  
    for (i=0; i<min(m,6); i++) {  
        for (j=0; j<min(n,6); j++) {  
            printf ("%12.5G", C[j+i*n]);  
        }  
        printf ("n");  
    }  
  
    printf ("n Deallocating memory nn");  
    mkl_free(A);  
    mkl_free(B);  
    mkl_free(C);  
  
    printf (" Example completed. nn");  
    return 0;  
}