## 2 Matrix Sign Function-Based Solvers
The use of several GPUs also increments the computational power of the platform and, thus, can potentially reduce the execution time. The implementations for this type of platform extend the hybrid CPU-GPU codes, while dealing with several difficulties, particularly, load balancing and data transfers. In general, the data layout exerts a relevant impact on load-balance and data transfers. However, the special properties of the GJE algorithm allow us to rely on a simple block cyclic data distribution, as shown in Figure 1. The matrix is partitioned there into blocks of columns, which are then distributed cyclically among the GPUs. Additionally, in the multi-GPU codes, some operations that are performed are merged, reducing the number of invocations to BLAS kernels and avoiding the computation of small to moderate matrix-matrix products that offer a poor performance on the GPUs. In particular, the update of \([A_{00}A_{02};A_{10}A_{12};A_{20}A_{22}]^{T}\) can be performed with a single matrix-matrix product. These and other features of this version are detailed in [19]. ### Parallel Matrix Equations Solvers

The optimization of the matrix equation solvers is not restricted to the matrix inversion routines and, by extension, the calculation of the matrix sign function. All the phases of Algorithm 2--GECLNC, Algorithm 3--GECRSG and Algorithm 4--ROS1 were carefully optimized for the different platforms. Figure 1: Data distribution between graphics processing units (GPUs) for the multi-GPU implementation (example with four GPUs). In the case of the Lyapunov solver, high performance linear algebra kernels from BLAS and LAPACK were employed to execute most of the additional computations that appear in Algorithm 2--GECLNC (QR factorization, matrix-matrix products, computation of the scaling factor, convergence test, _etc._). OpenMP directives were also employed to parallelize other minor computations, like matrix addition and matrix scaling. The ARE solvers build matrix \(H\) in parallel and employ OpenMP directives to efficiently compute the matrix addition and scaling required for the update of \(H_{k+1}\). Besides, multi-threaded routines from BLAS and LAPACK were used to solve the over-determined system in the last stage of the method. In the DRE solvers, the computation of the low rank factor for the solution, \(Z_{k+1}\), dramatically reduces the computational cost of all steps. Thus, the total cost of the algorithm basically boils down to that required by the solver of the Lyapunov equations. Once more, multi-threaded BLAS and LAPACK kernels and OpenMP directives were employed to parallelize some minor computations, like matrix additions and matrix norm computations. The updates of \(\tilde{A}_{k}\) and \(\tilde{N}_{k}\) require several matrix-matrix products, matrix additions and scalings that are computed on the CPU. Although matrix-matrix products of large matrices are suitable for the GPU architecture, the dimensions of \(B\) and \(Z^{T}\) are usually too small to amortize the cost of data transfers. ## 4 Experimental Results

We evaluate the performance of the implementations using two problems from the _Oberwolfach Model Reduction Benchmark Collection_[20]. For the STEEL problem, we employ two instances, steel\({}_{S}\) and steel\({}_{L}\). For both cases, \(m=7\) and \(p=6\). The order of the system is \(n=1,357\) for steel\({}_{S}\) and \(n=5,177\) for steel\({}_{L}\); for the FLOW_METER the dimensions are \(n=9,669\), \(m=1\), \(p=5\). Experiments are performed on two different platforms. The first platform, Kepler, is equipped with an Intel Sandy Bridge-E i7-3930K processor at 3.2 GHz, 24 GB of RAM, and is connected to an NVIDIA Tesla K20 via a PCI-e bus. The second platform, S2050, consists of two Intel Xeon QuadCore E5440 processors at 2.83 GHz, with 16 GB of RAM, connected to an NVIDIA Tesla S2050 (four NVIDIA M2050 GPUs) via two PCI-e buses. Platform Kepler is employed to evaluate the multicore and hybrid CPU-GPU kernels, while the multi-GPU kernels are evaluated in S2050.