Hacker News Reader

Comment by AnimalMuppet | Hacker News Reader

The argument seems to be that C isn't suited to SIMD because some functions (exp, for instance) are not. But it then talks about how exp is implemented in hardware, in a way that forces it to be a non-SIMD calculation. That's not a C issue, that's a CPU issue. No language can save you there (unless it uses a software-only function for exp).

Is this article really confused, or did I misunderstand it?

The thing that makes C/C++ a good language for SIMD is how easily it lets you control memory alignment.

View on HN

Not in the article, but I've read that the way C does pointers, a single address with length implicit, makes it hard for compilers to assert that 2 arrays don't overlap/alias and this is an obstacle for generating SIMD instructions.

Not really. Pointers from different allocations are guaranteed to never alias. If they do, it is UB.

Right, but the compiler in general has no way to know that pointers are from different allocations - especially several function calls away from the allocations, with some pointer arithmetic in between.

Yes, that is why restrict was introduced.

In general, the information as to whether or not the pointers will alias might not even be known, right? Like we could be compiling a library that will take in some pointers from the user.

In C++ there is the align_alloc specificier. https://en.cppreference.com/w/c/memory/aligned_alloc

Not sure for C

aligned_alloc is a C function. It doesn't help the compiler prove that two pointers can't alias. restrict is the keyword that does.

The compiler is allowed to assume the results of malloc and certain other allocation functions do not alias.

Restrict does too.

I thought you could simply add `restrict`[1] to indicate that only one pointer points to a specific object? Shouldn't this help?

[1]https://en.cppreference.com/w/c/language/restrict

You can add that qualifier, but nothing will enforce that objects don't overlap. So you replaced slow correct code with fast incorrect code.

If you believed the pointers didn't alias and they actually did then you replaced slow incorrect code with fast incorrect code

That doesn't follow? If you replaced every call of memcpy with memmove (to use an example from the standard library), then your program is no less correct than it was before (and possibly more). The converse is that adding "restrict" qualifiers to existing functions can only make your program less correct, not more.

In this context it follows. If you believe the pointers dont alias (and want to use simd) then there should be no data-dependencies inbetween indices. (If there are data-dependencies it is quite likely that you have designed your algorithm to accomodate that, i.e. make it suitable for vectorised access. parallel prefix-sum would be an example of this). And without data-dependencies vectorized code is exactly equal to unvectorized code in all circumstances. If however, you do have data-dependencies, then your algorithm is wrong. Vectorizing can surface this error, while keeping it unvectorized may hide it. In my previous example (parallel prefix-sum) if you had a direct data-dependency to the previous index your program would've been wrong - regardless of vectorization.

In short, while your statement is true in general, I believe that it is not applicable in the context of the discussion.

That's what the restrict type qualifier[0] is for.

[0] https://en.cppreference.com/w/c/language/restrict

When compiling a library, it is generally impossible for the compiler to know whether or not the pointers will be aliased, right? That decision is made after the library is already compiled.

The function declaration in the header file of the library can carry the required 'restrict'. It works for c++ invokers too, as most c++ compilers also support and check the __restrict__ for old plain pointer types.

> But it then talks about how exp is implemented in hardware, in a way that forces it to be a non-SIMD calculation.

The fact that exp is implemented in hardware is not the argument. The argument is that exp is a library function, compiled separately, and thus the compiler cannot inline the function and fuse it with an array-wide loop, to detect later on an opportunity to generate the SIMD instructions.

It is true however that exp is implemented in hardware in the X86 world, and to be fair, perhaps a C compiler only needs to represent that function as an intrinsic instead, to give itself a chance to later replace it either with the function call or some SIMD instruction; but I guess that the standard doesn't provide that possibility?

The compiler can (and does) know enough about the stdlib functions to treat them specially. The C standard doesn't require that libm exist as a separately linked library either, that's just an implementation detail from an earlier time that got enshrined in POSIX.

> perhaps a C compiler only needs to represent that function as an intrinsic instead, to give itself a chance to later replace it either with the function call or some SIMD instruction

I gesture to this in the blog post:

> In C, one writes a function, and it is exported in an object file. To appreciate why this is special, consider sum :: Num a => [a] -> a in Haskell. This function exists only in the context of GHC. > ... > perhaps there are more fluent methods for compilation (and better data structures for export à la object files).

