Inline vector operations in C

Question

I am using C to do some scientific computing, and need to move around a lot of vectors.

I've defined some basic Multiply and Add functions of the form

add(int size, double *out, double *x, double *y)

But for complex operations the code quickly becomes long and hard to read.

Is it possible to define inline operators (V1 + V2) ? Or any general best practices that will make it easier to spot check code for math errors? Perhaps some #define shenanigans?

Answer 1

At first I thought this was a C++ question for some reason!

If you can use C++ you might be able to use an STL array (where the size is a template parameter, not a stored value).

It might look something like this:

std::array<double, 7> a = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0};
std::array<double, 7> b = {0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0};
auto c = a + b;
auto d = a + b + c;

And you can define operator overload for operator+. E.g:

#include <array>
// ...
// (version 1) returns a new copy/temporary
template <class T, size_t size>
std::array<T, size> operator+(const std::array<T, size>& a, 
                              const std::array<T, size>& b)
{
    std::array<T, size> c;
    // add them however you want; simple version:
    for (size_t i = 0; i < size; ++i)
    {
        c[i] = a[i] + b[i];
    }
    return c;
}

// (version 2) no temporaries via rvalue ref/move semantic
template <class T, size_t size>
std::array<T, size>&& operator+(std::array<T, size>&& a, 
                                const std::array<T, size>& b)
{
    for (size_t i = 0; i < size; ++i)
    {
        a[i] += b[i];
    }
    return std::move(a);
}

so that:

auto c = a + b;       // invokes version 1
auto d = a + b + c;   // invokes version 1 for (b+c) 
                      // and then version 2 for a+(temp)

so you use at most one temporary in any chain of operations.

The memory layout of a std::array should be the same as a native array, so you should be able to "inject" this code quite easily into your existing program (by (ab)using casts, typedefs or preprocessor) if you want to touch a minimal amount of existing code.

Answer 2

If you don't mind a non-standard solution to the problem, the Intel compiler supports extensions for array notation in C. Here's a link:

http://software.intel.com/sites/products/documentation/studio/composer/en-us/2011/compiler_c/optaps/common/optaps_par_cean_prog.htm

Basically, you can write code like this:

int A[4] = { 1, 2, 3, 4 };
int B[4] = { 4, 3, 2, 1 };
int C[4];
C[:] = A[:] + B[:];

It also provides vector versions of many math functions, so you can do stuff like this:

B[:] = cosf(C[:]);

As an added bonus, if you write your code like this, it gets vectorized in SSE/AVX instructions.

Answer 3

This is called “operator overloading” and it is not a feature of C. It is a feature of C++ though. If C++ is an option for you, there are plenty of tutorials on how to overload operators. Some consider it an evil because:

1. When code is taken out of context it is difficult to tell what will actually happen. For example, you can see operators there but how do you know whether they are the built-in ones or overloaded ones?

2. It is possible to overload operators and apply meaning to them that is not common, for example (and this is something I have seen in actual code), you can overload the “/” operator (division operator) to act as a path separator for a class that handles filenames in an OS-independent way:

   Filename f = Filename::MainFileSystem() / "folder1" / "folder2" / "file.txt";
   

   Don't do this.

There are also benefits to using overloaded operators:

1. Code can be heavily simplified. Rather than multiple function calls, you can use

   result = vec1 + vec2 + vec3;
   

   Note that as above, it is difficult to tell just from looking at that line of code whether any functions will be called but whether that difficulty outweighs the benefit of simplification is subjective.

2. You can perform some nifty tricks, such as smart pointers, that give you added functionality without modifying the original semantics of the operators.

Answer 4

Is it possible to define inline operators (V1 + V2) ?

No. But it is in C++.