See what's going on with flipcode!




This section of the archives stores flipcode's complete Developer Toolbox collection, featuring a variety of mini-articles and source code contributions from our readers.

 

  Matrix Mult Tests
  Submitted by



I was at an AMD performance talk this past week, and they were showing how one of their perf tools. To demonstrate, they had 3 different matrix multiplication implementations which they ran several million times. One implementation was AMD's 3dnow implementation, another was C code, and the other D3DXMatrixMultiply. Out of the 3, obviously the 3dnow version was the fastest. However, the C code was faster than what I thought would be optimized assembly in D3DXMatrixMultiply.

I was skeptical so I wrote up an inline assembly implementation and a C version and compared them vs D3DXMatrixMultiply. Running the test, I got the following results (measuing using QueryPerformanceCounter over 5 million matrix mult operations): Asm: 1923883
C : 3669129
D3DX: 3755883 Which lead me to believe that the D3DX implementation is only vanilla C... which is crap compared to assembly =P I sent my program to the work computer the next day to test them vs the latest build of D3DX, dx8.1 (i'm currently working for microsoft/directx), and the D3DX implementation was now around 15% faster vs assembly on the same machine. I traced into D3DXMatrixMultiply and the code was completely different from what i saw on my machine at home (dx8.0a); 3dnow ops were being used. I found out later that day that the D3DX team has been playing around with processor specific optimizations, and D3DX detects your proc and uses functions specific to the type detected. (which means there are SSE/SSE2 optimizations as well if u have p3/4) Anyways, this all fits into COTD cause I have the assembly implementation attached (along with the tests), in case anyone wanted it.

Currently browsing [multmatrixasm.zip] (4,499 bytes) - [matmult.cpp] - (16,204 bytes)

// matmult.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include <memory.h>
#include <iostream>
#include <windows.h>
#include <D3dx8math.h>
#include <stdio.h>
#include <assert.h>

using namespace std;

// Load an identity matrix void IdentityMatrixAsm( float *pfM ) { __asm { push ebx mov eax, pfM mov ebx, 00000000h // 0.0 in float mov dword ptr[ eax + 4 ], ebx mov dword ptr[ eax + 8 ], ebx mov dword ptr[ eax + 12 ], ebx

mov dword ptr[ eax + 16 ], ebx mov dword ptr[ eax + 24 ], ebx mov dword ptr[ eax + 28 ], ebx

mov dword ptr[ eax + 32 ], ebx mov dword ptr[ eax + 36 ], ebx mov dword ptr[ eax + 44 ], ebx

mov dword ptr[ eax + 48 ], ebx mov dword ptr[ eax + 52 ], ebx mov dword ptr[ eax + 56 ], ebx

mov ebx, 3f800000h // 1.0 in float mov dword ptr[ eax ], ebx mov dword ptr[ eax + 20 ], ebx mov dword ptr[ eax + 40 ], ebx mov dword ptr[ eax + 60 ], ebx

pop ebx } }

// loop unrolled matrix multiplication // Basically, I load a single row from Matrix 1, push it onto the floating point stack // and then get the vector product each column from Matrix 2 to get a row of results. // The assembly is no longer "straightforward", I optimized by grouping reads so // the their use is delayed as long as possible. (prevents lost cycles) // Since the initial row values are cached in the fp stack, it IS legal to do: // MulMatrixAsm( M1, M1, M2 ); for matrix stack operations // MulMatrixAsm( M2, M1, M2 ); is not possible though. void MultMatrixAsm( float *pfResults, const float *pfM1, const float *pfM2 ) { // assertions dont show up in release builds, so it incurs no performance loss. assert( pfResults != pfM2 );

