Optimizing C++ equivalent of Matlab `filter` function

Minimum working example below. The Matlab filter function filters the input data x using a rational transfer function defined by the numerator and denominator coefficients b and a and initial conditions z. My C++ version of the code is way slower than the Matlab one. My signal x is a huge valarray, a and b are three elements long each and the initial conditions z is just three zeros.

#include 
#include
#include 
using std::rand;
#include 
using std::valarray;
#include 

#define D_SCL_SECURE_NO_WARNINGS 1

struct FilterResults
{
    valarray filteredSignal;
    valarray finalConditions;
};

FilterResults Filter(const valarray &b,
    const valarray &a, const valarray &x, const valarray &z)
{
    valarray y(x.size());
    // Pad z
    valarray zOut(a.size());
    std::copy(std::begin(z), std::end(z), std::begin(zOut));

    double Xm, Ym;
    size_t i, n = a.size();
    for (size_t m = 0; m  s = zOut[std::slice(0, zOut.size() - 1, 1)];
    FilterResults r;
    r.filteredSignal = std::move(y);
    r.finalConditions = std::move(s);
    return r;
}

int main() {
    std::srand(std::time(0));
    valarray b{ { (double)rand(), (double)rand(), (double)rand() } };
    valarray a{ { (double)rand(), (double)rand(), (double)rand() } };
    valarray z{ { 0, 0, 0 } };
    valarray x( 500000 );
    for (size_t i = 0; i < 500000; i++)
    {
        x[i] = (double)rand();
    }

    clock_t begin = clock();
    auto r = Filter(b, a, x, z);
    clock_t end = clock();
    double elapsed_secs = double(end - begin) / CLOCKS_PER_SEC;
    std::cout << "elapsed_secs: " << elapsed_secs << std::endl;
    system("pause"); // C++ is 0.184 seconds on my computer, Matlab is 0.006060.
}

Any idea on how to speed it up?

Few highlights:

std::copy(std::begin(z), std::end(z), std::begin(zOut));

To copy input values to overwrite them is just a waste of time. You will overwrite them then you do not need to copy old values, use z instead of zOut inside loops.

y.size() is possibly evaluated multiple times. It's a minor performance issue (and it may be inlined) but it depends on your actual implementation.

I don't know about valarray implementation but if this is performance critical you may want to use a plain array.

I hope FilterResults constructor doesn't allocate anything (what valarray default constructor does?) That said, this is C++, filter needs a state then you should design this using a Filter class with its private instance fields. You return just filtered data and you do not create/copy/move status. I imagine a Filter class like this:

class Filter
{
public:
    Filter(const valarray& b, const valarray& b)
    {
        // Validate parameters...
        this->a = a;
        this->b = b;

        // I omit z, I think there are very few chances you want
        // to set initial filter state...
    }

    valarray& Process(const valarray& input)
    {
        // Your code here, when you're done just return computed
        // result. "Final conditions" (filter status?) has not to
        // be returned, it's a private field.
    }

private:
    valarray b, a, z;
};

As you can note some calculations/allocations will be performed in constructor. It won't change anything if each time you filter a different chunk of samples, using different filter, but in the other cases it will allow a big speed-up. Next step in design (to do something in C++) is to abstract away filter and its design: you may have a pure abstract base class Filter, an IirFilter implementation and one or more classes derived from FilterDesigner (let's say ChebyshevIirFilterDesigner to calculate filter coefficients). In general to mimic MATLAB approach may result in a sub-optimal implementation (because languages are different and because it will lead to a less than perfect OOP design.) If this class is a use-once function and you want to keep a function instead of a class you may drop finalConditions, do you actually use it? If you're implementing an on-line filter (for a continuous stream) then you'd better save it as private field (it's the z parameter you're using now)

You should never have using namespace std; You already use std:: consistently then you just add a possible source of conflicts.

Few more notes about performance:

In this case (I guess they do not use any esoteric math theory to perform an optimization) C++ performance should be pretty similar or even better than MATLAB. MATLAB may (I do not know what current implementation does) use SIMD/AVX instructions. Do not think you can have same throughput if you don't.

To perform a benchmark is more complex than this! At very least you should setup your machine for that and to perform calculations few thousands times to have an average and a distribution then inspect results to understand caching issues...

StackExchange Code Review Q#141798, answer score: 6