【发布时间】:2019-10-31 16:00:34
【问题描述】:
我已经实现了一个基于向量的 C# Log 近似值。它包括不安全的代码。它在许多环境中运行良好,但在最近的部署中失败了。该实现通过 System.Numerics.Vectors 库使用 SIMD。
很遗憾,我无法在软件无法运行的系统上进行测试。但是,我想知道我对图书馆所做的哪些假设是无效的:
- Vector
.Count 是否总是返回 2 的幂? - Vector
.Count == Vector.Count * 2 吗? - 我是否可以使用 Unsafe.AsPointer 获取向量的指针,然后执行标准操作,就好像它是内存中的 N 个压缩数字一样?
- 是否有运行 dotNet 4 的处理器具有不同的字节序或不以 IEEE754 格式存储浮点数?
代码如下:
const float invLn2 = 1.44269504089f; // 1 / ln(2)
const float pow2_126 = 8.5070592e+37f; //2^126
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public static Vector<float> QuickLog2(Vector<float> vecOrig)
{
//32 bit Float specification:
//Leftmost bit is sign bit.
//Next 8 bits are exponent
//Next 23 bits are mantissa
unsafe
{
var ints = Vector.AsVectorUInt32(vecOrig);
var exponents = Vector.BitwiseAnd(ints, new Vector<uint>(0x7F800000));
BitshiftVector23(Unsafe.AsPointer(ref exponents));
var unsignedExponents = exponents - new Vector<uint>(127);
var signedExponents = Vector.AsVectorInt32(unsignedExponents);
var localMantissBitmask = Vector.AsVectorSingle(new Vector<UInt32>(0x807FFFFF));
var maskedMantissas = Vector.BitwiseAnd(vecOrig, localMantissBitmask);
var mantissas = maskedMantissas * new Vector<float>(pow2_126);
var mantissasLogged = LogPolynomialFunction2(mantissas) * new Vector<float>(invLn2);
Vector<float> floatExponents;
#if false
floatExponents = Vector.ConvertToSingle(signedExponents);
#else
ConvertIntToFloatInPace(Unsafe.AsPointer(ref signedExponents));
floatExponents = Vector.AsVectorSingle(signedExponents);
#endif
return mantissasLogged + floatExponents;
}
}
const float log10_2 = 0.30102999566398119521373889472449f;
/// <summary>
/// A vectorized implementation of Log10(N). Uses bitshift, bitmasks, and unsafe code.
/// Does not have the same safety as Math.Log10: Behaviour for infities, zero, negative numbers are undefined.
/// </summary>
/// <param name="vec">The vector to take the log of</param>
/// <returns>The log, to the base 10, of the vector</returns>
/// <remarks>
/// Accurate to about 10^-7, which is the limit of a 32 bit float anyway.
/// In my (BS) tests, takes about twice as long to run on as Math.Log10(...), but operates on 8 numbers,
/// so 4x faster.
/// Reverts to Math.Log10(...) if vectors are not hardware accelerated.
/// Given the extra memory copies required, that will be much slower than using scalar code.
/// It'll be nice once intrinsics make it into dotNet and we can replace this with a single instruction...
/// </remarks>
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public static Vector<float> QuickLog10(Vector<float> vec)
{
if (Vector.IsHardwareAccelerated)
return QuickLog2(vec) * new Vector<float>(log10_2);
else
{
float[] tmp = new float[Vector<float>.Count];
vec.CopyTo(tmp);
for (int i = 0; i < Vector<float>.Count; i++)
tmp[i] = (float)Math.Log10(tmp[i]);
return new Vector<float>(tmp);
}
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static unsafe void BitshiftVector23(void* vector)
{
UInt64* asUlong = (UInt64*)vector;
if (Vector<UInt64>.Count == 4)
{
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
}
else if (Vector<UInt64>.Count == 8)
{
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
asUlong++;
*asUlong = *asUlong >> 23;
}
else
for (int i = 0; i < Vector<UInt64>.Count; i++)
asUlong[i] = asUlong[i] >> 23;
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static unsafe void ConvertIntToFloatInPace(void* vector)
{
int* asInt = (int*)vector;
if (Vector<int>.Count == 8)
{
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
}
else if (Vector<UInt64>.Count == 16)
{
for (int i = 0; i < 2; i++)
{
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
*(float*)asInt = *asInt;
asInt++;
}
}
else
for (int i = 0; i < Vector<UInt64>.Count; i++)
{
*(float*)asInt = *asInt;
asInt++;
}
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static Vector<float> LogPolynomialFunction2(Vector<float> mantissas)
{
var zm1 = mantissas;
var zp1 = mantissas + new Vector<float>(2f);
var zm1Divzp1 = Vector.Divide(zm1, zp1);
var squared = zm1Divzp1 * zm1Divzp1;
var cur = zm1Divzp1;
//Manual loop unwinding:
#if false
var mantissasLogged = Vector<float>.Zero;
for (float i = 0; i < 4; i++)
{
var fac = 2f / (2f * i + 1f);
mantissasLogged += cur * new Vector<float>(fac);
cur *= squared;
}
#else
//i = 0;
const float fac0 = 2f / (2 * 0 + 1);
var mantissasLogged = cur * new Vector<float>(fac0);
cur *= squared;
//i = 1;
const float fac1 = 2f / (2 * 1 + 1);
mantissasLogged += cur * new Vector<float>(fac1);
cur *= squared;
//i = 2;
const float fac2 = 2f / (2 * 2 + 1);
mantissasLogged += cur * new Vector<float>(fac2);
cur *= squared;
//i = 3;
const float fac3 = 2f / (2 * 3 + 1);
mantissasLogged += cur * new Vector<float>(fac3);
cur *= squared;
//i = 4;
const float fac4 = 2f / (2 * 4 + 1);
mantissasLogged += cur * new Vector<float>(fac4);
#endif
return mantissasLogged;
}
编辑:我在启动时对程序进行了一些简单的测试。 Vector.IsHardwareAccelerated == true; Vector.Count == 4;这个向量化的 Log 给出了前两个输入的正确答案,但后两个输入不正确。也许 Unsafe.AsPointer(Vector) 给我一个指向向量元素的指针作为四个连续浮点数的假设是不正确的。
日志输出:
DEBUG Vector.IsHardwareAccelerated: True
DEBUG Vector<float>.Count: 4
DEBUG Vector<Uint64>.Count: 2
DEBUG MathUtils test input data: 5.967E+009,1.072E+006,9.521E+017,4.726E+000
DEBUG MathUtils required output: 9.776,6.030,17.979,0.674
DEBUG MathUtils actual output: 9.776,6.030,0.218,0.072
(还有机会检查位模式...)
【问题讨论】:
-
如果目标支持 AVX-512 那么它就不可能是一个奇怪的架构
-
@harold:与 AArch64 SVE 或类似 Cray 风格的 vector machine 相比,您将如何具体检测 AVX512?即使 C#
float始终是 IEEE binary32,在理论上,字节序仍然可能是未来架构/未来 C# 实现的一个问题。 -
@PeterCordes 代码给出了错误的结果现在,假设的期货并没有导致这个问题
-
@harold:哦,我明白了,您说的是导致 OP 提出这个一般性问题的具体问题。一般不是标题问题。
标签: c# floating-point endianness simd unsafe