__asm { push ebx push edx

mov eax, pfM1 mov ebx, pfM2 mov edx, pfResults

fld dword ptr[ eax ] fld dword ptr[ eax + 4 ] fld dword ptr[ eax + 8 ] fld dword ptr[ eax + 12 ] fld dword ptr[ ebx ] fmul st(0), st(4) fld dword ptr[ ebx + 4 ] fmul st(0), st(5) fld dword ptr[ ebx + 16 ] fld dword ptr[ ebx + 20 ] fmul st(0), st(6) faddp st(2), st(0) fmul st(0), st(5) faddp st(2), st(0) fld dword ptr[ ebx + 32 ] fld dword ptr[ ebx + 36 ] fmul st(0), st(5) faddp st(2), st(0) fmul st(0), st(4) faddp st(2), st(0) fld dword ptr[ ebx + 48 ] fld dword ptr[ ebx + 52 ] fmul st(0), st(4) faddp st(2), st(0) fmul st(0), st(3) faddp st(2), st(0) fstp dword ptr[ edx + 4 ] fstp dword ptr[ edx ]

fld dword ptr[ ebx + 12 ] fmul st(0), st(4) fld dword ptr[ ebx + 8 ] fmulp st(5), st(0) fld dword ptr[ ebx + 24 ] fld dword ptr[ ebx + 28 ] fmul st(0), st(5) faddp st(2), st(0) fmulp st(4), st(0) fld dword ptr[ ebx + 40 ] fld dword ptr[ ebx + 44 ] fmul st(0), st(4) faddp st(2), st(0) fmulp st(3), st(0) fld dword ptr[ ebx + 56 ] fld dword ptr[ ebx + 60 ] fmul st(0), st(3) faddp st(2), st(0) fmulp st(2), st(0) fstp dword ptr[ edx + 12 ] faddp st(1), st(0) faddp st(1), st(0) faddp st(1), st(0) fstp dword ptr[ edx + 8 ]





fld dword ptr[ eax + 16 ] fld dword ptr[ eax + 20 ] fld dword ptr[ eax + 24 ] fld dword ptr[ eax + 28 ] fld dword ptr[ ebx ] fmul st(0), st(4) fld dword ptr[ ebx + 4 ] fmul st(0), st(5) fld dword ptr[ ebx + 16 ] fld dword ptr[ ebx + 20 ] fmul st(0), st(6) faddp st(2), st(0) fmul st(0), st(5) faddp st(2), st(0) fld dword ptr[ ebx + 32 ] fld dword ptr[ ebx + 36 ] fmul st(0), st(5) faddp st(2), st(0) fmul st(0), st(4) faddp st(2), st(0) fld dword ptr[ ebx + 48 ] fld dword ptr[ ebx + 52 ] fmul st(0), st(4) faddp st(2), st(0) fmul st(0), st(3) faddp st(2), st(0) fstp dword ptr[ edx + 20 ] fstp dword ptr[ edx + 16 ]

fld dword ptr[ ebx + 12 ] fmul st(0), st(4) fld dword ptr[ ebx + 8 ] fmulp st(5), st(0) fld dword ptr[ ebx + 24 ] fld dword ptr[ ebx + 28 ] fmul st(0), st(5) faddp st(2), st(0) fmulp st(4), st(0) fld dword ptr[ ebx + 40 ] fld dword ptr[ ebx + 44 ] fmul st(0), st(4) faddp st(2), st(0) fmulp st(3), st(0) fld dword ptr[ ebx + 56 ] fld dword ptr[ ebx + 60 ] fmul st(0), st(3) faddp st(2), st(0) fmulp st(2), st(0) fstp dword ptr[ edx + 28 ] faddp st(1), st(0) faddp st(1), st(0) faddp st(1), st(0) fstp dword ptr[ edx + 24 ]





fld dword ptr[ eax + 32 ] fld dword ptr[ eax + 36 ] fld dword ptr[ eax + 40 ] fld dword ptr[ eax + 44 ] fld dword ptr[ ebx ] fmul st(0), st(4) fld dword ptr[ ebx + 4 ] fmul st(0), st(5) fld dword ptr[ ebx + 16 ] fld dword ptr[ ebx + 20 ] fmul st(0), st(6) faddp st(2), st(0) fmul st(0), st(5) faddp st(2), st(0) fld dword ptr[ ebx + 32 ] fld dword ptr[ ebx + 36 ] fmul st(0), st(5) faddp st(2), st(0) fmul st(0), st(4) faddp st(2), st(0) fld dword ptr[ ebx + 48 ] fld dword ptr[ ebx + 52 ] fmul st(0), st(4) faddp st(2), st(0) fmul st(0), st(3) faddp st(2), st(0) fstp dword ptr[ edx + 36 ] fstp dword ptr[ edx + 32 ]

fld dword ptr[ ebx + 12 ] fmul st(0), st(4) fld dword ptr[ ebx + 8 ] fmulp st(5), st(0) fld dword ptr[ ebx + 24 ] fld dword ptr[ ebx + 28 ] fmul st(0), st(5) faddp st(2), st(0) fmulp st(4), st(0) fld dword ptr[ ebx + 40 ] fld dword ptr[ ebx + 44 ] fmul st(0), st(4) faddp st(2), st(0) fmulp st(3), st(0) fld dword ptr[ ebx + 56 ] fld dword ptr[ ebx + 60 ] fmul st(0), st(3) faddp st(2), st(0) fmulp st(2), st(0) fstp dword ptr[ edx + 44 ] faddp st(1), st(0) faddp st(1), st(0) faddp st(1), st(0) fstp dword ptr[ edx + 40 ]





fld dword ptr[ eax + 48 ] fld dword ptr[ eax + 52 ] fld dword ptr[ eax + 56 ] fld dword ptr[ eax + 60 ] fld dword ptr[ ebx ] fmul st(0), st(4) fld dword ptr[ ebx + 4 ] fmul st(0), st(5) fld dword ptr[ ebx + 16 ] fld dword ptr[ ebx + 20 ] fmul st(0), st(6) faddp st(2), st(0) fmul st(0), st(5) faddp st(2), st(0) fld dword ptr[ ebx + 32 ] fld dword ptr[ ebx + 36 ] fmul st(0), st(5) faddp st(2), st(0) fmul st(0), st(4) faddp st(2), st(0) fld dword ptr[ ebx + 48 ] fld dword ptr[ ebx + 52 ] fmul st(0), st(4) faddp st(2), st(0) fmul st(0), st(3) faddp st(2), st(0) fstp dword ptr[ edx + 52 ] fstp dword ptr[ edx + 48 ]

fld dword ptr[ ebx + 12 ] fmul st(0), st(4) fld dword ptr[ ebx + 8 ] fmulp st(5), st(0) fld dword ptr[ ebx + 24 ] fld dword ptr[ ebx + 28 ] fmul st(0), st(5) faddp st(2), st(0) fmulp st(4), st(0) fld dword ptr[ ebx + 40 ] fld dword ptr[ ebx + 44 ] fmul st(0), st(4) faddp st(2), st(0) fmulp st(3), st(0) fld dword ptr[ ebx + 56 ] fld dword ptr[ ebx + 60 ] fmul st(0), st(3) faddp st(2), st(0) fmulp st(2), st(0) fstp dword ptr[ edx + 60 ] faddp st(1), st(0) faddp st(1), st(0) faddp st(1), st(0) fstp dword ptr[ edx + 56 ]



pop edx pop ebx } }

