Comments 74
В статье очень обрывистое и несвязное объяснение, которое все сводится к одному:
Используя векторные команды процессора SSE можно ускорить произведение матриц. Вот тривиальный код для случая 4x4. Мои тесты показывают прирост в 10 раз.
Никаких откровений и трюков тут нет. Но в качестве пособия для ознакомления эта статья тоже не годится, ибо там надо гораздо основательнее, подробнее и медленнее объяснять что это такое, зачем и почему. Также надо уделить внимание транспонированию матриц. Вы точно этот код запускали и он точно вам один и тот же ответ выдает с наивным методом?
спасибо, извините, насчет порядка замера и хода расчета в статье указан результат, насчет платформы указано в статье, извините что не предоставил Транспонирование вы можете посмотреть как делать транспонирование в соседней статье по расчету инверсной матрицы, эта статья конкретно про эти 2 метода и о том где эти файлики и что можно начать по-тихоньку пробовать начиная с перемножения матриц через тесты интегрировать 128 регистры в расчеты я считаю это может кому-то помочь, код протестирован на реализации граффики и граффика вздохнула с новой скоростью, насчет тонких нюансов математики тоже нет описаний в этой статье.
ускорение х10 может помочь на старте, дальше просто читаем файл тестируем и пр.
Вам имело смысл соревноваться не с Си кодом, а с библиотеками, вот в Intel MKL вроде это есть, ну и там gemm() из LAPACK, по-моему и то и другое под никсами живёт (хотя я тут не спец). Вам имело смысл попробовать реализовать это дело для матриц произвольного размера, как одинарной так и двойной точности и если бы Вы обогнали по скорости существующие библиотеки, то получилась бы хорошая статья. Кроме того, у Вас своеобразный стиль изложения (Вы дайте кому-нибубдь прочитать перед публикацией), ну и по коду можно придраться, ну вот зачем столько файлов включать для интрисиков, там же один immintrin.h всего нужен, смысл верчения циклов "с ускорением" от меня тоже ускользает, и т.д.
спасибо тут не про сравнения с lapack, а вводная как начать через перемножения матриц пользоваться ускоряющимися регистрами, там есть еще некоторые моменты сходу когда я писал статью тоже думал ну и вставил пока так, дело в том что там по файликам разбросаны некоторые возможно удобные функции, и поэтому я выделил их и вставил, иминтрин файл как головной от проекта по ускорениям регистров, как таблица, и каждый новый файл дополняет части функционала, которые нужны возможно в математике, но вот в 3д я даже както смог глазом заметить, что получше, хотя это конечно замер глазом не обьективно, но побыстрее работает вообщем, библиотека написана в рамках 3д, не со смыслом гонятся с другими библиотеками, смысл моей библиотеки быть стабильной, и не иметь развития, потомучто библиотеки развиваются, был выбран курс lts - для этого пришлось всё реализовать
как начать через перемножения матриц пользоваться ускоряющимися регистрами
Вы термин "ускоряющихся регистров" сами придумали? Вот я бы порекомендовал раздобыть книжку Даниэля Куссвюрма "Профессиональное программирование на ассемблере x64 с расширениями AVX, AVX2 и AVX-512" и прочесть, хотя бы ради терминологии. Второе издание переведено на руский, третье пока нет, вроде. Лучше возьмите обе, на английском тоже "Kusswurm Daniel — Modern X86 Assembly Language Programming" (там в русском таки есть косяки перевода). Хотя эти книжки и не идеальны, но в них есть ровно всё то, что Вы выше изложили, причём весь код из книги издательством выложен на Гитхабе, включая использование FMA и произвольные размеры матриц. Если на интрисиках, то "Modern Parallel Programming with C++ and Assembly Language" того же автора, а код к ней вот здесь.
вообще, статья умышленно написана таким образом спасибо что вы вот написали это, дело в том что не везде есть поддержка bmi avx-подобных и в идеале придётся проверять, поэтому начало с имминтрин, я указал минимально в ресурсах где можно войти в курс дела на вики тоже описано что это, потом далее, придётся проверять если это касается запуска на других платформах такого кода с расширениями, процессор архитектуру и возможность расширения, это не указано в статье, прошу простите,
насчет того что я не професионально выразился простите, когда мы приходим в 3д ассемблер не на первом месте мотивация происходит из визуала первый вопрос у новичка почему медленно, и далее он уже разбирается если поставлена задача не глм, вообще я с вами солидарен, спасибо
не везде есть поддержка bmi avx-подобных и в идеале придётся проверять, поэтому начало с имминтрин,
Ну почитайте же книжки или хотя бы документацию. этот immintrin - это intrinsics - специальные функции компилятора, которые перобразуются в векторные инструкции процессоров. Eсли целевая архитектура их поддерживает. А если не поддерживает, то они заменяются обычными, не векторными инструкциями.
почитать и применить разное я как раз в своих ответах пытаюсь натолкнуть на такую нить, что знать и применить разное, придётся так и так проверить, я с вами солидарен
в случае этой статьи ради интереса можно проверить математику, указан ресурс с хорошими описаниями и картинками, этого будет достаточно в простейшем калькуляторе, матрицы view projection и шейдерная составляющая дотяните например на learnopengl, попробуйте найти или понять или попробовать эти моменты (view, projection поищите в интернете) попробуйте погрузитесь в математику, читать книги по интрисинкам есть смысл, но сначала надо удостовериться нужно ли точно? прикинуть примеры, и вообще сформировать свой подход к получению ответа
новичек будет мотивирован визуальной составляющей, если новичек догадался что ему нет смысла торопиться в анриал первая рекомендация будет SSE, ее проще включить и более вероятно что она есть, а AVX придётся проверить
по части графики до реализации разбиения пространства я могу ответить если вам интересно, а так же о наблюдениях на платформах, и паре библиотек, движками не пользовался извините, не пользовался так же реактивной разработкой, зато заходил даже в Lisp( лисп на уровне захотел написал, лучший язык для прототипа - на любителя )
А если не поддерживает
То будет ошибка компиляции, никакой эмуляции на скалярах нет:
g++ foo.cpp
foo.cpp: In function ‘int main()’:
foo.cpp:8:5: warning: AVX512F vector return without AVX512F enabled changes the ABI [-Wpsabi]
8 | res = _mm512_mask_add_ps(src, k, a, b);
| ~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
In file included from /usr/lib/gcc/x86_64-linux-gnu/9/include/immintrin.h:55,
from foo.cpp:1:
/usr/lib/gcc/x86_64-linux-gnu/9/include/avx512fintrin.h:12301:1: error: inlining failed in call to always_inline ‘__m512 _mm512_mask_add_ps(__m512, __mmask16, __m512, __m512)’: target specific option mismatch
12301 | _mm512_mask_add_ps (__m512 __W, __mmask16 __U, __m512 __A, __m512 __B)
| ^~~~~~~~~~~~~~~~~~
foo.cpp:8:25: note: called from here
8 | res = _mm512_mask_add_ps(src, k, a, b);
| ~~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~
Для компиляции обязательно указывать соответствующий -march
сравните скорость с этим если вам интересно
void atTheTestMatMul(const float *a,const float *b)
{
/* vec4 a1=(vec4){2,9,40,5};vec4 b1=(vec4){50 ,85,89,99}; */
/* vec4 a2=(vec4){8,6,5,6};vec4 b2=(vec4){30,100,65,45}; */
/* vec4 a3=(vec4){8,9,7,4};vec4 b3=(vec4){88,0,50,14}; */
/* vec4 a4=(vec4){7,5,3,10};vec4 b4=(vec4){70,10,60,80}; */
vec4 c1;
for(int i=0;i<4;++i)
{
c1=(vec4){0,0,0,0};
float r=*(b+i);
c1=Mulv4((vec4){*(a),*(a+1),*(a+2),*(a+3)}, (vec4){r,r,r,r});
r=*(b+i)+4;
c1=Addv4(c1,Mulv4((vec4){*(a+4),*(a+5),*(a+6),*(a+7)}, (vec4){r,r,r,r}));
r=*(b+i)+8;
c1=Addv4(c1,Mulv4((vec4){*(a+8),*(a+9),*(a+10),*(a+11)}, (vec4){r,r,r,r}));
r=*(b+i)+12;
c1=Addv4(c1,Mulv4((vec4){*(a+12),*(a+13),*(a+14),*(a+15)}, (vec4){r,r,r,r}));
printf("%f %f %f %f\n",c1.x,c1.y,c1.z,c1.w);
}
}
это прям по канонам С всё на указателях и сдвиг от первого адреса массива по сути и есть ускорение надо попасть просто в теже регистры
на данный момент тут ошибка в шифте, но работает быстро
А реализация Addv4/Mulv4 где?
Скрытый текст
#pragma GCC optimize("Ofast,unroll-loops,-ffast-math")
#pragma GCC target("sse")
#include <xmmintrin.h>
#include <wmmintrin.h>
#include <smmintrin.h>
#include <nmmintrin.h>
#include <immintrin.h>
#include <stdio.h>
#include <unistd.h>
#include <time.h>
#include <errno.h>
#include <string.h>
typedef float f32;
struct vec4_t
{
f32 x;
f32 y;
f32 z;
f32 w;
};
typedef struct vec4_t vec4;
struct mat4_t
{
f32 v[16];
};
typedef struct mat4_t mat4;
vec4 Mulv4(vec4 a, vec4 b)
{
return (vec4){a.x * b.x, a.y * b.y, a.z * b.z, a.w * b.w};
}
vec4 Addv4(const vec4 a, const vec4 b)
{
return (vec4){a.x + b.x, a.y + b.y, a.z + b.z, a.w + b.w};
}
void atTheTestMatMul(const float *a,const float *b)
{
/* vec4 a1=(vec4){2,9,40,5};vec4 b1=(vec4){50 ,85,89,99}; */
/* vec4 a2=(vec4){8,6,5,6};vec4 b2=(vec4){30,100,65,45}; */
/* vec4 a3=(vec4){8,9,7,4};vec4 b3=(vec4){88,0,50,14}; */
/* vec4 a4=(vec4){7,5,3,10};vec4 b4=(vec4){70,10,60,80}; */
for(int i=0;i<4;i++)
{
vec4 c1=(vec4){0,0,0,0};
int ShapeShifter=i*4;
float r=*(b+ShapeShifter);
c1=Mulv4((vec4){*(a),*(a+1),*(a+2),*(a+3)}, (vec4){r,r,r,r});
r=*(b+ShapeShifter+1);
c1=Addv4(c1,Mulv4((vec4){*(a+4),*(a+5),*(a+6),*(a+7)}, (vec4){r,r,r,r}));
r=*(b+ShapeShifter+2);
c1=Addv4(c1,Mulv4((vec4){*(a+8),*(a+9),*(a+10),*(a+11)}, (vec4){r,r,r,r}));
r=*(b+ShapeShifter+3);
c1=Addv4(c1,Mulv4((vec4){*(a+12),*(a+13),*(a+14),*(a+15)}, (vec4){r,r,r,r}));
printf("%f %f %f %f\n",c1.x,c1.y,c1.z,c1.w);
}
}
int main()
{
mat4 a=(mat4){
2,9,40,5,
8,6,5,6,
8,9,7,4,
7,5,3,10
};
mat4 b=(mat4){
50,30,88,70,
85,100,0,10,
89,65,50,60,
99,45,14,80
};
struct timespec Ta,Tb;
float GlobalTime=0;
float localTime=0;
clock_gettime(CLOCK_MONOTONIC,&Ta);
atTheTestMatMul(&a.v[0],&b.v[0]);
clock_gettime(CLOCK_MONOTONIC,&Tb);
localTime=(float)(Tb.tv_nsec - Ta.tv_nsec) / 1000000 +(float)(Tb.tv_sec - Ta.tv_sec);
printf("%f\n",localTime);
return 0;
}
странно в год болте не то а у меня на компе то, вот пофиксил
эх вот на моем локальном пк эта реализация быстро тоже работает вроде
int ShapeShifter=i<<2; так будет стабильнее вроде
вообще интересно конечно она время показывает ппц просто
#pragma GCC optimize("Ofast,unroll-loops,-ffast-math")
#pragma GCC target("sse")
#include <xmmintrin.h>
#include <wmmintrin.h>
#include <smmintrin.h>
#include <nmmintrin.h>
#include <immintrin.h>
#include <stdio.h>
#include <unistd.h>
#include <time.h>
#include <errno.h>
#include <string.h>
#pragma pack(1)
typedef float f32;
struct vec4_t
{
f32 x;
f32 y;
f32 z;
f32 w;
};
typedef struct vec4_t vec4;
struct mat4_t
{
f32 v[16];
};
typedef struct mat4_t mat4;
vec4 Mulv4(vec4 a, vec4 b)
{
return (vec4){a.x * b.x, a.y * b.y, a.z * b.z, a.w * b.w};
}
vec4 Addv4(const vec4 a, const vec4 b)
{
return (vec4){a.x + b.x, a.y + b.y, a.z + b.z, a.w + b.w};
}
vec4 Processing(vec4 c1,vec4 a,vec4 b,int shift)
{
return Addv4(c1,Mulv4(a,b));
}
vec4 TakeIting(vec4 c1,const float *a,const float *b,int shiftA,int ShapeShifter)
{
return Processing(c1,(vec4){*(a+shiftA),*(a+shiftA+1),*(a+shiftA+2),*(a+shiftA+3)},(vec4){*(b+ShapeShifter),*(b+ShapeShifter),*(b+ShapeShifter),*(b+ShapeShifter)},ShapeShifter);
}
float* Seting(float *a,int shift,vec4 b)
{
*(a+shift)=b.x;
*(a+shift+1)=b.y;
*(a+shift+2)=b.z;
*(a+shift+3)=b.w;
return a;
}
void atTheTestMatMul(float *result,const float *a,const float *b)
{
vec4 c1=(vec4){0,0,0,0};
int ShapeShifter=0<<2;
c1=TakeIting(c1,a,b,0,ShapeShifter+0);
c1=TakeIting(c1,a,b,4,ShapeShifter+1);
c1=TakeIting(c1,a,b,8,ShapeShifter+2);
c1=TakeIting(c1,a,b,12,ShapeShifter+3);
Seting(result,0*4,c1);
ShapeShifter=1<<2;
c1=TakeIting(c1,a,b,0,ShapeShifter+0);
c1=TakeIting(c1,a,b,4,ShapeShifter+1);
c1=TakeIting(c1,a,b,8,ShapeShifter+2);
c1=TakeIting(c1,a,b,12,ShapeShifter+3);
Seting(result,1*4,c1);
c1=(vec4){0,0,0,0};
ShapeShifter=2<<2;
c1=TakeIting(c1,a,b,0,ShapeShifter+0);
c1=TakeIting(c1,a,b,4,ShapeShifter+1);
c1=TakeIting(c1,a,b,8,ShapeShifter+2);
c1=TakeIting(c1,a,b,12,ShapeShifter+3);
Seting(result,2*4,c1);
c1=(vec4){0,0,0,0};
ShapeShifter=3<<2;
c1=TakeIting(c1,a,b,0,ShapeShifter+0);
c1=TakeIting(c1,a,b,4,ShapeShifter+1);
c1=TakeIting(c1,a,b,8,ShapeShifter+2);
c1=TakeIting(c1,a,b,12,ShapeShifter+3);
Seting(result,3*4,c1);
c1=(vec4){0,0,0,0};
}
int main()
{
mat4 a=(mat4){
2,9,40,5,
8,6,5,6,
8,9,7,4,
7,5,3,10
};
mat4 b=(mat4){
50,30,88,70,
85,100,0,10,
89,65,50,60,
99,45,14,80
};
struct timespec Ta,Tb;
float GlobalTime=0;
float localTime=0;
mat4 c;
clock_gettime(CLOCK_MONOTONIC,&Ta);
atTheTestMatMul(&c.v[0],&a.v[0],&b.v[0]);
clock_gettime(CLOCK_MONOTONIC,&Tb);
localTime=(float)(Tb.tv_nsec - Ta.tv_nsec) / 1000000 +(float)(Tb.tv_sec - Ta.tv_sec);
printf("%f\n",localTime);
// printf(
// "%f %f %f %f \n%f %f %f %f \n%f %f %f %f \n%f %f %f %f \n",
// c.v[0], c.v[1], c.v[2], c.v[3],
// c.v[4], c.v[5], c.v[6], c.v[7],
// c.v[8], c.v[9], c.v[10], c.v[11],
// c.v[12], c.v[13], c.v[14], c.v[15]
// );
return 0;
}
так самое лучшее вроде
Посмотрел на годболте что последний gcc выдаёт - это же ужас какой то, кроме одной реально нужной инструкции mulps/addps целая портянка кода. Всё таки ABI не заточено под передачу пользовательских структур через векторные регистры.
знаете что я заметил, я замерил более скурпулёзно с топиковым расчетом более тщательно на сколько это возможно, и расчет матриц с SSE не на много но быстрее, тоесть она уводит расчет в такой формат
0.000089 SSE медиана clock_gettime()
0.000130 - 0.000200 - медиана исходной матрицы
0.000100-0.000200 - медина расчета на указателях, на моём пк расчет указателей гонится в код SSE он буквально чуть чуть даст ускорения по отношению с исходной, но только SSE уведет его порой в скорость ниже 100
тоесть я запускал 1 расчет с 1 замером
а тот замер что я привел надо будет еще раз подумать там в нем чтото не так я ошибся с тем замером какой в статье указан
код на указателях повторяет логику SSE как видите он не уводит скорость ниже 100 а если уводит то реже чем SSE
верно - было замечено что ассемблер сделал 2 операции отсюда незамысловатые 4 строки превратились в 4 функции и там по-лучше чем в исходной ито ассемблером тоесть за счет того как ассемблер примет, видимо не предназначено, я тоже крайние разы смотрел и по годболту и на локальном
Скрытый текст
#pragma GCC optimize("Ofast,unroll-loops,-ffast-math")
#pragma GCC target("sse")
#include <xmmintrin.h>
#include <wmmintrin.h>
#include <smmintrin.h>
#include <nmmintrin.h>
#include <immintrin.h>
#include <stdio.h>
#include <unistd.h>
#include <time.h>
#include <errno.h>
#include <string.h>
#pragma pack(1)
typedef float f32;
struct vec4_t
{
f32 x;
f32 y;
f32 z;
f32 w;
};
typedef struct vec4_t vec4;
struct mat4_t
{
f32 v[16];
};
typedef struct mat4_t mat4;
vec4 Mulv4(vec4 a, vec4 b)
{
return (vec4){a.x * b.x, a.y * b.y, a.z * b.z, a.w * b.w};
}
vec4 Addv4(const vec4 a, const vec4 b)
{
return (vec4){a.x + b.x, a.y + b.y, a.z + b.z, a.w + b.w};
}
vec4 Processing(vec4 c1,vec4 a,vec4 b,int shift)
{
return Addv4(c1,Mulv4(a,b));
}
vec4 TakeIting(vec4 c1,const float *a,const float *b,int shiftA,int ShapeShifter)
{
return Processing(c1,(vec4){*(a+shiftA),*(a+shiftA+1),*(a+shiftA+2),*(a+shiftA+3)},(vec4){*(b+ShapeShifter),*(b+ShapeShifter),*(b+ShapeShifter),*(b+ShapeShifter)},ShapeShifter);
}
float* Seting(float *a,int shift,vec4 b)
{
*(a+shift)=b.x;
*(a+shift+1)=b.y;
*(a+shift+2)=b.z;
*(a+shift+3)=b.w;
return a;
}
vec4 GetingSum(vec4 c1,const float *a,const float *b,int ShapeShifter)
{
c1=TakeIting(c1,a,b,0<<2,ShapeShifter+0);
c1=TakeIting(c1,a,b,1<<2,ShapeShifter+1);
c1=TakeIting(c1,a,b,2<<2,ShapeShifter+2);
c1=TakeIting(c1,a,b,3<<2,ShapeShifter+3);
return c1;
}
vec4 Proc2(vec4 c1,float *result,const float *a,const float *b,int shiftA,int ShapeShifter)
{
c1=GetingSum((vec4){0,0,0,0},a,b,ShapeShifter);
Seting(result,shiftA,c1);
return c1;
}
void atTheTestMatMul(float *result,const float *a,const float *b)
{
vec4 c1;
c1=Proc2(c1,result,a,b,0<<2,0<<2);
c1=Proc2(c1,result,a,b,1<<2,1<<2);
c1=Proc2(c1,result,a,b,2<<2,2<<2);
c1=Proc2(c1,result,a,b,3<<2,3<<2);
}
mat4 Mulm4(const mat4 a, const mat4 b)
{
mat4 result;
result.v[0] = a.v[0] * b.v[0] + a.v[4] * b.v[1] + a.v[8] * b.v[2] + a.v[12] * b.v[3];
result.v[1] = a.v[1] * b.v[0] + a.v[5] * b.v[1] + a.v[9] * b.v[2] + a.v[13] * b.v[3];
result.v[2] = a.v[2] * b.v[0] + a.v[6] * b.v[1] + a.v[10] * b.v[2] + a.v[14] * b.v[3];
result.v[3] = a.v[3] * b.v[0] + a.v[7] * b.v[1] + a.v[11] * b.v[2] + a.v[15] * b.v[3];
result.v[4] = a.v[0] * b.v[4] + a.v[4] * b.v[5] + a.v[8] * b.v[6] + a.v[12] * b.v[7];
result.v[5] = a.v[1] * b.v[4] + a.v[5] * b.v[5] + a.v[9] * b.v[6] + a.v[13] * b.v[7];
result.v[6] = a.v[2] * b.v[4] + a.v[6] * b.v[5] + a.v[10] * b.v[6] + a.v[14] * b.v[7];
result.v[7] = a.v[3] * b.v[4] + a.v[7] * b.v[5] + a.v[11] * b.v[6] + a.v[15] * b.v[7];
result.v[8] = a.v[0] * b.v[8] + a.v[4] * b.v[9] + a.v[8] * b.v[10] + a.v[12] * b.v[11];
result.v[9] = a.v[1] * b.v[8] + a.v[5] * b.v[9] + a.v[9] * b.v[10] + a.v[13] * b.v[11];
result.v[10] = a.v[2] * b.v[8] + a.v[6] * b.v[9] + a.v[10] * b.v[10] + a.v[14] * b.v[11];
result.v[11] = a.v[3] * b.v[8] + a.v[7] * b.v[9] + a.v[11] * b.v[10] + a.v[15] * b.v[11];
result.v[12] = a.v[0] * b.v[12] + a.v[4] * b.v[13] + a.v[8] * b.v[14] + a.v[12] * b.v[15];
result.v[13] = a.v[1] * b.v[12] + a.v[5] * b.v[13] + a.v[9] * b.v[14] + a.v[13] * b.v[15];
result.v[14] = a.v[2] * b.v[12] + a.v[6] * b.v[13] + a.v[10] * b.v[14] + a.v[14] * b.v[15];
result.v[15] = a.v[3] * b.v[12] + a.v[7] * b.v[13] + a.v[11] * b.v[14] + a.v[15] * b.v[15];
return result;
}
void mulm4VVV(const float *result, const float *a, const float *b)
{
__m128 row0 = _mm_loadu_ps(&a[0]);
__m128 row1 = _mm_loadu_ps(&a[4]);
__m128 row2 = _mm_loadu_ps(&a[8]);
__m128 row3 = _mm_loadu_ps(&a[12]);
__m128 newRow0 = _mm_mul_ps(row0, _mm_set1_ps(b[0]));
newRow0 = _mm_add_ps(newRow0, _mm_mul_ps(row1, _mm_set1_ps(b[1])));
newRow0 = _mm_add_ps(newRow0, _mm_mul_ps(row2, _mm_set1_ps(b[2])));
newRow0 = _mm_add_ps(newRow0, _mm_mul_ps(row3, _mm_set1_ps(b[3])));
__m128 newRow1 = _mm_mul_ps(row0, _mm_set1_ps(b[4]));
newRow1 = _mm_add_ps(newRow1, _mm_mul_ps(row1, _mm_set1_ps(b[5])));
newRow1 = _mm_add_ps(newRow1, _mm_mul_ps(row2, _mm_set1_ps(b[6])));
newRow1 = _mm_add_ps(newRow1, _mm_mul_ps(row3, _mm_set1_ps(b[7])));
__m128 newRow2 = _mm_mul_ps(row0, _mm_set1_ps(b[8]));
newRow2 = _mm_add_ps(newRow2, _mm_mul_ps(row1, _mm_set1_ps(b[9])));
newRow2 = _mm_add_ps(newRow2, _mm_mul_ps(row2, _mm_set1_ps(b[10])));
newRow2 = _mm_add_ps(newRow2, _mm_mul_ps(row3, _mm_set1_ps(b[11])));
__m128 newRow3 = _mm_mul_ps(row0, _mm_set1_ps(b[12]));
newRow3 = _mm_add_ps(newRow3, _mm_mul_ps(row1, _mm_set1_ps(b[13])));
newRow3 = _mm_add_ps(newRow3, _mm_mul_ps(row2, _mm_set1_ps(b[14])));
newRow3 = _mm_add_ps(newRow3, _mm_mul_ps(row3, _mm_set1_ps(b[15])));
_mm_storeu_ps(&result[0], newRow0);
_mm_storeu_ps(&result[4], newRow1);
_mm_storeu_ps(&result[8], newRow2);
_mm_storeu_ps(&result[12], newRow3);
}
int main()
{
mat4 a=(mat4){
2,9,40,5,
8,6,5,6,
8,9,7,4,
7,5,3,10
};
mat4 b=(mat4){
50,30,88,70,
85,100,0,10,
89,65,50,60,
99,45,14,80
};
mat4 c;
struct timespec Ta,Tb;
float GlobalTime=0;
float localTime=0;
clock_gettime(CLOCK_MONOTONIC,&Ta);
//c=Mulm4(a,b);
//atTheTestMatMul(&c.v[0],&a.v[0],&b.v[0]);
//mulm4VVV(&c.v[0],&a.v[0],&b.v[0]);//moment time
clock_gettime(CLOCK_MONOTONIC,&Tb);
localTime=(float)(Tb.tv_nsec - Ta.tv_nsec) / 1000000 +(float)(Tb.tv_sec - Ta.tv_sec);
printf("%f\n",localTime);
printf(
"%f %f %f %f \n%f %f %f %f \n%f %f %f %f \n%f %f %f %f \n",
c.v[0], c.v[1], c.v[2], c.v[3],
c.v[4], c.v[5], c.v[6], c.v[7],
c.v[8], c.v[9], c.v[10], c.v[11],
c.v[12], c.v[13], c.v[14], c.v[15]
);
return 0;
}
вот на годболте например гарантировано только в 0.00080 удаленно сработать может только SSE. остальное зависимо от состояний
тоесть я запускал 1 расчет с 1 замером
Во первых, чтобы цифрам можно было доверять, надо бы завернуть код в цикл с выполнением хотя бы в секунду. Во вторых, важнее даже не насколько вы ускорили исходный код, а насколько приблизились к пиковой производительности на вашем процессоре - сможете её оценить теоретически?
по пику я упираюсь только в 3д по нагрузке если я отпускаю держатель времени, и он отсекается драйвером в 75 фпс, если включать мангохуд, нагрузка начинается только если сложная сцена и нет кешурующих нюансов, если всё более менее в пределах видимости будет всё ок с моей математикой, тут кватернионы 1 команда на ОСЬ-ориентир в нейже получается поворот камеры, все остальные обьекты по тому же принципу, ну и дальше умножения, ну и плохонький лист с ним я смотрю и туда внедрял в тот тест консольный с костями SSE и там я увидел как раз производительность, потомучто много расчетов, ну и в 3д внедрил, больше по упору нагрузки в предел ПК я пока ничего не могу сказать
всё равно так или иначе нюансы кеша, кватернион, и умножение матриц - пока это были ключевые моменты до сих пор в моих подходах, кватернионы решены, кеши это типо шрифт, умножение вот как раз будет SSE
эх на данный момент наверно пока не смогу
0,000080 дабы избежать путаницы в порядках
-Ofast -ffast-math -msse4.1 -msse4.2 -lm -lc -ldl
or
-Ofast -ffast-math -fopenmp -msse4.1 -msse4.2 -lm -lc -ldl
вобщем да, +- как я написал выше теперь и в годболте если все везде прописать
Всё таки - какой теоретический максимум и сколько процентов от него дают реализации на интринсиках и без?
по наблюдениям учитывая что системные вызовы тоесть наивное перемножение занимает общий ресурс, а SSE выделенный то стабильно процентов 30 даст, у SSE как я могу понять стабильный выделенный канал, вся система работает в другом канале, хотя словом канал пользоваться может быть не професионально
тоесть есть общие ресурсы, Юзерленд исполняется от ресурсов, латенси удается добится потомучто ресурс ктото даст использовать на время тоесть ресурс будет свободен, я пока так представляю, и вот ресурс SSE это как выделенный канал
и так как это ресурсы там происходят проверки так называемые зависимые запросы на состояния свободности и вот тут латенси наивных расчетов может скакать или быть локально быстрым а на другом пк по какойто причине иным, тоесть тут предельные значения, общего верхнего отклика на верхний запрос
Много слов и мало цифр )
понял, я покажу стабильную работу на фрибсд, а на линуксе может на цпу быть нагрузка
с цифрами да пока сложнее
Я ж не про нагрузку, а про то, сколько тактов процессора требует идеальная реализация умножения матриц и насколько близки к ней ваши варианты.
спасибо, пока я такого не знаю и пока не знаю как посмотреть кодом, всё что я знаю, только то что попробовал
я пользуюсь clock_gettime() пока только такой расчет у меня
CLOCK_MONOTONIC ближе всех поидее потомучто 90 тактов более реально на ССЕ чем 1 тысяча
я вот щас посмотрел CLOCK_PROCESS_CPUTIME_ID и там уже другие цифры и порядки
а вот по монотонику более понятно что первый самый медленный, второй чуть лучше это на указателях и тот который к 90 стремится это самый фаст, а тут не всё очевидно они почти все в 1 пределе такого не должно быть, по крайней мере я себе по другому представлял, поэтому я пока не знаю как циклы посчитать
если CLOCK_PROCESS_CPUTIME_ID это оно то от 483 до 900 = SSE. второе третье место остальные реализации при этом прок на ускорение редко появлялось я видел, но суть в том что ССЕ даёт стабильный доступ мне кажется, это ключевой момент, скорость должна быть стабильной
хотя я пока до конца не уверен с циклами процессора и временем
хм а я понял тогда скажу так
пилотная скорость у всех одинаковая
криты у SSE выше по циклам, тоесть она прокнуть может и до 483,
криты я редко видел на указателях там просто пилотный предел в пределе,
и бывает критует неоптимизированный код, но это реже, и пилотная скорость не такая как у SSE
но мне почемуто кажется CLOCK_PROCESS_CPUTIME_ID это не совсем то мне кажется там по другому меряют еще как-то
Можно теоретически посчитать сколько сложений и умножений нужно в функции, можно посмотреть сколько процессор по максимуму может исполнить за такт (с векторами подлиннее - у вас это скорее 256-битный AVX, и желательно используя FMA, а не отдельные ADD/MUL) - дальше поделить одно на другое.
если просто скомпилировать мою реализацию с mavx512f кроме msse... https://dpaste.com/ETX2EMBUT насчет подсчетов, ок попробую так, в 512f не запускается на моей тачке уже примтивный счет матрицы (Illegal instruction. dump не смотрел, скорее всего какаято инструкция не поддерживается)
C avx512 поосторожнее - на новых десктопных процессорах их вообще вырезали, на некоторых десктопных они были, но прироста скорости не давали. Реально стоит использовать на серверных, и то лучше Cascadelake или новее - на SKX неосторожным использованием 512-битных инструкций можно было замедлить код.

