Webseite von Michael Rennecke

#include <iostream>
#include <stdlib.h>
#include <time.h>
 
using namespace std;
 
void mul(double *A, double *B, double *C, int dim)
{
    int i, j, k;
    double s;
 
    for (i = 0; i < dim; i++){
        for (j = 0; j < dim; j++) {
            s = 0.0;
            for (k = 0; k < dim; ++k)
                s += A[dim*i +k] * B[dim*k+j];
 
            C[dim*i+j] = s;
        }
    }
}
 
void transpose(double *A, int dim){
 
        int i,j;
        double tmp;
 
    for (i = 0; i < dim; i++){
        for (j = i; j < dim; j++){
            tmp = A[dim*i+j];
            A[dim*i+j] =A[dim*j+i];
            A[dim*j+i] = tmp;
        }
    }
 
}
 
void mulfast(double *A, double *B, double *C, int dim)
{
    int i, j, k;
    double s;
 
    transpose(B, dim);
 
    for (i = 0; i < dim; ++i){
        for (j = 0; j < dim; ++j) {
            s = 0.0;
            for (k = 0; k < dim; ++k)
                // B is transposed !!
                s += A[dim*i +k] * B[dim*j+k];
 
            C[dim*i+j] = s;
        }
    }
 
    transpose(B, dim);
}
 
int main(){
 
    int dim = 1024;
    clock_t start, end;
 
    double* A = (double*) calloc(dim* dim, sizeof(double));
    double* B = (double*) calloc(dim* dim, sizeof(double));
    double* C = (double*) calloc(dim* dim, sizeof(double));
 
    start = clock();
    mul(A,B,C, dim);
    end = clock();
    cout << "normal implementation: " << (end-start) / CLOCKS_PER_SEC << " seconds" << endl;
 
    start = clock();
    mulfast(A,B,C, dim);
    end = clock();
    cout << "fast implementation: " << (end-start) / CLOCKS_PER_SEC << " seconds" << endl;
 
    return 0;
}