/* // This is the loop version of the above... __forceinline void MultMatrixAsm( float *pfResults, const float *pfM1, const float *pfM2 ) { // assertions dont show up in release builds, so it incurs no performance loss. assert( Results != M2 );

__asm { push ebx push ecx push edx

mov eax, pfM1 mov ebx, pfM2 mov edx, pfResults mov ecx, 48



row: fld dword ptr[ eax + ecx ] fld dword ptr[ eax + ecx + 4 ] fld dword ptr[ eax + ecx + 8 ] fld dword ptr[ eax + ecx + 12 ]

fld dword ptr[ ebx ] fmul st(0), st(4) fld dword ptr[ ebx + 4 ] fmul st(0), st(5) fld dword ptr[ ebx + 16 ] fld dword ptr[ ebx + 20 ] fmul st(0), st(6) faddp st(2), st(0) fmul st(0), st(5) faddp st(2), st(0) fld dword ptr[ ebx + 32 ] fld dword ptr[ ebx + 36 ] fmul st(0), st(5) faddp st(2), st(0) fmul st(0), st(4) faddp st(2), st(0) fld dword ptr[ ebx + 48 ] fld dword ptr[ ebx + 52 ] fmul st(0), st(4) faddp st(2), st(0) fmul st(0), st(3) faddp st(2), st(0) fstp dword ptr[ edx + ecx + 4 ] fstp dword ptr[ edx + ecx ]

fld dword ptr[ ebx + 12 ] fmul st(0), st(4) fld dword ptr[ ebx + 8 ] fmulp st(5), st(0) fld dword ptr[ ebx + 24 ] fld dword ptr[ ebx + 28 ] fmul st(0), st(5) faddp st(2), st(0) fmulp st(4), st(0) fld dword ptr[ ebx + 40 ] fld dword ptr[ ebx + 44 ] fmul st(0), st(4) faddp st(2), st(0) fmulp st(3), st(0) fld dword ptr[ ebx + 56 ] fld dword ptr[ ebx + 60 ] fmul st(0), st(3) faddp st(2), st(0) fmulp st(2), st(0) fstp dword ptr[ edx + ecx + 12 ] faddp st(1), st(0) faddp st(1), st(0) faddp st(1), st(0) fstp dword ptr[ edx + ecx + 8 ]

sub ecx, 16 jge row pop edx pop ecx pop ebx } } */


// C Version.. VC7 can loop unroll this, but its not that much faster still.. // M1, M2, and Results must be different pointers, a version could be written void CMultMatrix( float *pfResults, float *pfM1, float *pfM2 ) { // assertions dont show up in release builds, so it incurs no performance loss. assert( pfResults != pfM1 ); assert( pfResults != pfM2 );