замер gettime(CLOCK_REALTIME)clock_gettime
1 запуск по 1000 счетов на 1 тип 3 запуска на 1 тип
выше представлено
1 SSE
2 указатели
3 примтивная где прописаны массивы.
по крайней мере у меня так сейчас
Всё таки в секундах сложно понимать, насколько это хорошо или плохо, проще в тактах на конкретный кусок ассемблерного кода.
хорошо посмотрю, спасибо
наивное avx2 будет таким после 128
Скрытый текст
#include <immintrin.h>
#include <avxintrin.h>
#include <time.h>
void mulm4VVV(float *result, const float *a, const float *b)//256
{
__m256 row0 = _mm256_loadu_ps(&a[0]);
__m256 row1 = _mm256_loadu_ps(&a[4]);
__m256 row2 = _mm256_loadu_ps(&a[8]);
__m256 row3 = _mm256_loadu_ps(&a[12]);
__m256 nR0= _mm256_mul_ps(row0,_mm256_set1_ps(b[0]));
nR0 = _mm256_add_ps(nR0,_mm256_mul_ps(row1,_mm256_set1_ps(b[1])));
nR0 = _mm256_add_ps(nR0,_mm256_mul_ps(row2,_mm256_set1_ps(b[2])));
nR0 = _mm256_add_ps(nR0,_mm256_mul_ps(row3,_mm256_set1_ps(b[3])));
__m256 nR1= _mm256_mul_ps(row0,_mm256_set1_ps(b[4]));
nR1 = _mm256_add_ps(nR1,_mm256_mul_ps(row1,_mm256_set1_ps(b[5])));
nR1 = _mm256_add_ps(nR1,_mm256_mul_ps(row2,_mm256_set1_ps(b[6])));
nR1 = _mm256_add_ps(nR1,_mm256_mul_ps(row3,_mm256_set1_ps(b[7])));
__m256 nR2= _mm256_mul_ps(row0,_mm256_set1_ps(b[8]));
nR2 = _mm256_add_ps(nR2,_mm256_mul_ps(row1,_mm256_set1_ps(b[9])));
nR2 = _mm256_add_ps(nR2,_mm256_mul_ps(row2,_mm256_set1_ps(b[10])));
nR2 = _mm256_add_ps(nR2,_mm256_mul_ps(row3,_mm256_set1_ps(b[11])));
__m256 nR3= _mm256_mul_ps(row0,_mm256_set1_ps(b[12]));
nR3 = _mm256_add_ps(nR3,_mm256_mul_ps(row1,_mm256_set1_ps(b[13])));
nR3 = _mm256_add_ps(nR3,_mm256_mul_ps(row2,_mm256_set1_ps(b[14])));
nR3 = _mm256_add_ps(nR3,_mm256_mul_ps(row3,_mm256_set1_ps(b[15])));
_mm256_storeu_ps(&result[0], nR0);
_mm256_storeu_ps(&result[4], nR1);
_mm256_storeu_ps(&result[8], nR2);
_mm256_storeu_ps(&result[12], nR3);
}
как измерять не знаю
минимальная вводная AVX
открыл файлик <avxintrin.h> там коментарии и описание
на указателях примерный алгоритм чтоб понять что происходит и почему не дот продакты, в кратце там суммирование тоесть именно приближено к тому интегралу какой в вики
4 суммы
можно со сдвигами делать команд ассемблера чуть побольше будет *(a+(0<<2))
простите но на С++ есть нюансы, я не говорю что bad, но я просто тонко скажу нюансы есть, на С не так работает( у меня не хватит технического опыта описать почему 3д на С и С++ по разному работает это империческое наблюдение 2 теста грубо говоря, хотя при определенном подходе на С++ будет удобнее и возможно так же или быстрее ) в 3д в общем масштабе так сказать, и добавлю всё проверяйте, я могу проверить только на своей тачке и тех парах штук которые у меня имеются, еще добавлю просто в самой математике тоже есть нюансы, ктото пользуется кватернионами а ктото поворотными матрицами, в блендере например на сколько я могу понять теперь кватернионы

