romeric
commented
4 years ago I can see that you have implemented a LinearSolver in your code which seems like a standard LU solver with partial pivoting. May I ask what are the size of your LHS matrices?
You can by-pass all that code and use one of Fastor's solvers if you want to get a bit more performance
// inversion based solve - very fast
Tensor<T,N> x = solve(A,b);
// or equivalent of what you are doing but faster using template recursion
Tensor<T,N> x = solve<SolveCompType::BlockLUPiv>(A,b);
// or exactly equivalent of what you are doing I guess
Tensor<T,N> x = solve<SolveCompType::SimpleLUPiv>(A,b); You will gain significant speed up if your matrix sizes are less than 32x32 using inversion based or block LU decomposition. For matrices of bigger sizes you will still get the speed-up but given that Fastor implements most things using compile time recursion rather than for loops, the compilation time might blow-up quite quickly, specially with GCC that has nearly quadratic template instantiation cost. Clang should be fine.
These solvers are not heavily tested with Visual Studio but compiling with clang Fastor solvers are faster than Eigen, MKL and straight-forward for loop style implementations (such as yours) by quite a big margin. Here are the results from a comparison I did a while ago on my macbook with Julia 1.4 backslash solve and Fastor BlockLU solve compiled with clang using -mavx2 -mfma -O3 if it is of any interest
| Size | Julia | Julia StaticArrays | Fastor |
|---|---|---|---|
| 8x8 | 942ns | 350ns | 95ns |
| 16x16 | 2.831 μs | 3.036 μs | 807ns |
| 32x32 | 25.851 μs | 25.517 μs | 4.11538 µs |
| 64x64 | 95.105 μs | 104.029 μs | 25.8607 µs |
martinholters
ghost
Hi Martin,
First of all, thank you very much for all your scientific works, and especially for ACME that is a great project.
Since a moment, I attempt to port ACME to C++ using Armadillo for Linear Algebra. Unfortunatly, the Sparse Matrix of Armadillo do not cover the potential of Julia's Sparse Arrays standard library, in this case the maintain of structural zero. So i've write a C++ version of the Julia Sparse Arrays library.
After that, I've follow the ACME code design using latest standards of C++, heavy usage of meta programming to keep in place with the flexibility of the Julia code. The discretization of the circuit is, for my own configuration more faster in C++ than the Julia code. I mean, computation of incidence and topology matrices, the non-linear decomposition and so on. But of course, no possibility to reproduce the symbolic trick for the build of the circuit.
For example, the build of the SuperOver circuit become:
For a build and design simplification, the DiscreteModel class follow the PIMPL idiom and hide the Model Runner :
Now, I can run the circuit using a 30 seconds unprocessed guitar excerpt. No problem, I get the same result than Julia ACME process ... but in 57 seconds where the excerpt take 28 seconds with the Julia code !!! What ?!?!
My fault, by following the ACME design using lambda function (for genericity) make it unusable in such of process. Optimization with Intel MKL and cblas call do not have effect. Poor performance.
So after reading this : https://medium.com/@romanpoya/a-look-at-the-performance-of-expression-templates-in-c-eigen-vs-blaze-vs-fastor-vs-armadillo-vs-2474ed38d982 , I take a look to the famous Fastor library.
After some code change, I'm able to produce on-the-fly a static header only file that is a freezed version of the DiscreteModel of a circuit using Fastor at backend for the matrices and linear algebra. No dependencies except Fastor that is header only too. There's some changes on the design to allow some generic code, for example LinearSolver, SimpleSolver, HomotopySolver become generic (template based), ParametricNonLinEq become abstract and template base to allow multiple solvers instance (same as ACME core solvers module).
So, in a subnamespace details, you will find reusable code :
The others classes are produced programmaticaly including compile time static versions of the precomputed matrices as Fastor's Tensor. Because the possibility of multiple solvers for the non-linear parts, the code produce overrided classes of ParametricNonLinEq base class and for the SuperOver with fixed potentiometer values case, only one :
Ofcourse, this is not yet fully optimized, but this class can be used by the generic Solver classes and be runned with a class that is a mixed version of the DiscreteModel and ModelRunner of ACME :
As you can see, for the moment the SuperOver class can be instanciated with a compile-time chunck size for the I/O, fixed by default at 1024 columns. A static method return the sampling rate used for the discretization of the circuit.
SuperOver-Fastor The full code (single file header-only, no dependencies except Fastor) is available here for the proof of concept : https://codeshare.io/5QEE4q
Now this is the report for the performance with my configuration, all the tests are running 10 times on a 30 seconds audio sample file at 44100 Hz with the SuperOver circuit discretized with the 3 fixed potentiometers at full value (1.0) :
Julia 32bits and 64 bits
min : 23.407315 seconds max : 29.414382 seconds mean : 25.026536 seconds
21.512428 seconds (430.89 k allocations: 21.301 MiB) 21.094247 seconds (430.88 k allocations: 21.295 MiB) 20.971705 seconds (430.88 k allocations: 21.295 MiB) 21.836396 seconds (430.88 k allocations: 21.295 MiB) 21.668114 seconds (430.88 k allocations: 21.295 MiB) 20.907806 seconds (430.88 k allocations: 21.295 MiB) 20.982056 seconds (430.88 k allocations: 21.295 MiB) 21.016102 seconds (430.88 k allocations: 21.295 MiB) 20.938037 seconds (430.88 k allocations: 21.295 MiB) 21.450474 seconds (430.88 k allocations: 21.295 MiB)
min : 20.907806 seconds max : 21.836396 seconds mean : 21.237736 seconds
17.2277 seconds 17.4744 seconds 17.4055 seconds 17.3357 seconds 17.2751 seconds 17.3368 seconds 17.4383 seconds 17.3296 seconds 17.3059 seconds 17.2541 seconds
min : 17.2277 seconds max : 17.4744 seconds mean : 17.3383 seconds
13.9769 seconds 13.9896 seconds 14.1125 seconds 13.9835 seconds 13.9530 seconds 14.7414 seconds 14.2672 seconds 14.1663 seconds 14.7137 seconds 14.3044 seconds
min : 13.9769 seconds max : 14.7414 seconds mean : 14.2208 seconds
9.24264 seconds 9.26479 seconds 9.35169 seconds 9.35638 seconds 9.51796 seconds 9.44793 seconds 9.34992 seconds 9.52753 seconds 9.87442 seconds 9.74665 seconds
min : 9.24264 seconds max : 9.87442 seconds mean : 9.46799 seconds
8.62546 seconds 8.29256 seconds 8.29380 seconds 8.30038 seconds 8.30197 seconds 8.32373 seconds 8.51325 seconds 8.30997 seconds 8.41125 seconds 8.31047 seconds
min : 8.29380 seconds max : 8.62546 seconds mean : 8.36828 seconds