for( int i = 0; i < 16; i+=4 ) { for ( int j = 0; j < 4; j++ ) pfResults[ i + j ] = pfM1[ i ] * pfM2[ j ] + pfM1[ i + 1 ] * pfM2[ j + 4 ] + pfM1[ i + 2 ] * pfM2[ j + 8 ] + pfM1[ i + 3 ] * pfM2[ j + 12 ]; } }

void PrintMatrix( float * M ) { cout << endl; for ( int i = 0; i < 16; i+=4 ) { for ( int j = 0; j < 4; j++ ) cout << M[ i + j ] << "\t"; cout << endl; } }

int __cdecl main(int argc, char* argv[]) { float M1[16]; float M2[16]; float Results[16]; int i;

for( i = 0; i < 16; i++ ) { M1[ i ] = 16.0f - i; M2[ i ] = i; } LARGE_INTEGER CTestStart, CTestFinish; LARGE_INTEGER D3DTestStart, D3DTestFinish; LARGE_INTEGER AsmTestStart, AsmTestFinish;

// Run CMultMatrix test memset( Results, 0, sizeof( float ) * 16 ); QueryPerformanceCounter( &CTestStart ); for ( i = 0; i < 5000000; i++ ) CMultMatrix( Results, M1, M2 ); QueryPerformanceCounter( &CTestFinish ); PrintMatrix( Results ); LARGE_INTEGER CTestTime; CTestTime.QuadPart = CTestFinish.QuadPart - CTestStart.QuadPart; // number too big to print assert( CTestTime.HighPart == 0 ); cout << "C :\t" << CTestTime.LowPart << endl;

// Run D3DXMatrixMultiply test D3DXMATRIX D3DResults; D3DXMATRIX D3DM1(M1); D3DXMATRIX D3DM2(M2); memset( D3DResults, 0, sizeof( float ) * 16 ); QueryPerformanceCounter( &D3DTestStart ); for ( i = 0; i < 5000000; i++ ) D3DXMatrixMultiply( &D3DResults, &D3DM1, &D3DM2 ); QueryPerformanceCounter( &D3DTestFinish ); PrintMatrix( D3DResults ); LARGE_INTEGER D3DTestTime; D3DTestTime.QuadPart = D3DTestFinish.QuadPart - D3DTestStart.QuadPart; // number too big to print assert( D3DTestTime.HighPart == 0 ); cout << "D3D:\t" << D3DTestTime.LowPart << endl;

// Run MultMatrixAsm test memset( Results, 0, sizeof( float ) * 16 ); QueryPerformanceCounter( &AsmTestStart ); for ( i = 0; i < 5000000; i++ ) MultMatrixAsm( Results, M1, M2 ); QueryPerformanceCounter( &AsmTestFinish ); PrintMatrix( Results ); LARGE_INTEGER AsmTestTime; AsmTestTime.QuadPart = AsmTestFinish.QuadPart - AsmTestStart.QuadPart; // number too big to print assert( AsmTestTime.HighPart == 0 ); cout << "Asm:\t" << AsmTestTime.LowPart << endl;

getchar();

return 0; }

Currently browsing [multmatrixasm.zip] (4,499 bytes) - [StdAfx.cpp] - (294 bytes)

// stdafx.cpp : source file that includes just the standard includes
//	matmult.pch will be the pre-compiled header
//	stdafx.obj will contain the pre-compiled type information

#include "stdafx.h"

// TODO: reference any additional headers you need in STDAFX.H // and not in this file

Currently browsing [multmatrixasm.zip] (4,499 bytes) - [StdAfx.h] - (769 bytes)

// stdafx.h : include file for standard system include files,
//  or project specific include files that are used frequently, but
//      are changed infrequently
//

#if !defined(AFX_STDAFX_H__8F4EF30C_35AE_46C4_B626_F3D0050F8467__INCLUDED_)
#define AFX_STDAFX_H__8F4EF30C_35AE_46C4_B626_F3D0050F8467__INCLUDED_

#if _MSC_VER > 1000 #pragma once #endif // _MSC_VER > 1000 #define WIN32_LEAN_AND_MEAN // Exclude rarely-used stuff from Windows headers #include <stdio.h>

// TODO: reference additional headers your program requires here //{{AFX_INSERT_LOCATION}} // Microsoft Visual C++ will insert additional declarations immediately before the previous line. #endif // !defined(AFX_STDAFX_H__8F4EF30C_35AE_46C4_B626_F3D0050F8467__INCLUDED_)

The zip file viewer built into the Developer Toolbox made use of the zlib library, as well as the zlibdll source additions.

 

Copyright 1999-2008 (C) FLIPCODE.COM and/or the original content author(s). All rights reserved.
Please read our Terms, Conditions, and Privacy information.