на Расте будет либо так же по скорости либо быстрее аналогично, судя по моим наблюдениям наблюдал в пару видео, лично дел не имел
tan tan показал свой опыт с 3д на расте и рассказал о своём пути
так же есть разница в платформах, (виндовс/linux/unix/other)+gpu+cpu, на линуксе наивная анимация нагружает(25 бон, гуманоид, тривиальная математика), когда на юниксе уже другое отобразится в диспетчере и наваливать не будет имейте ввиду(на юниксе состояния процессора фиксируются иным образом, другая архитектура, на линуксе, не такая архитектура как на Юниксе нету селекта и прочее и рендер процесса(контекст составляющей ОС не 3д(процессы, архитектура вот это вот всё)) другой ), на виндовсе не смотрел с текущими знаниями, которые приобрёл по пути вникания в процесс математика+3д
и даже тут есть нюансы - конфигурации ПК, вот, но в целом вот как есть :( потомучто нагружать может на линуксе только на моём ПК и только на том дистрибутиве на котором я запускал-смотрел
простите но на С++ есть нюансы,
Ну при чём тут плюсплюсы-то? Я Вас уверяю, интринсики будут одинаковы что там что на Си. Вот, смотрите, код из пятой главы - это ровно тоже что и у Вас, только с FMA до кучи. Тут как раз нет вообще никаких нюансов.
на Расте будет либо так же по скорости либо быстрее ... так же есть разница в платформах,...
Вы мешаете вообще всё в одну кучу.
да смешал, всё в одну кучу, чтобы вот сейчас ответить, конкретно инстрисинков будет одинаково, но фон общего масштаба как рисуется и прочее будет разный, просто я когда открыл вашу ссылку обратил внимание на С++, и решил уточнить, эх код как там, я не смотрел никуда, только описание математики и применил наблюдение из воксельной составляющей, я лично только листал книжку совсем другую
Computer Graphics, C Version (2nd Edition) я её смотрел
Computation Geometry in C а эту щас нашел пока очередной раз искал верхнюю
советую, тут мало возможно того что сегодня, но касательно 3д имеет место быть я точно edtition не скажу, надоело искать по названию
(просто бесчетное количество попыток и гайдов разного рода по типо learnopengl)
точка 1 воксельная часть хранится в упаковке(сетка или точки(минимальный воксель это наличие в точке вокселя)) плоской и всё 3д пространство мы прыгаем смещениями по упаковке получая нужный кубик в каждом кубике стороны итдитп, на ютубе есть конференция по воксельной составляющей и сотни тысяч примеров как делали можно посмотреть, там математика интересная и увлекательная,
точка 2 делаем ради теста замер математика плюс такое пространство далее делаем выводы или не делаем
сравниваем если требуется С++ и С
(!) а этот подход у меня очередной в нём получилось пройти дальше, отладить математику(избавится от глм обновлений), привести всё в рабочее состояние, сгенерировать поверхности, выйти на новый уровень - тоесть оперировать более большими частями пространства поверхностей и они стали не маленькими, в общем, базовое представление 3д, так же понял свет, и прочее (сравниваем)
(зная из первых тестов, которые были нечаяно найдены имперически тоесть я наблюдал еффективность работы обоих на своём пк, я сначала прикинул как бы лучше и посмотрел опять С++ и зная С и как работает ушел просто в С)
как видно за 1 день не догодаться!
до FMA надо будет как раз входить уже в тонкости, на простейшем этапе через воксельное понимание, я первым делом смотрел что у меня работает и первая рекомендация по ускорению встретится не FMA а SSE далее уже тонкости, которые надо будет дотягивать теорией и описаниями технического характера, еще раз скажу конкретно в такой последовательности даже ускорения х10 достаточно, остальное придётся разбираться
вот тут можно почитать
what-are-the-benefits-of-quaternion-rotation-over-eulerian-rotation
и в вики много пояснений
del
округление опасно FreeBSD предлагает округлить лично у меня на одном из тестов когда начал править по ворнингам всё упало(приходилось откатывать), всё проверять надо, эта реализация с linux(сначала я тестил на Линуксе и думал что на Линуксе буду), и соблюдает все расчеты по стабильности в т.ч. того оптимального достаточного округления
(какоето из sqrtf или fabs может положить расчеты)
Вы точно этот код запускали и он точно вам один и тот же ответ выдает с наивным методом?
Я не поленился проверить, да, он выдаёт тот же результат, что и референсный dgemm(). Правда я перебросил его на double:
Скрытый текст
MATMULDLL_API void mulm4VVVd(double* result, double* a, double* b)
{
__m256d row0 = _mm256_loadu_pd(&b[0]);
__m256d row1 = _mm256_loadu_pd(&b[4]);
__m256d row2 = _mm256_loadu_pd(&b[8]);
__m256d row3 = _mm256_loadu_pd(&b[12]);
__m256d newRow0 = _mm256_mul_pd(row0, _mm256_set1_pd(a[0]));
newRow0 = _mm256_add_pd(newRow0, _mm256_mul_pd(row1, _mm256_set1_pd(a[1])));
newRow0 = _mm256_add_pd(newRow0, _mm256_mul_pd(row2, _mm256_set1_pd(a[2])));
newRow0 = _mm256_add_pd(newRow0, _mm256_mul_pd(row3, _mm256_set1_pd(a[3])));
__m256d newRow1 = _mm256_mul_pd(row0, _mm256_set1_pd(a[4]));
newRow1 = _mm256_add_pd(newRow1, _mm256_mul_pd(row1, _mm256_set1_pd(a[5])));
newRow1 = _mm256_add_pd(newRow1, _mm256_mul_pd(row2, _mm256_set1_pd(a[6])));
newRow1 = _mm256_add_pd(newRow1, _mm256_mul_pd(row3, _mm256_set1_pd(a[7])));
__m256d newRow2 = _mm256_mul_pd(row0, _mm256_set1_pd(a[8]));
newRow2 = _mm256_add_pd(newRow2, _mm256_mul_pd(row1, _mm256_set1_pd(a[9])));
newRow2 = _mm256_add_pd(newRow2, _mm256_mul_pd(row2, _mm256_set1_pd(a[10])));
newRow2 = _mm256_add_pd(newRow2, _mm256_mul_pd(row3, _mm256_set1_pd(a[11])));
__m256d newRow3 = _mm256_mul_pd(row0, _mm256_set1_pd(a[12]));
newRow3 = _mm256_add_pd(newRow3, _mm256_mul_pd(row1, _mm256_set1_pd(a[13])));
newRow3 = _mm256_add_pd(newRow3, _mm256_mul_pd(row2, _mm256_set1_pd(a[14])));
newRow3 = _mm256_add_pd(newRow3, _mm256_mul_pd(row3, _mm256_set1_pd(a[15])));
_mm256_storeu_pd(&result[0], newRow0);
_mm256_storeu_pd(&result[4], newRow1);
_mm256_storeu_pd(&result[8], newRow2);
_mm256_storeu_pd(&result[12], newRow3);
}
И он, кстати, в данном частном случае 4х4 обгоняет по производительности реализацию LabVIEW. Ну, может и пригодится когда-нибудь. Хотя сейчас такой код практически сразу ИИ сходу набрасывает.
это всё интересно если вы из другой области, автор этим пользуется в 3д
она может это набрасывать, потомучто есть определенные закономерности плоских наблюдений на сколько мне известно, посмотрите внимательно на Воксели и его математику, на оптимизации, greedy_mesh, и вторая базовая где отрисовкой играемся ради производительности типо meshing
так же есть какието библиотеки с flat историей я не смотрел что это такое
тут решается вопрос, когда уже решен главный момент, как видите нету шифта, транспонирования тут нету, это решение общее для трансформации, не тратье время на шифт, на первых приближениях изучения оно не нужно, и если есть SSE ассемблер сам как видите сделает шифт, а второй способ, делает тоже самое без единого цикла в том что я показал
немного конкретики если вдаваться в подробности есть еще третий(и еще 1) способ решения,
Скрытый текст

нас интерересует вторая матрица, только тогда у меня встречный вопрос, а зачем это надо если итак уже это делается без цикла, чтобы доказать что библиотека сама себя посчитать может? так может
напомню https://en.wikipedia.org/wiki/Matrix_multiplication мы можем идти в цикле, а можем не в цикле!
правильность делания в цикле или без я не знаю, мне без цикла пока нравится
итого выходит около 5 решений, ускоренное, то которое я показал исходное, наивное сложение + шифт, и дот продакты, тут извините, формула указана, там решение в любом гайде в цикле бежим считаем, вобщем, да, вот так, еще раз извините, под трансформы и просто рисовать можно не ломать голову разными решениями
Скрытый текст
void atTheTest()
{
vec4 t5=(vec4){0.0f,0.0f,0.0f,0.0f};
vec4 a1=(vec4){2,8,8,7};vec4 b1=(vec4){50 ,30,88,70};
vec4 a2=(vec4){9,6,9,5};vec4 b2=(vec4){85,100,0,10};
vec4 a3=(vec4){40,5,7,3};vec4 b3=(vec4){89,65,50,60};
vec4 a4=(vec4){5,6,4,10};vec4 b4=(vec4){99,45,14,80};
vec4 c1;
c1=Mulv4(a1,(vec4){b1.x,b1.y,b1.z,b1.w});
vec4 c2;
c1=Addv4(c1,Mulv4(a2,(vec4){b2.x,b2.y,b2.z,b2.w}));
vec4 c3;
c1=Addv4(c1,Mulv4(a3,(vec4){b3.x,b3.y,b3.z,b3.w}));
vec4 c4;
c1=Addv4(c1,Mulv4(a4,(vec4){b4.x,b4.y,b4.z,b4.w}));
printf("%f %f %f %f\n",Dotv4(a1,b1),Dotv4(a2,b1),Dotv4(a3,b1),Dotv4(a4,b1));
printf("%f %f %f %f\n",Dotv4(a1,b2),Dotv4(a2,b2),Dotv4(a3,b2),Dotv4(a4,b2));
printf("%f %f %f %f\n",Dotv4(a1,b3),Dotv4(a2,b3),Dotv4(a3,b3),Dotv4(a4,b3));
printf("%f %f %f %f\n",Dotv4(a1,b4),Dotv4(a2,b4),Dotv4(a3,b4),Dotv4(a4,b3));
}
проблема в том что я храню значения матриц не в векторах, вот, а так самопроверка работает.
а то что указано выше в примере то что считается не вывода, для него надо наоборот проходить по матрице чтобы попасть в текущую ситуацию счета, и там нужен шафл или шифт(по игрику матрицы), как раз то что делает ассемблер (наоборот это значит в матрице а строка, а в матрице б колонка)
всё теперь никто не запутается
Скрытый текст
void atTheTest2()//0,1,2,3
{
vec4 a1=(vec4){2,9,40,5};vec4 b1=(vec4){50 ,85,89,99};
vec4 a2=(vec4){8,6,5,6};vec4 b2=(vec4){30,100,65,45};
vec4 a3=(vec4){8,9,7,4};vec4 b3=(vec4){88,0,50,14};
vec4 a4=(vec4){7,5,3,10};vec4 b4=(vec4){70,10,60,80};
//Mulv4 умножение векторов add сложение
vec4 c1; //x,y,z,w
c1=Mulv4( a1,(vec4){b1.z,b1.z,b1.z,b1.z}); //a row
c1=Addv4(c1,Mulv4(a2,(vec4){b2.z,b2.z,b2.z,b2.z}));
c1=Addv4(c1,Mulv4(a3,(vec4){b3.z,b3.z,b3.z,b3.z}));
c1=Addv4(c1,Mulv4(a4,(vec4){b4.z,b4.z,b4.z,b4.z}));
/* printf("%f %f %f %f\n",Dotv4(a1,b1),Dotv4(a2,b1),Dotv4(a3,b1),Dotv4(a4,b1)); */ == /* result.v[0] = a.v[0] * b.v[0] + a.v[4] * b.v[1] + a.v[8] * b.v[2] + a.v[12] * b.v[3]; */
/* printf("%f %f %f %f\n",Dotv4(a1,b2),Dotv4(a2,b2),Dotv4(a3,b2),Dotv4(a4,b2)); */
/* printf("%f %f %f %f\n",Dotv4(a1,b3),Dotv4(a2,b3),Dotv4(a3,b3),Dotv4(a4,b3)); */
/* printf("%f %f %f %f\n",Dotv4(a1,b4),Dotv4(a2,b4),Dotv4(a3,b4),Dotv4(a4,b3)); */
/* result.v[0] = a.v[0] * b.v[0] + a.v[4] * b.v[1] + a.v[8] * b.v[2] + a.v[12] * b.v[3]; *///a collumn
/* result.v[1] = a.v[1] * b.v[0] + a.v[5] * b.v[1] + a.v[9] * b.v[2] + a.v[13] * b.v[3]; */
/* result.v[2] = a.v[2] * b.v[0] + a.v[6] * b.v[1] + a.v[10] * b.v[2] + a.v[14] * b.v[3]; */
/* result.v[3] = a.v[3] * b.v[0] + a.v[7] * b.v[1] + a.v[11] * b.v[2] + a.v[15] * b.v[3]; */
printf("%f %f %f %f\n",c1.x,c1.y,c1.z,c1.w);
}
заключительный пример тут симулируется шафл/шифт
ситуация следующая, то что b его компонента шифтуется вот теперь всё, для такого расчета используется представленная в примере матрица а и б (и походу интринсик так и работает, наверно)
аффтар. ты ффто куришь?
извините, в каком смысле
Похоже, это тот же самый автор, что и Трампа хотел замочить.
Обращаюсь к комментаторам этой ветки - ваши подколки вас совершенно не красят.
Автор искренне увлекается компьютерной графикой и по его постам видно, что не всегда у него получается в слаженное повествование и его может занести. Но это не школьник/студент, рекламирующий с помощью LLM свой телеграм-канал.
Обратите внимание, что остальные пользователи дают советы и замечания по теме (причем, без перехода на личности), что полезно как автору, так и остальным читателям.
А почему "граффика" с двумя "ф". Это такой авторский стиль?
Зачем вводится f32 и потом в коде мешается с float?
Зачем вводится mat4_t и mat4?
Почему getN() возвращает int и потом везде кастится во float (при касте в mat4)?
Зачем localTime объявлено как float, а справа выражение double?
Зачем там же два раза кастить в double, когда достаточно добавить ".0" к миллиону?
Почему в некоторых функциях входные float параметры передаются по указателю?
Ну и названия переменных и функций - просто пример, как не надо называть.
И суффикс _t
зарезервирован, и сторонний код не должен использовать его.
Ах, это ещё и для стороннего пользования... А почему только t? Почему бы не зарезервировать остальные буквы? А что насчёт остальных моих пунктов? А вы знаете, что float не гарантирует 4 байта?
Ты чё разошёлся, утро уже, не проспался?
библиотека реализована в рамках 3д, а пример же работает
float не гарантирует 4 байта(потомучто плавающая точка?), что тогда делать? мои примеры по 3д и с которыми я имел дело, все на флоат, в шейдере тоже флоат, кроме где не флоат. (вектора флоат, матрицы флоат)
если читать сериализацию из блендера вам поможет только сдвиг 4 байта( в части флоат ), + от платформы зависит, соответственно следующие вводные данные в матрицах ( они будут тоже флоат ) View Projection, соотв флоат измерение соотв вектора флоат и пошел расчет моделей, а если нету сериализации а до этого всё было определено как флоат то будет работать, потомучто всё на флоат, более тонкий момент не гарантирования 4 байтов, не рассматривался мной
(я знаю что в этих темах есть много вопросов на которые нету ответов, и тем неменее с присутствием этих вопросов у меня работает
https://stackoverflow.com/questions/8746443/float-and-int-both-4-bytes-how-come
)
у меня даже высота расчитывается, по IEEE 754 этой точности и тоже флоат, там есть нюансы но вы просто поймите как расчеты производить(простой 3д на простом мат аппарате сегодня укладывается в точность, раньше были нюансы, сейчас не знаю как, "но проблем как раньше не встречал на всём пути - построения мат аппарата с типом флоат, сегодняшней реализации С а не например какогото прошедшего года"), а вот раньше да у меня на расчете высоты вылетало, откуда это берётся я не знаю, откуда взялся этот вопрос тоже, с сериализацией лично имел дело в наши дни, а по мат аппарату наивно просто прошелся не встретил вашего вопроса, ни на Vulkan ни на OpenGL, было бы интересно увидеть примеры, что вы имеете ввиду, я бы смог дать более детальный ответ

в .3ds формате столкнулся, понятно что можно записать скриптом модельку в свой формат
потомучто мой вектор так реализован на момент написания статьи и при тестах я не во все добавления ввёл свой тип
потомучто mat_4 и mat4 у меня так используеются и определен mat4(соотв матрица у меня определенного вида ) как тип я пришел к такому соглашению лично
потомучто в примере rand()
другие замеры не работали с там есть пример функции с clock_t
чтобы передать адрес первого элемента(я пользуюсь плоским представлением матрицы если вы заметили, это было проверено вообще в игре, типо кубик крафт ) при вызове функции
название удобное у меня другие названия и такое название выбивает эти определения из общей библиотеки
Не знаю, как вы ухитрились получить 10 раз, возможно, компилятор не понял Ofast и не оптимизировал.
Я бы предположил, что большого ускорения не будет, потому что умножение матриц компиляторы давно векторизуют сами, причём под любую систему команд, включая FMA и AVX512, про которые выше писали. У gcc прям красивый ассемблер получается, clang почему-то чуть хуже: https://gcc.godbolt.org/z/91b9M4GeM
Правда, пол-функции - это перетасовки, чтобы использовать полную длину AVX-регистров. Не факт, что это быстро на самом деле, ну напишите -msse4 вместо -march=haswell, будет примерно копия ваших sse-интринсиков.
Попиарю тогда уж и свою статью, там как раз про затачивание кода под автовекторизатор:
https://habr.com/ru/articles/685228/
не не не, гцц имеет место, в линуксе гцц и наваливает на анимациях, ок можно предположить код не ок, но на кланге на фряхе не наваливает, далее в этой статье не векторизированный тоже векторизируется только компилятором и как бы прирост есть(SSE мультик матриц) даже в самом 3д
Насчет AVX не знаю не пользовался , насчет ФМА пока спорно ибо он может выйти боком(конкретно в таких расчетах, в каких-то других проектах где не будет ничего падать он может быть будет классный)
ФМА придется плотно тестировать в целом округление не нужно и абсолютность пока тоже не нужна какая-то сверх той какая уже есть
по-сути надо просто операции перекинуть в SSE/AVX но не стараться с округлениями наверно, просто прокинуть операции в такой вид, и я думаю этого будет достаточно на первом приближении (мультик матриц интересен в любом ускорении, дело в том что помимо позиционирования он участвует очень активно в анимации скелетной, там может быть много вызовов), и тут SSE даст мягкость, а не ускоренный как-либо мультик матриц просто покажет нагрузку если убрать чуть задежку отпустить
вот посмотрите
Скрытый текст
void atTheTestMatMul(const float *a,const float *b)
{
/* vec4 a1=(vec4){2,9,40,5};vec4 b1=(vec4){50 ,85,89,99}; */
/* vec4 a2=(vec4){8,6,5,6};vec4 b2=(vec4){30,100,65,45}; */
/* vec4 a3=(vec4){8,9,7,4};vec4 b3=(vec4){88,0,50,14}; */
/* vec4 a4=(vec4){7,5,3,10};vec4 b4=(vec4){70,10,60,80}; */
vec4 c1;
for(int i=0;i<4;++i)
{
c1=(vec4){0,0,0,0};
float r=*(b+i);
c1=Mulv4((vec4){*(a),*(a+1),*(a+2),*(a+3)}, (vec4){r,r,r,r});
r=*(b+i)+4;
c1=Addv4(c1,Mulv4((vec4){*(a+4),*(a+5),*(a+6),*(a+7)}, (vec4){r,r,r,r}));
r=*(b+i)+8;
c1=Addv4(c1,Mulv4((vec4){*(a+8),*(a+9),*(a+10),*(a+11)}, (vec4){r,r,r,r}));
r=*(b+i)+12;
c1=Addv4(c1,Mulv4((vec4){*(a+12),*(a+13),*(a+14),*(a+15)}, (vec4){r,r,r,r}));
printf("%f %f %f %f\n",c1.x,c1.y,c1.z,c1.w);
}
}
тут пока ошибка сделаете шифт правильный будет почти так же как в статье, только надо будет смотреть асемблер как не крути так как в примере на скринах из емакса используются другие регистры
Начальное ускорение математики