I now realize that I somewhat rephrased what was said in the blog. How strange!

> But it then talks about how exp is implemented in hardware, in a way that forces it to be a non-SIMD calculation.

This is especially surprising given that very few architectures have an `exp` in hardware. It's almost always done in software.

I think the interpretation is more conceptual about functions, reusability and composability.

The function exp has an abstract mathematical definition as mapping from numbers to numbers, and you could implement that generically if the language allows for it, but in C you cannot because it's bound to the signature `double exp(double)` and fixed like that at compile time. You cannot use this function in a different context and pass e.g. __m256d.

Basically because C does not have function dispatch, it's ill-suited for generic programming. You cannot write an algorithm on scalar, and now pass arrays instead.

> You cannot write an algorithm on scalar, and now pass arrays instead.

It is limited to cases where you know all overloads beforehand and limited because of C’s weak type system (can’t have an overload for arrays of int and one for pointers to int, for example), but you can. https://en.cppreference.com/w/c/language/generic:

  #include <math.h>
  #include <stdio.h>
 
  // Possible implementation of the tgmath.h macro cbrt
  #define cbrt(X) _Generic((X),     \
              long double: cbrtl, \
                  default: cbrt,  \
                    float: cbrtf  \
              )(X)
 
  int main(void)
  {
    double x = 8.0;
    const float y = 3.375;
    printf("cbrt(8.0) = %f\n", cbrt(x));    // selects the default cbrt
    printf("cbrtf(3.375) = %f\n", cbrt(y)); // converts const float to float,
                                            // then selects cbrtf
  }

It also, IMO, is a bit ugly, but that fits in nice with the rest of C, as seen through modern eyes.

There also is a drawback for lots of cases that you cannot overload operators. Implementing vector addition as ‘+’ isn’t possible, for example.

Many compilers partially support that as a language extension, though, see for example https://clang.llvm.org/docs/LanguageExtensions.html#vectors-....

I am not sure what you mean by "weak" type system. C's type system is not weak IMHO and also differentiates between arrays and pointers:

   _Generic((typeof(x), int*: 1, int[]: 2)

\_Generic with type argument is a new feature, it is a bit less elegant without it:

   _Generic((typeof(x)*)0, typeof(int*)*: 1, typeof(int[])\*: 2)

Any pointer converts to void* without warning. That's pretty weak.

This is the point of void*, one does not have to use it.

void* malloc(size_t size);

void free(void *ptr);

You use it a lot more than you think ;)

This is not a type safety issue though. For malloc the size has to agree. But in any case, you could wrap it into type-aware allocation demonstrating that this is not a limitation of the type system.

Perhaps not in the sense of weak vs strong, but rather in weak vs powerful. C type system isn't very powerful, if you compare it to languages from the scientific community (Python is an exception here imho). It is simple and probably good enough for what it was designed for (system programming), but it's not very good for scientific computing.

Why? I use it for scientific computing and find it very useful for exactly this purpose. But maybe you could be more precise and explain what exactly about the type system is not powerful enough.

Imho, a good type system for scientific computing will let you check that your computations are well formed. For instance, if one tries to apply a matrix to a vector of the wrong size, the type system should report it, or rather one should be able to express this sort of constraint using the type system to the same effect.

You are however correct, that doesn't mean that C cannot be useful for scientific computing, but I think that there are better alternatives. Although, one can write the core functionalities in C for speed and efficiency, and use another language as a driver.

You can definitely implement generic algorithms in C and there are different ways to do this. For example, you create an abstract data type "ring" with addition, multiplication, and multiplication with a scalar, and pass objects to this type and functions for these three operations to a generic exponential function.

i think the article is confused. i write quite a bit of SIMD code for HPC and c is my goto for low level computing because of the low level memory allocation in c. In fact that tends to be the lion share of work to implement simd operations where memory isn't the bottleneck. if only it were as simple as calling vfmadd132ps on two datatypes.

Even high level languages like c# allow you to express such a thing

https://stackoverflow.com/questions/73269625/create-memory-a...

Also structlayout and fieldoffset

Reply on news.ycombinator.com