#include <windows.h>
#include <malloc.h>
#include <math.h>
#include <stdio.h>
#include <assert.h>
#include "..\inc\ezsimd.h"
#include "..\inc\checksimd.h"
// (a) 通常方式
// out[i] = a[i] * b[i]の計算をする
void mul_float(float *out, float *a, float *b, int len)
{
int i;
for(i=0;i<len;i++){
out[i] = a[i] * b[i];
}
}
// (b) SIMD使用
// SIMDを使用して記述。a,b配列から4つづつXMM0,XMM1レジスタに値を転送し
// mulps (multiply packed single-fp) 命令で4float同時に掛け算する
// 結果をXMM0からoutputに転送する。
// lenの値は４の倍数でなければならない
void mul_float4_simd_u(float *output, float *a, float *b, int len)
{
int i;
assert(len % 4 == 0); //lenは４の倍数でなければならない
// mmxレジスターは128bit= 32*4の大きさを持ち一度にfloatを４つ処理できる
// ループ回数は1/4
len = len / 4;
for(i=0;i<len;i++){
_asm mov edi, a; // 次のmovups( )は直接C言語変数を引数にとれないので一旦レジスターに値移す
movups(XMM0, P_EDI); // xmm0 = [edi]; xmm0にaの示すfloat(4word)を転送する
a+=4;
_asm mov edi, b;
movups(XMM1, P_EDI); // xmm1 = [edi]; xmm1にaの示すfloat(4word)を転送する
b+=4;
mulps(XMM0, XMM1); // xmm0 = mmx0 * mmx1; 4word同時に計算される
_asm mov edi, output;
movups(P_EDI, XMM0); // oの示すアドレスにmmx0のfloat(4word)を転送する
output += 4;
}
}
// (c) SIMD使用
// XMMレジスターとメモリーとの転送にmovaps(move aligned four packed single-fp)
// を使う。movupsにくらべて高速に転送できるが転送元・先のアドレスが16byte
// 境界（アドレス値が16の倍数）でなければならない
// ポインターをレジスターに置き換え高速化
void mul_float4_simd_a(float *output, float *a, float *b, int len)
{
assert(len % 4 == 0); //lenは４の倍数でなければならない
// movapsに使うポインターは16byte境界にアライメントされていなければならない
assert(((unsigned int)a % 16) == 0);
assert(((unsigned int)b % 16) == 0);
assert(((unsigned int)output % 16) == 0);
_asm {
mov ecx, len; // ecx = len
shr ecx, 2; // ecx = ecx/4
mov ebx, a; // ebx = a
mov edx, b; // edx = b
mov edi, output; // edi = output
loop_start:
movaps(XMM0, P_EBX); // mmx0 = [ebx]
add ebx, 16; // ebx = ebx + sizeof(float) * 4
movaps(XMM1, P_EDX); // mmx1 = [edx]
add edx, 16; // ebx = edx + sizeof(float) * 4
mulps(XMM0, XMM1); // mmx0 = mmx0 * mmx1
movaps(P_EDI, XMM0); // [edi] = mmx0
add edi, 16; // edi = edi + sizeof(float) * 4
dec ecx; // ecx = ecx -1;
jnz loop_start; // if ecx!=0 goto loop_start
}
}
void dump_float(float *data, int len)
{
int i;
for(i=0;i<len;i++){
if(i % 4 ==0) printf("%d : ",i);
printf(" %10f",data[i]);
if(i % 4 ==3) printf("\n");
}
printf("\n");
}
void main(void)
{
float a[] = {1,2,3,4,5,6,7,8};
float b[] = {2,3,4,5,6,7,8,9};
float out1[8], out2[8];
float *mem1, *mem2, *mem3;
float *a16, *b16, *out3;
int status;
status = check_simd_support();
if( !(status & STATUS_CPU_SUPPORT) ){
printf("このCPUはSIMD命令をサポートしていません\n");
return;
}
if( !(status & STATUS_OS_SUPPORT) ){
printf("このOSはSIMD命令をサポートしていません\n");
return;
}
// 通常演算
mul_float(out1, a, b, 8);
printf("normal function results\n");
dump_float(out1, 8);
// SIMD使用（１）
mul_float4_simd_u(out2, a, b, 8);
printf("simd function results\n");
dump_float(out2, 8);
// SIMD使用（２）の実行には16byte境界にアライメントされた
// メモリーが必要になる。それを用意する
// 16 byte境界に切り上げる
#define ALIGNEFLOAT16(p) (float*) ((((unsigned int)(p))+15) &0xfffffff0);
// 16byte境界に切り上げしても容量が不足しないように15byte多く確保する
mem1 = (float *)malloc(sizeof(float) * 8 +15);
mem2 = (float *)malloc(sizeof(float) * 8 +15);
mem3 = (float *)malloc(sizeof(float) * 8 +15);
a16 = ALIGNEFLOAT16(mem1);
b16 = ALIGNEFLOAT16(mem2);
out3 = ALIGNEFLOAT16(mem3);
memcpy(a16, a, sizeof(float) * 8);
memcpy(b16, b, sizeof(float) * 8);
// SIMD使用（２）
mul_float4_simd_a(out3, a16, b16, 8);
printf("simd2 function results\n");
dump_float(out3, 8);
free(mem1);// 間違えてfree(a16)としないように注意
free(mem2);
free(mem3);
}
