I’m a researcher, and I have been working on different integration methods. I wanted to benchmark one of my integration methods against the acados integrators, as these are the fastest most accurate integration methods that I have found so far. I’m using the Matlab interface to create a simulation model using the ‘acados_sim’ function for an IRK integrator. I want to compare this integrator to one that I have programmed in plain matlab (consisting of just basic matrix operations like addition, subtraction, multiplication, indexing). At the moment, the IRK integrator is about 2x faster, but I was wondering if the comparison is apples to apples. My system is a moderate size (I usually go with 77 states), and my understanding is that acados may be doing something under the hood that is more efficient for such system sizes that Matlab may not be doing.

To summarize, I want to know if it is fair to compare an acados integrator with one written in plain Matlab. If it is not, I would like some guidance on how to make the comparison fair. I’m afraid of wasting a lot of time on this issue, only to discover some trivial fact that makes this whole discussion a moot point. Any guidance would be much appreciated.

That sounds interesting.
I would have expected actually a bigger speedup of acados IRK with respect to a plain Matlab IRK implementation.
However, you system size is larger than what I typically tested.
In acados, we use BLASFEO linear algebra, which is fast for two reasons.

BLASFEO uses a dedicated matrix format (panel major).

BLASFEO has optimized routines common advanced hardware architectures. You should make sure acados is compiled with the dedicated target (which should automatically happen by default).

However, as long as you mention these differences in implementation, I think you can make a reasonable benchmark comparison.
Some things that come to mind:

do you use the same discretization options in your other integrator (Butcher tableau, number of steps, number of newton iterations, tolerance for the implicit system of equations?

do you also look into sensitivity propagation? or just a pure simulation?

I do want to add that my matlab code is not an IRK, but rather an exponential integrator. I want to compare performance of my exponential integrator with the fastest existing solutions. In this case, acados IRK is the fastest that I have found for my application. In particular, I’m concerned about states and sensitivities, so I do use the .get method to generate those as well. At the end of the day I’m trying to solve an optimal control problem, and for my system the clear limiting factor is the integration method. What I would really like is to have both methods implemented in an equivalent way, such that the dominant difference is the method itself. To me, the best manner to do this is to optimize my method to match the level of optimization of the acados IRK.

Naturally I have selected the same step sizes, and generally the same accuracy (though they have different characteristics in this regard). I have seen the GNSF IRK, but my system does not benefit from this formulation.

I was hoping that there was some kind of automatic method to translate matlab functions with only very simple matrix operations into an efficient BLASFEO implementation, but I have not found such functionality.

So simulation and forward sensitivities?
Then just make sure to turn adjoint and hessian sensitivities off.

I guess it depends on how important an efficient implementation is for you.
Do you just need it for the comparison or actually want to solve your OCP as fast as possible?

Such a functionality does not exist.
Matlab does offer some way to generate “efficient” C code with some toolbox, but I dont know much about that. Did you look into that already?

I guess it depends how much code your integrator is and the questions above if a reimplementation with BLASFEO is worth it for you.

Yes my mistake, specifically simulation and forward sensitivities. Indeed, I have the adjoint and hessian sensitivities turned off. I want a comparison, but also yes, I want to solve the OCP as fast as possible (it is convenient for demonstration purposes on the somewhat modest hardware I have).

I also don’t know too much about matlab code generation. From what I have seen, there is little speed benefit for built in matrix functions, see here: c++ - MATLAB Mex function isn't faster than regular function - Stack Overflow. My code is already calling a bunch of built in matrix functions, so code generation is not likely to be significantly better, if at all.

If I can expect at least 2x faster execution time, then I would consider a reimplementation. Anything less would make no difference. At 4x, it would be a no-brainer to do the reimplementation. To this end, my integrator is not very long (~50 lines).

That makes sense!
I can refer you to the benchmarks in https://publications.syscop.de/Frison2018.pdf
where Matrix-Matrix multiplications (dgemm) are compared to OpenBLAS, MKL and Eigen on Intel Haswell architecture.

I do not know how well the Matlab linear algebra performs.

But given it is a very compact code, I think a reimplementation would be worth it.