remarkRemark hypothesisHypothesis claimClaim Coded QR DecompositionQuang Minh
Nguyen, Iain Weissburg and Haewon Jeong ex_supplement

# Coded Computing for Fault-Tolerant Parallel QR Decomposition ††thanks: A
preliminary version of this paper was presented at IEEE International
Symposium on Information Theory 2020.

Quang Minh Nguyen  Massachusetts Institute of Technology<EMAIL_ADDRESS>Iain Weissburg  University of California, Santa Barbara<EMAIL_ADDRESS>Haewon Jeong  University of California, Santa Barbara<EMAIL_ADDRESS>
###### Abstract

QR decomposition is an essential operation for solving linear equations and
obtaining least-squares solutions. In high-performance computing systems,
large-scale parallel QR decomposition often faces node faults. We address this
issue by proposing a fault-tolerant algorithm that incorporates ‘coded
computing’ into the parallel Gram-Schmidt method, commonly used for QR
decomposition. Coded computing introduces error-correcting codes into
computational processes to enhance resilience against intermediate failures.
While traditional coding strategies cannot preserve the orthogonality of $Q$,
recent work has proven a _post-orthogonalization condition_ that allows low-
cost restoration of the degraded orthogonality. In this paper, we construct a
checksum-generator matrix for multiple-node failures that satisfies the _post-
orthogonalization condition_ and prove that our code satisfies the maximum-
distance separable (MDS) property with high probability. Furthermore, we
consider _in-node checksum storage_ setting where checksums are stored in
original nodes. We obtain the minimal number of checksums required to be
resilient to any $f$ failures under the in-node checksum storage, and also
propose an in-node systematic MDS coding strategy that achieves the lower
bound. Extensive experiments validate our theories and showcase the negligible
overhead of our coded computing framework for fault-tolerant QR decomposition.

###### keywords:

High-performance computing, Matrix factorization, Coding theory, Algorithm-
based fault-tolerance

15A23, 94B05, 65Y05

##  1 Introduction

Motivated by the widespread use of machine learning and scientific computing
applications that require large-scale distributed computations, _coded
computing_ has been a burgeoning area of research [27, 14, 44, 53, 15]. The
aim of coded computing is to add systematically-designed redundancies in the
computing algorithm to build fault resilience that is more efficient than
replication or checkpointing. In this work, we consider protecting the QR
decomposition algorithm from failures in high-performance computing (HPC)
systems. Building a reliable supercomputer has been a long-standing problem,
but the emergence of exascale computing poses new challenges that require a
new and innovative solution that goes beyond traditional reliability
techniques. More than 20% of the computing capacity in today’s HPC systems is
wasted due to failures and ensuing recovery [18], and this wastage is only
expected to grow as the system size grows. To reduce the overhead of fault
tolerance in upcoming HPC systems, _algorithm-based fault-tolerance (ABFT)_
for HPC has been suggested [6], the core idea of which is essentially the same
as coded computing: adding encoded redundancy tailored to a given numerical
algorithm.

QR decomposition factors a matrix into a product of an orthogonal matrix ($Q$)
and an upper triangular matrix ($R$). It is an essential building block of
linear algebraic computations as it provides a numerically stable method for
solving linear equations and linear least squares problem (e.g. linear
regression). As it is an important computation primitive, ABFT for parallel QR
decomposition has been studied in the HPC literature [12, 50, 31]. While these
existing works have shown that ABFT can be a more efficient resilience
solution than periodic checkpointing, they proved that coding-based fault-
tolerance can only be a half of the solution: while the $R$-factor can be
protected with error-correcting codes (ECCs), the same does not apply to the
$Q$-factor as encoding destroys the orthognal structure of the $Q$-factor [12,
Theorem 5.1]. Our prior work [32] introduced a novel strategy that can
circumvent this fundamental theorem, by using a specialized code design that
can restore the $Q$-factor’s orthogonality with minimal overhead.However, this
approach was limited to single-node failure scenarios, lacking a clear
extension to multiple-node fault tolerance.

In this paper, we make critical contributions in generalizing the coded QR
decomposition with low-cost post-orthogonalization to a wider variety of
settings, as summarized below:

  * •

We propose a novel coding strategy for parallel QR decomposition that
tolerates multiple-node failures with minimal overhead. Standard coding
designs adequately protect the $R$-factor. However, protecting the $Q$-factor
poses a challenge as its orthogonality is compromised by conventional linear
coding [12]. We address this by devising a post-processing framework that
restores $Q$’s orthogonality. Specifically, we derive a _post-
orthogonalization condition_ which states that if a checksum-generator matrix
follows certain structure, there exists a sparse matrix $G_{0}$ that can
linearly transform the coded $Q$ matrix into an orthogonal matrix. To this
end, we design a semi-random checksum generator matrix satisfying the post-
orthogonalization condition, while achieving the maximum-distance separable
(MDS) property, i.e. code optimality, with high probability. The proposed
algorithm can be applied to solving a full-rank square system of linear
equations in HPC.

  * •

We further consider the emerging setting of _in-node checksum storage_ where
the coded data is stored in original systematic processors instead of being
distributed over extra checksum processors for fault-tolerance. We prove a
lower bound on the number of checksums required to be resilient to any $f$
failures and also demonstrate coding schemes for both $Q$-factor and
$R$-factor protection that meet the bound. The proposed in-node checksum
coding strategy is optimal, thereby improving over the results in [12, 32].

Finally, we note that this paper is an extension of our conference paper [32].
For both contributions, $Q$-factor protection codes and in-node checksum
storage codes, we had a result for the single node failure case in [32]. In
this paper, we generalize those results to the multiple-node failure scenario.
The results and constructions in [32] can be regarded as a special case within
the broader framework developed in this paper, tailored specifically for
single-node failures.

##  2 System Model and Problem Formulation

We first describe our notations throughout the paper, and then present the
system model and problem formulation within interests.

Notations: We use upper-case letters to denote matrices and bold lower-case
letters to denote vectors. For matrix $A$, we use $a_{i,j}$ to denote the
$(i,j)$-th entry of the matrix. We use $I_{m}$ to denote the identity matrix
of dimension $m\times m$, and write $I$ when the dimension is clear from the
context.

Computation of interest: We consider solving a system of linear equations
$A\bm{x}=\bm{b}$ through QR decomposition, where $A$ is a square and full-rank
matrix of dimension $n\times n$. Solving square and non-singular system of
linear equations is a fundamental building block for many applications in HPC
[11, 2, 42], and QR decomposition is a popular choice in practice due to its
guaranteed stability and computational efficiency [22].

Distributed computing model: We assume that $A$ is a large matrix and the
computation has to be distributed to multiple processors (nodes)111We use
_nodes_ and _processors_ interchangeably in this paper.. We assume that we
have a total of $P$ nodes aligned as a $p_{r}\times p_{c}$ grid
($p_{r}p_{c}=P$). We denote the $(i,j)$-th node on the grid as $\Pi(i,j)$. We
consider _2D block distribution_ of data where the matrix $A$ is split into
$p_{r}\times p_{c}$ blocks and each block $A_{ij}$
($i\in[0,p_{r}-1],j\in[0,p_{c}-1]$) is distributed to processor $\Pi(i,j)$.
The $j^{th}$ block-column is the ”column” of data blocks
$A_{0j},A_{1j},...,A_{(p_{r}-1)j}$. Block-rows are defined similarly. For
simplicity, we assume that $p_{r}=p_{c}$ and consequently set the dimension of
each block $A_{i,j}$ to be $b=\frac{n}{p_{r}}=\frac{n}{p_{c}}$. Note that
$n\gg P$, and thus $n$ is much greater than $p_{r}$ or $p_{c}$.

Failure model: We assume that nodes are subject to failures, especially to
_“fail-stop errors.”_ In this model, a failure corresponds to a node that
completely stops responding, and loses all of its data, including its inputs
and computation results. This is a realistic failure model in HPC [4] and also
similar to the _“straggler model”_ commonly used in the coded computing
literature. At any iteration of the QR decomposition, we assume that at most
$f$ failures can occur per column or per row, and that the identities of the
processors that fail are provided by some external source (e.g. Message
Passing Interface (MPI) library [23]). Further, we assume that $f$ is small
compared to $p_{r},p_{c}$. Specifically, in our proposed scheme, we require
$f\leq\frac{1}{2}\min\\{p_{r},p_{c}\\}$.

Checksums: In this work, we use checksums to recover from failures. Before we
perform QR factorization, we encode the $n\times n$ matrix $A$ is as:

(1) $\displaystyle\widetilde{A}=\begin{bmatrix}A&AG_{h}\\\
G_{v}A&G_{v}AG_{h}\end{bmatrix},$

where $\widetilde{A}$ is of dimension $(n+c)\times(n+d)$ and $G_{h}$ and
$G_{v}$ are checksum-generator matrices of size $n\times d$ and $c\times n$.
We call $G_{v}A$ _“vertical checksums”_ or _“checksum rows”_ and $AG_{h}$
_“horizontal checksums”_ or _“checksum columns”_. After computation, vertical
checksums protect $Q$ and horizontal checksums protect $R$ [32, Lemma 1 and
Corollary 1.1]. We assume that $c$ and $d$ are small to keep the overhead
negligible.

Conventionally, to store the checksum part, we introduce additional nodes
called _“checksum nodes.”_ We call this setting _out-of-node checksum
storage_. In this setting, the original matrix $A$ is stored in _“systematic
nodes”_ and all the checksum entries are stored in the _checksum nodes._
Assuming that we use the same encoding coefficients at each node for the
entire data block, $G_{v}$ and $G_{h}$ can be represented as:

(2) $\displaystyle G_{v}=\widetilde{G}_{v}\otimes
I_{n/p_{r}},\;\;G_{h}=\widetilde{G}_{h}\otimes I_{n/p_{c}},$

where $\otimes$ is a Kronecker product, $\widetilde{G}_{v}$ and
$\widetilde{G}_{h}$ are respectively $m_{r}\times p_{r}$ and $p_{c}\times
m_{c}$ matrices, and $m_{r}$ and $m_{c}$ are respectively the number of
checksum nodes per row and column.

In this paper, we will also consider an _in-node checksum storage_ setting,
where one stores the checksums within the original set of nodes instead of
introducing a separate set of checksum nodes. This will increase the size of
the computing task per node, but it is better suited for a scenario where a
user has a fixed number of available nodes.

###  2.1 Problem Statement

Our goal is to develop a coding scheme that can add effective checksums for
parallel QR decomposition, which can protect both Q-factor and R-factor from
any $f$ fail-stop errors ($f\geq 1$). This builds up on our previous work [32]
that provided a coding scheme for a single fail-stop error ($f=1$) and
generalizes to multiple errors and applies to both out-of-node and in-node
checksum storage settings (see Table 1). We will provide theoretical and
experimental evidences to demonstrate the small overhead of the proposed
coding scheme.

| Coding
---
Scheme
QR Algorithm |  $Q$-factor |  $R$-factor | $f>1$ | Out-of-node | In-node
[12] | Householder | ✗ | ✓ | ✗ | ✗ | ✓
[26] | Householder | ✗ | ✓ | ✗ | ✓ | ✗
[31] | Givens Rotation | ✗ | ✓ | ✗ | ✓ | ✗
[32] | Gram-Schmidt | ✓ | ✓ | ✗ | ✓ | ✓
This work | Gram-Schmidt | ✓ | ✓ | ✓ | ✓ | ✓
Table 1: Summary of work that utilize checksum-protection for QR
decomposition. Those work [12, 26, 31] without coding scheme for $Q$-factor
protection rely on costly checkpointing for fault-tolerance.

##  3 Preliminaries and Related Work

The core part of this work is incorporating error-correcting codes into a
parallel QR decomposition algorithm for efficient checksum generation and
failure recovery. For the QR factorization algorithm, we use parallel block
modified Gram-Schmidt (PBMGS) algorithm. For error-correcting codes, we adapt
the well-known maximum distance separable (MDS) codes. In this section, we
describe essential preliminaries for these.

Input:
$\widetilde{A}=\begin{bmatrix}\bm{\widetilde{a}_{1}}&\bm{\widetilde{a}_{2}}&\cdots&\bm{\widetilde{a}_{n+d}}\end{bmatrix}$
is $(n+c)\times(n+d)$ matrix

Output:
$Q=\begin{bmatrix}\bm{q_{1}}&\bm{q_{2}}&\dots&\bm{q_{n+d}}\end{bmatrix}$ is
$(n+c)\times(n+d)$ matrix , and $R=(r_{ij})$ is $(n+d)\times(n+d)$ matrix

Result: $\widetilde{A}=QR$, where $Q^{T}Q=I$

1 for  _$i=1$ to $n$ _ do

2 $\bm{u_{i}}=\bm{\widetilde{a}_{i}}$

3 end for

4for  _$i=1$ to $n+d$ _ do

5 $r_{ii}=||\bm{u_{i}}||_{2}$

6 $\bm{q_{i}}=\bm{u_{i}}/r_{ii}$

7 for  _$j=i+1$ to $n+d$ _ do

8 $r_{ij}=\bm{q_{i}}^{T}\bm{u_{j}}$

9 $\bm{u_{j}}=\bm{u_{j}}-r_{ij}\bm{q_{i}}$

10 end for

11

12 end for

Algorithm 1 Modified Gram-Schmidt

###  3.1 Parallel block modified Gram-Schmidt (PBMGS) algorithm

There are three classes of commonly used algorithms for QR decomposition:
Gram-Schmidt (GS) [43, 39, 34], Householder Transformation [38, 3], and Givens
Rotation [17]. In this work, we consider parallel block modified Gram-Schmidt
(PBMGS), which performs QR factorization on a grid of distributed nodes. For
extreme-scale QR factorization, PBMGS is a preferred choice as it has a low
computational cost and is easy to implement. In the Gram-Schmidt algorithm, we
iteratively obtain a set of orthonormal vectors by subtracting the projections
of each column onto the space spanned by orthogonal vectors computed in the
previous iterations. Modified Gram-Schmidt (MGS) performs this operation in a
slightly different order to achieve better numerical stability (see Algorithm
1). Block MGS performs the same operation in column-blocks instead of column
by column. How this can be parallelized in the PBMGS algorithm is depicted in
Fig. 1 (see Appendix A for the full algorithm and its overhead analysis).

Figure 1: The parallel block modifed Gram-Schmidt (PBMGS) algorithm on a 3x3
grid.

###  3.2 Maximum distance separable (MDS) codes

What is the minimal number of checksum nodes we need to recover the data when
$f$ nodes fail? This is a classic question in the coding theory literature.
The information-theoretic optimum is $f$, given by the Singleton bound [25].
MDS codes are the class of error-correcting codes that meet this bound [29].
MDS codes are linear block codes, meaning that each checksum is a linear
combination of original entries and checksum generation can be represented
with matrices as we did in (1). For codes to have the MDS property, checksum
generator matrices should satisfy certain rank conditions. For our problem
setup, these conditions can be translated as below:

######  Theorem 3.1 (Modified Theorem 8, p. 321 in [29]).

$\widetilde{G}_{v}$ (resp. $\widetilde{G}_{h}$) has the MDS property if and
only if every square submatrix of $\widetilde{G}_{v}$ (resp.
$\widetilde{G}_{h}$) is full-rank.

###  3.3 Checksum preservation

We now review the desirable checksum preservation property of the PBMGS
algorithm. When we perform QR decomposition on the encoded matrix
$\widetilde{A}$, the PBMGS algorithm updates the $Q$ and $R$ factors in $T$
iterations. At the end of each iteration the algorithm maintains the
intermediate $Q$-factor $Q^{(t)}$ of size $(n+c)\times(n+d)$ and $R$-factor
$R^{(t)}$ of size $(n+d)\times(n+d)$. The initial values of the $Q$-factor and
$R$-factor are set to be $Q^{(0)}=\widetilde{A}$ and
$R^{(0)}=0_{(n+d)\times(n+d)}$. After the $T$-th iteration, we retrieve the
orthogonal $Q=Q^{(T)}$ and the upper triangular $R=R^{(T)}$ as the output for
the factorization $\widetilde{A}=QR$. We further consider the following
submatrices of $Q^{(t)}$ and $R^{(t)}$:

$\displaystyle Q^{(t)}=\begin{bmatrix}Q_{1}^{(t)}\\\
Q_{2}^{(t)}\end{bmatrix},R^{(t)}=\begin{bmatrix}R_{1}^{(t)}&R_{2}^{(t)}\end{bmatrix}$

where $Q_{1}^{(t)}$: $n\times(n+d)$, $Q_{2}^{(t)}$: $c\times(n+d)$,
$R_{1}^{(t)}$: $(n+d)\times n$, and $R_{2}^{(t)}$: $(n+d)\times d$. Lemma 3.2
shows that the original vertical checksums $G_{v}A$ and horizontal checksums
$AG_{h}$ are preserved after each iteration as $Q_{2}^{(t)}=G_{v}Q_{1}^{(t)}$
for $Q$-protection and $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$ for $R$-protection.

######  Lemma 3.2.

For both sequential and parallel MGS, the following checksum relations hold
for $t\in[0,T]$:

(3) $\displaystyle Q_{2}^{(t)}=G_{v}Q_{1}^{(t)},$ (4) $\displaystyle
R_{2}^{(t)}=R_{1}^{(t)}G_{h}.$

######  Corollary 3.3.

At the end of the QR decomposition, the factorization of $\widetilde{A}$ has
the final form:

(5) $\displaystyle\widetilde{A}=\begin{bmatrix}Q_{1}\\\
Q_{2}\end{bmatrix}\begin{bmatrix}R_{1}&R_{2}\end{bmatrix}=\begin{bmatrix}Q_{1}\\\
G_{v}Q_{1}\end{bmatrix}\begin{bmatrix}R_{1}&R_{1}G_{h}\end{bmatrix}$

where $Q_{i}=Q_{i}^{(T)}$ and $R_{i}=R_{i}^{(T)}$.

The proofs of Lemma 3.2 and Corollary 3.3 can be found in Appendix B.

######  Remark 1.

All the prevalent QR decomposition methods—Householder [12, 50], Givens
Rotations [31], and MGS—preserve the $R$-factor checksums. However, only MGS
permits $Q$-factor checksum-preservation (i.e., satisfies (3)) since it uses
linear projection of the input matrix’s column vectors for orthogonalization.
On the other hand, Householder and Givens Rotations rely on non-linear
operations such as reflections and rotations, which do not preserve the
linearity of $Q$-factor checksums.

###  3.4 Related work

Algorithm-based fault tolerance (ABFT) has been proposed to alleviate the high
overhead of checkpoint-restart (C/R) schemes in exascale high-performance
applications [9, 5, 51, 12, 13, 50, 28, 35, 21, 10, 31]. In the information
theory community, a similar idea called _coded computing_ has been extensively
studied, with more focus on renovating the code design under a theoretical
computing model [53, 16, 14, 45, 36, 37, 24, 52, 1, 20, 19, 30, 49, 33, 41,
7]. There exist several works on ABFT for QR decomposition [12, 50, 31, 26],
and they are summarized in Table 1. For the out-of-node checksum storage
setting, in [26, 50], the authors develop ABFT schemes for the Householder
algorithm and in [31], they study ABFT for the Givens Rotation algorithm. The
ABFT scheme for Householder is further adapted to the in-node checksum storage
setting in [12]. However, all of them only protect the $R$ factor with
checksums and rely on replication or checkpointing for the $Q$-factor. In our
prior work [32], we proposed a coded MGS algorithm that can provide full
protection for both $R$ and $Q$ factors when _there is only one failure during
computation ( $f=1$)_. This work extends [32] and provide a generalized scheme
that can have resilience for multiple failures.

Figure 2: Coded QR decomposition framework.

###  3.5 Overview of the Proposed Coded QR Decomposition Scheme

We present the overall framework of coded QR decomposition in Figure 2. It
comprises of encoding, (coded) QR factorization, failure recovery (decoding),
and post-orthogonalization steps. In the encoding phase, the input matrix $A$
is encoded into $\widetilde{A}$ as given in (1). We then perform QR
decomposition of the encoded matrix $\widetilde{A}$. If a node fails during
the computation, we use vertical checksums to recover the intermediate
$Q$-factor from the failed node and horizontal checksums to recover the
$R$-factor. The original input matrix $A$ can be restored by either type of
checksums. At the end of the computation, we retrieve $A=Q_{1}R$. As we will
show in Section 5, $Q_{1}$ is not orthogonal, and hence we deploy the post-
orthogonalization technique which transforms $Q_{1}$ into $G_{0}Q_{1}$ and
$\mathbf{b}$ into $G_{0}\mathbf{b}$, which are a very low-cost operations.
After post-orthogonalization, we can solve the original problem:
$A\mathbf{x}=\mathbf{b}$.

##  4 Checksum Generation for $R$-factor protection

In this section, we show how coded computing for $R$-factor protection can be
extended from the out-of-node checksum storage setting to the in-node checksum
setting. For the out-of-node setting, the coding process is straightforward.
To protect from node failures for the $R$ factor, we can do row-wise encoding
and add horizontal checksums $RG_{h}$ to checksum nodes as depicted in Fig. 2.
Following Lemma 3.2, any linear coding schemes can preserve the checksums of
$R$ over the iterations in the PBMGS algorithm:
$R_{2}^{(t)}=R_{1}^{(t)}G_{h}$. Furthermore, $R_{1}^{(t)}$ still maintains the
upper-triangular shape. Thus, there is no need for further processing if there
is no failure. For optimal redundancy, we can use the checksum generator
matrix of any MDS codes (e.g., Reed Solomon codes). We also refer to [12, 26]
for more details of out-of-node $R$ protection.

The challenge we address in this paper is to extend the out-of-node checksum
storage to the in-node checksum storage setting. In the coded computing
literature, a common assumption is out-of-node checksum storage, which means
that we add “coded” nodes that contain the encoded inputs. This follows the
framework in coding theory where we have systematic symbols and parity-check
symbols. However, in the HPC literature considers the in-node checksum storage
setting for practical reasons. Many HPC algorithms are designed to perform
optimally when the number of nodes is a power-of-two. If a fault-tolerant
strategy requires changing the number of nodes, it can be a less appealing
solution to practitioners. While practically appealing, in-node checksum
storage imposes challenges on the design of checksums. Specifically, since the
failure of a node now corresponds to the loss of not only input data, i.e.
systematic symbols, but also checksums, i.e. parity-check symbols, the setting
of which differs from traditional coding theory, thereby hindering direct
adoption of existing MDS coding schemes. To this end, we dedicate this Section
to derive a new coding scheme for $R$-factor protection under the in-node
checksum storage, and further show that it achieves the optimal amount of
required redundancies.

In [32], we proved a lower bound on the number of checksum blocks and a scheme
that achieves the lower bound for single failure recovery. This improved the
existing strategy for in-node checksum storage [12] by $\sim 2$x. In this
paper, we generalize this result for any number of failures $f\geq 1$. We
first show a lower bound on the number of checksum blocks $K$ and then propose
_in-node systematic MDS coding_ scheme that achieves the lower bound.

Throughout this Section, we assume that nodes are load-balanced, and we make
this assumption explicit for the checksums in the following.

######  Definition 4.1 (Load-balanced In-node Checksum Distribution).

Let $P$ be the total number of nodes and let $K$ be the number of checksum
blocks. Then, we say the checksum distribution is load-balanced if:

  * •

$\overline{a}\triangleq K-P\lfloor\frac{K}{P}\rfloor$ of the nodes have
$\lceil\frac{K}{P}\rceil$ checksum blocks, and

  * •

$\overline{b}\triangleq P-K+P\lfloor\frac{K}{P}\rfloor$ of the nodes have
$\lfloor\frac{K}{P}\rfloor$ checksum blocks.

We now show a lower bound on $K$ in the next Theorem, whose proof can be found
in Appendix D.1.

######  Theorem 4.2.

Let $L$ be the number of data blocks and assume that the total number of nodes
$P$ divides $L$. Under the load-balanced in-node checksum distribution, to
recover from any $r$ failed nodes out of $P$ nodes,

(6) $K\geq\frac{fL}{P}+f\lceil\frac{fL}{(P-f)P}\rceil.$

In the following, we show that a simple adaptation of MDS coding with $K$
satisfying (6) can tolerate any $f$ failures.

######  Construction 1 (In-Node Systematic MDS Coding).

Let $D_{1},\ldots,D_{L}$ be data blocks and $C_{1},\ldots C_{K}$ be checksum
blocks, and let $P$ be the number of nodes. We assume that $P$ divides $L$ and
we assume the load-balanced setting, i.e, each node has $L/P$ data blocks and
the checksum distribution follows the load-balanced in-node checksum
distribution. Let $\mathcal{I}_{p}$ be the set of data block indices in the
$p$-th node and let $\mathsf{Proc}(D_{i})$ or $\mathsf{Proc}(C_{i})$ denote
the processor index the data or checksum block belongs to. Let
$\widetilde{G}\in\mathbb{R}^{K\times L}$ be a checksum-generating matrix of a
$(K+L,L)$ systematic MDS code. Then, In-Node Systematic MDS Coding encodes the
$i$-th checksum block $C_{i}$ as follows:

(7)
$C_{i}=\sum_{j\notin\mathcal{I}_{\mathsf{Proc}(C_{i})}}\widetilde{G}_{i,j}D_{j}.$

######  Theorem 4.3.

In-node systematic MDS coding given in Construction 1 with $K$ that satisfies
(6) can tolerate any $f$ failed nodes out of $P$ nodes.

The proof of Theorem 4.3 is given in Appendix D.2.

######  Remark 2.

We note that directly adopting the conventional MDS coding, i.e. encoding
$C_{i}=\sum_{j=1}^{L}\widetilde{G}_{i,j}D_{j}$ instead of (7), also provides
tolerance to $f$ failures following the similar argument as the proof of
Theorem 4.3. Nevertheless, this would incur additional local computation,
though negligible, for failure recovery due to the in-node structure.
Particularly, we note that
$C_{i}-\sum_{j\in\mathcal{I}_{\mathsf{Proc}(C_{i})}}\widetilde{G}_{i,j}D_{j}=\sum_{j\notin\mathcal{I}_{\mathsf{Proc}(C_{i})}}\widetilde{G}_{i,j}D_{j}$.
Thus, if $\mathsf{Proc}(C_{i})$ is a surviving node, then it must locally
compute the quantity
$C_{i}-\sum_{j\in\mathcal{I}_{\mathsf{Proc}(C_{i})}}\widetilde{G}_{i,j}D_{j}$,
which is the numerical value of the weighted sum of other processors’ data
blocks, in the process of recovering failed nodes. To this end, our checksum
construction (7) alleviates such unnecessary step of local computation.

##  5 $Q$-Factor Protection

We have demonstrated in Section 4 that, as the resulting encoded matrix
naturally preserves the upper triangular structure of the $R$ factor, any
existing MDS checksum-generator matrix and a well-designed MDS checksum-
generator matrix would respectively suffice for the design of horizontal
checksums $AG_{h}$ (1) of out-of-node and in-node checksum storage. On the
other hand, if we apply off-the-shelf MDS encoding schemes for the vertical
checksums $G_{v}A$, it does not preserve the orthogonality of the $Q$ factor
[12, Theorem 5.1]. The main challenge in this paper is thus designing a
vertical checksum-generator matrix $G_{v}$ that can exhibit the MDS property
and preserve orthogonality at the same time.

In this section, we discuss how coding can be applied to the parallel Gram-
Schmidt algorithm to protect the left $Q$ factor (an orthogonal matrix). From
(1), we have both horizontal and vertical checksums, but in this section we
will only consider the vertical checksum for simpler presentation. When we use
$A$, one can regard it as $\begin{bmatrix}A&AG_{h}\end{bmatrix}$.

In [12, Theorem 5.1], it was concluded: _“ $Q$ in Householder QR factorization
cannot be protected by performing factorization along with the vertical
checksum.” _ The basis of this claim was that in the output retrieval
$A=Q_{1}R$ after the QR factorization of the vertically encoded matrix
$\widetilde{A}$:

(8) $\displaystyle\widetilde{A}=\begin{bmatrix}A\\\
GA\end{bmatrix}=\begin{bmatrix}Q_{1}\\\ Q_{2}\end{bmatrix}R,$

$Q_{1}$ is not orthogonal, i.e. $Q_{1}^{T}Q_{1}\neq I$. Thus, $Q_{1}R$ is not
the correct QR factorization of $A$. While the Theorem statement was proven
for the Householder algorithm, the reasoning is general and not limited to a
specific algorithm.

A major contribution of our work is showing that we can convert $Q_{1}$ into
an orthogonal matrix with a small amount of overhead (same as the overhead of
encoding in scaling sense). In Section 5.1, we prove that if the checksum-
generator matrix satisfies certain conditions (given in Theorem 5.1), there
exists a low-cost linear transform that orthogonalizes $Q_{1}$. In Section
5.2, we propose a checksum-generator matrix construction for out-of-node
checksum storage that satisfies the conditions given in Theorem 5.1 while
providing resilience to multiple-node failure. Then, in Section 5.3, we adapt
the checksum-generator matrix construction to the setting of in-node checksum
storage.

###  5.1 Low-cost post-orthogonalization

How “non-orthogonal” is $Q_{1}$? Can we still utilize $Q_{1}$ to recover the
original QR factorization? We show that with a low-cost linear transform,
$Q_{1}$ can be transformed into an orthogonal matrix, if the checksum-
generator matrix $G$ satisfies certain conditions.

######  Theorem 5.1.

Let $G_{1},V$ be submatrices of the vertical checksum-generator matrix $G$ as
follows:

(9) $G=\begin{bmatrix}G_{1}&V\end{bmatrix},$

where $G_{1}$ and $V$ have dimensions $c\times c$ and $c\times(n-c)$,
respectively. Let $G_{0}$ be an $n\times n$ by matrix as follows:

(10) $G_{0}=\begin{bmatrix}I_{c}+G_{1}&V\\\ V^{T}&-I_{n-c}\end{bmatrix}.$

If $G$ satisfies the following condition:

(11) $\displaystyle G_{1}=-\frac{1}{2}VV^{T},$

we can prove the following:

Claim 1: $G_{0}Q_{1}$ is orthogonal, i.e. $(G_{0}Q_{1})^{T}(G_{0}Q_{1})=I$.

Claim 2: $G_{0}$ is invertible.

The proof of Theorem 5.1 is given in Appendix D.3. Based on Theorem 5.1, we
next show how a checksum generator matrix satisfying (11) can be used to
retrieve the orthogonality of the original $Q$-factor. In Section 5.2, we will
design a semi-random checksum generator matrix (Construction 2) that jointly
satisfies (11) and is MDS with high probability (Theorem 5.4). Now, notice
that the matrix $G_{0}$ is very sparse as the bottom-right submatrix is simply
an $(n-c)\times(n-c)$ identity matrix, and $c$ is negligible. Claim 1 in
Theorem 5.1 suggests that by multiplying this sparse matrix $G_{0}$, we can
convert $Q_{1}$ into an orthogonal matrix. We now demonstrate how we can use
$G_{0}Q_{1}$ in place of the original $Q$ in solving a full-rank square system
of linear equations $A\bm{x}=\bm{b}$. Let $A^{\prime}=G_{0}A$. Then the
following factorization is a valid QR factorization of $A^{\prime}$:

(12) $A^{\prime}=(G_{0}Q_{1})R,$

with the left factor $(G_{0}Q_{1})$ and the right factor $R$. As $G_{0}$ is
invertible by Claim 2 in Theorem 5.1,

(13) $A\bm{x}=\bm{b}\iff(G_{0}A)\bm{x}=G_{0}\bm{b}.$

The linear system on the right side can be solved using the QR factorization
of $A^{\prime}$ given in (12). i.e.,

(14) $\displaystyle(G_{0}Q_{1})R\bm{x}$ $\displaystyle=G_{0}\bm{b},$ (15)
$\displaystyle(G_{0}Q_{1})^{T}(G_{0}Q_{1})R\bm{x}$
$\displaystyle=(G_{0}Q_{1})^{T}(G_{0}\bm{b}),$ (16) $\displaystyle R\bm{x}$
$\displaystyle=(G_{0}Q_{1})^{T}(G_{0}\bm{b}).$

Then, we can perform triangular solve to get the final answer $\bm{x}$.
Remember that we already have $Q_{1}$ and $R$ from the QR factorization of the
encoded matrix $\widetilde{A}$. Hence, all we need to perform in the above
steps is computing post-orthogonalization, $G_{0}Q_{1}$ and $G_{0}\bm{b}$. We
show in Theorem 6.1 that the overhead of post-orthogonalization is negligible.

###  5.2 Checksum-Generator Matrix for Out-of-Node Checksum Storage

The low-cost post-orthogonalization scheme exists under the constraint (11) on
the checksum-generator matrix $G$. One crucial question is whether we can
construct $G$ that has good error correction/detection capability while
satisfying (11). In this subsection, we present one such construction of $G$
for $f$ node failures recovery. Throughout this section, we assume $p_{r}$
divides $n$ for simplicity, but results generalize to any $p_{r}$ and $n$.

Recall from (2) that we can write $G=\widetilde{G}\otimes
I_{\frac{n}{p_{r}}}$. Let $\tilde{G_{1}},\tilde{V}$ be submatrices of
$\widetilde{G}$ such that
$\widetilde{G}=\begin{bmatrix}\widetilde{G}_{1}&\widetilde{V}\end{bmatrix}$,
where $\widetilde{G}_{1}$ and $\widetilde{V}$ have dimensions $m_{r}\times
m_{r}$ and $m_{r}\times(p_{r}-m_{r})$, respectively. Since
${G}=\begin{bmatrix}{G}_{1}&{V}\end{bmatrix}$ where ${G}_{1}$ and ${V}$ have
dimensions $c\times c$ and $c\times(n-c)$, we have:
$G_{1}=\widetilde{G}_{1}\otimes I_{n/p_{r}}$ and $V=\widetilde{V}\otimes
I_{n/p_{r}}$. The following Lemma establishes the equivalence of the post-
orthogonalization condition being imposed on $G$ versus $\widetilde{G}$.

######  Lemma 5.2 (Equivalence of post-orthogonalization condition).

The following relation holds:

(17) $\displaystyle
G_{1}=-\frac{1}{2}VV^{T}\iff\widetilde{G_{1}}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}.$

######  Proof 5.3.

We have:

$\displaystyle G_{1}=-\frac{1}{2}VV^{T}$
$\displaystyle\iff\widetilde{G_{1}}\otimes
I_{n/p_{r}}=-\frac{1}{2}(\widetilde{V}\otimes
I_{n/p_{r}})(\widetilde{V}^{T}\otimes I_{n/p_{r}})$
$\displaystyle\iff\widetilde{G_{1}}\otimes
I_{n/p_{r}}=(-\frac{1}{2}\widetilde{V}\widetilde{V}^{T})\otimes I_{n/p_{r}}$
$\displaystyle\iff\widetilde{G_{1}}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}.$

######  Construction 2 ($Q$-factor Checksum-Generator Matrix for Multiple-Node
Failure Recovery).

Let
$\widetilde{G}=\begin{bmatrix}\widetilde{G}_{1}&\widetilde{V}\end{bmatrix}=\\{\widetilde{g}_{i,j}\\}$
be a $f\times p_{r}$ matrix, where $\widetilde{V}$ is a $f\times(p_{r}-f)$
matrix whose entries are chosen randomly from $\mathrm{Unif}(0,1)$ and
$\widetilde{G}_{1}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}$ is an $f\times
f$ matrix.

The following $\frac{fn}{p_{r}}\times n$ checksum-generator matrix $G$
satisfies the restriction (11), and, as established in Theorem 5.4, guarantees
multiple-node fault tolerance for out-of-node checksum storage with
probability 1:

(18) $\displaystyle G=\widetilde{G}\otimes
I_{n/p_{r}}=\begin{bmatrix}\widetilde{g}_{0,0}I_{\frac{n}{p_{r}}}&\widetilde{g}_{0,1}I_{\frac{n}{p_{r}}}&\cdots&\widetilde{g}_{0,p_{r}-1}I_{\frac{n}{p_{r}}}\\\
\vdots&\ddots&&\vdots\\\
\widetilde{g}_{f-1,0}I_{\frac{n}{p_{r}}}&\widetilde{g}_{f-1,1}I_{\frac{n}{p_{r}}}&\cdots&\widetilde{g}_{f-1,p_{r}-1}I_{\frac{n}{p_{r}}}\end{bmatrix}.$

Under the _out-of-node checksum storage_ , the vertical checksums $C=GA$ are
distributed into $f\times p_{r}$ checksum nodes, and each
$\frac{n}{p_{r}}\times\frac{n}{p_{c}}$ block $C_{i,j}$ is owned by the
checksum processor $\Pi_{C}(i,j)$. Furthermore, by (18), each checksum
$C_{i,j}$ can be computed as:

(19) $\displaystyle C_{i,j}=\sum_{t=0}^{p_{r}-1}\widetilde{g}_{i,t}A_{t,j}.$

The following Theorem shows that the random checksum generator matrix
constructed in Construction 2 satisfies the MDS condition with probability 1.

######  Theorem 5.4.

Every square submatrix of the matrix $\widetilde{G}$ constructed following
Construction 2 is full rank with probability 1.

The proof of Theorem 5.4 and an example for the code construction are
respectively deferred to Appendix D.4 and C. We further empirically validate
Theorem 5.4, by comparing the smallest determinant and highest condition
number of all square submatrices of
$\widetilde{G}=\begin{bmatrix}\widetilde{G}_{1}&\widetilde{V}\end{bmatrix}$
following Construction 2 (i.e. random $\widetilde{V}$ and
$\widetilde{G}_{1}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}$ ) versus that
of the completely random matrix
$\widetilde{G}^{\prime}=\begin{bmatrix}\widetilde{G}_{1}^{\prime}&\widetilde{V}\end{bmatrix}$
where entries of $\widetilde{G}_{1}^{\prime}$ are also drawn iid from
$\mathrm{Unif}(0,1)$. We compare the least determinants and highest condition
numbers of these 2 matrices for $f=2,\ldots,6$ and $p_{r}=16,32,64$ in Figure
3. The results in the plot are the average of 100 trials222Just for
$(f=6,p_{r}=64)$, the number of trials was 10 due to the long running time..
The results show that the least determinant of square submatrices of the
structured matrix $\widetilde{G}$ is comparable to that of the completely
random matrix $\widetilde{G}^{\prime}$, and the highest condition number of
square submatrices of the structured matrix $\widetilde{G}$ is at most a
factor of $3$ from that of the completely random matrix
$\widetilde{G}^{\prime}$.

(a) Log of the least determinant of square submatrices.

(b) Log of the highest condition number of square submatrices.

Figure 3:  Comparison between the structured matrix $\widetilde{G}$ and the
completely random matrix $\widetilde{G}^{\prime}$, both of which have size
$f\times p_{r}$, in terms of the least determinant of square submatrices
(Figure 3(a)) and the highest condition number of square submatrices (Figure
3(b)). We report the results for $f=2,\ldots,6$ and $p_{r}=16,32,64$. Each
color (blue, green, red) corresponds to each value of $p_{r}$ ($=16,32,64$)
being used for comparing $\widetilde{G}$ versus $\widetilde{G}^{\prime}$.
Solid lines are for the structured matrices and dotted lines are for
completely random matrices.

###  5.3 Checksum-Generator Matrix for In-Node Checksum Storage

We highlight that the above results hold for the general in-node checksum
storage setting. Next, we demonstrate the application of our in-node
systematic MDS coding to the design of coded QR decomposition. For clarity, we
now turn back to differentiating between the vertical checksum generator
matrix $G_{v}$ and the horizontal checksum generator matrix $G_{h}$ as in (1),
and recall that the encoded input matrix $\widetilde{A}$ can be written as:

$\widetilde{A}=\begin{bmatrix}A&AG_{h}\\\ G_{v}A&G_{v}AG_{h}\end{bmatrix},$

in which, by Lemma 3.2, $G_{v}$ and $G_{h}$ respectively translate to
$Q$-factor protection and $R$-factor protection. Adapting Theorem 4.2 to our
setting of 2D block distribution, we have $L=p_{r}$ data blocks and
$P=p_{r}=p_{c}$, and consequently consider $K=f+f\lceil\frac{f}{P-f}\rceil$
that achieves equality in (6). Then the design of horizontal checksum
generator matrix $G_{h}$ can follow immediately Construction 1. On the other
hand, such direct adoption for $Q$-factor protection is not feasible, since as
discussed in Section 5.1 the vertical checksum generator matrix $G_{v}$ must
satisfy the post-orthogonalization condition (11) to circumvent the degraded
orthogonality in the decoding phase. Leveraging our idea of semi-random
checksum design in Construction 2 and the fact that the MDS coding is
applicable to in-node setting with negligible additional cost (Remark 2), we
propose the following vertical checksum generator matrix $G_{v}$ that with
high probability is MDS and thus tolerable to $f$ node failures.

######  Construction 3 (In-Node $Q$-factor Checksum-Generator Matrix for
Multiple-Node Failure Recovery).

Let
$\widetilde{G}=\begin{bmatrix}\widetilde{G}_{1}&\widetilde{V}\end{bmatrix}=\\{\widetilde{g}_{i,j}\\}$
be a $K\times p_{r}$ matrix, where $\widetilde{V}$ is a $K\times(p_{r}-K)$
matrix whose entries are chosen randomly from $\mathrm{Unif}(0,1)$ and
$\widetilde{G}_{1}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}$ is an $f\times
f$ matrix. The vertical checksum generator matrix is constructed as
$G_{v}=\widetilde{G}\otimes I_{n/p_{r}}$.

######  Theorem 5.5.

Every square submatrix of the matrix $\widetilde{G}$ constructed following
Construction 3 is full rank with probability 1.

The proof of Theorem 5.5 can be done similarly to that of Theorem 5.4.

##  6 Overhead Analysis of the $Q$-factor Protection Scheme

As discussed in Section 4, while there have been optimal coding designs and
overhead analysis for $R$-factor protection schemes, the literature for
$Q$-factor protection schemes via coding theory remains nascent, as
demonstrated in the summary of previous work in Table 1. This motivated our
proposal of our novel coding scheme for $Q$-factor protection in Section 5.
Nevertheless, even if a coding scheme can provide a great failure resilience,
if it has a prohibitive computational cost, it would not be a practical
solution. In this section, we show through a rigorous analysis that the
overhead of coding for $Q$-factor protection is negligible compared to the
cost of original computation. We analyze the overhead in terms of time
complexity that accounts for both communication and computation. A detailed
cost model is given in the following section.

###  6.1 Cost Model and MPI collective communication operations

We use the $\alpha$-$\beta$-$\gamma$ model [8, 46] to analyze communication
and computation costs:

(20) $\displaystyle T=\alpha C_{1}+\beta C_{2}+\gamma C_{3},$

where $C_{1}$ is the number of communication rounds, $C_{2}$ is the number of
bytes communicated on the critical path, and $C_{3}$ is the number of floating
point operations (flops). The $\alpha$ and $\beta$ terms model the
communication latency and bandwidth, respectively, and the $\gamma$ term is
the computation cost.

Both parallel QR decomposition and our coding strategy are implemented via MPI
collective operations, which are common communication-computation patterns in
parallel computing. Our overhead analysis depends on the costs of these
operations including broadcast, reduce and all-reduce. We denote by
$T_{broadcast}(p,w),T_{reduce}(p,w)$ and $T_{allreduce}(p,w)$ respectively the
cost of MPI_Broadcast($p,w$), MPI_Reduce($p,w$), and MPI_Allreduce($p,w$) for
$p$ processors to transfer an array of $w$ words 333It can be other units,
such as bytes, but then the term $\beta$ is the bandwidth per bytes. . We
consider the algorithms given in [47] for broadcasts, and [46] for reduce and
all-reduce, which are optimized for communication of long messages.
Illustrations and costs of these MPI operations are given in Table 2.

Table 2: MPI operations  MPI_Broadcast($p,w$) |  MPI_Reduce($p,w$) |  MPI_Allreduce($p,w$)
---|---|---
|  |
$T_{broadcast}(p,w)=(\sqrt{\alpha(\lceil\log p\rceil-1)}+\sqrt{\beta w})^{2}$ |  $T_{reduce}(p,w)=2\alpha\log p+2\beta\frac{p-1}{p}w+\gamma\frac{p-1}{p}w$ |  $T_{allreduce}(p,w)=2\alpha\log p+2\beta\frac{p-1}{p}w+\gamma\frac{p-1}{p}w$

###  6.2 Overhead Analysis

Two types of overhead are considered: the total overhead of coding
$T_{coding}$ and the overhead for recovering multiple-nodes failure
$T_{recov}$. The coding overhead $T_{coding}$ is modeled as:

(21) $\displaystyle T_{coding}=T_{enc}+T_{post}+T_{comp}$

where $T_{enc}$, $T_{post}$ and $T_{comp}$ are the overhead for encoding,
post-orthogonalization, and increased computation cost for QR factorization of
the encoded matrix. We compare $T_{coding}$ with the cost $T_{QR}$ of QR
factorization without coding for an $n\times n$ matrix.

The following Theorem, whose proof is given in Appendix D.5, summarizes the
overhead analysis:

######  Theorem 6.1.

Let $T_{QR}$ be the overall cost of the the MGS algorithm [39] for factorizing
an $n$-by-$n$ matrix on a $p_{r}$-by-$p_{c}$ processor grid. Then, the
overhead of the $Q$-factor protection scheme given in Construction 2 is given
as:

$\displaystyle T_{coding}=O\left(\frac{f}{p_{r}}+\frac{f}{p_{c}}\right)\cdot
T_{QR},\;\;\;T_{recov}=O\left(\frac{f}{p_{c}}\right)\cdot T_{QR},$

where $f$ is the number of failures. Especially,

(22) $\displaystyle T^{\alpha}_{coding}=O\left(\frac{f}{n}\right)\cdot
T^{\alpha}_{QR}.$

This Theorem shows that as $p_{r}$ and $p_{c}$ grow, the overhead of coding
becomes negligible compared to $T_{\text{QR}}$, as we assume that $f$ is a
small fixed constant. Also, (22) shows that the $\alpha$ cost, i.e., the
latency cost which is usually orders of magnitude higher than others, scales
as $O(\frac{1}{n})$.

##  7 Experimental Evaluation

Figure 4: Reverse-diagonal failure pattern for $f=2$.

In this section, we empirically demonstrate the negligible overhead of our
coded computing framework for fault-tolerant QR decomposition.

###  7.1 Implementation and Experimental Setting

For practicality, we have a minor deviation from the PBMGS (Algorithm 2)
considered in Section 3.1. Specifically, we replace the ICGS algorithm therein
with the MGS algorithm for further practical acceleration with negligible
error, on the order of $10^{-9}$ on average, as measured by $||A-QR||_{2}$.
Moreover, we make the following considerations in our implementation and cost
measurement of the post-orthogonalization step. First, as solving either
$R\bm{x}=Q_{1}^{T}\bm{b}$ or $R\bm{x}=(G_{0}Q_{1})^{T}(G_{0}\bm{b})$ (as done
in our coded computing framework in Section 5.1) requires one call of the
triangular solve operation on the original problem size and does not reflect
additional cost due to the coding framework, we neglect this step from the
cost measurement of the time $T_{QR}$ for a full QR decomposition. Yet, we
still take into account the overhead for computing $(G_{0}Q_{1})$ and
$(G_{0}\bm{b})$ in the post-orthogonalization overhead $T_{post}$. Second, for
the computation of $(G_{0}Q_{1})$ and $(G_{0}\bm{b})$, instead of deploying
the optimized procedure in Appendix D.5 that utilizes the sparsity of $G_{0}$,
we directly use the parallel matrix-matrix multiplication routine and
empirically demonstrate that the total overhead remains negligible.

All the experiments were performed on the CSC Pod cluster. Each batch used six
nodes containing dual 20-core Intel 6148s, 160 GB RAM, and an OmniPath
interconnect. At the software level, Intel MKL 2022.0.2 was used for BLAS and
LAPACKE routines, random number generation, and memory management. Intel MPI
2021.5.1 was used for message passing, and batch jobs were submitted through
SLURM.

Run-times are measured per step and do not include MPI initialization, input
matrix generation, or floating point error calculations. We report the total
overhead $T_{coding}$ of our coding framework and the recovery overhead
$T_{recov}$ in the form of the percentage increase with respect to the total
running time $T_{QR}$ of parallel QR decomposition on a $n\times n$ matrix.
All the results are averaged over 10 trials444The plots also include standard
deviation of every point, which may not be visible due to negligible
variation..

Setups and failure modeling: In all tests, we fix $n=24000$ and vary the
number of processors via $p=p_{r}=p_{c}\in\\{6,8,10,12\\}$ as well as the
number of tolerable failures (per row/column) $f\in\\{0,1,2,3\\}$. Note that
$f=0$ corresponds to a vanilla run of parallel QR decomposition without the
coding scheme. For failure modeling, we inject fail-stop errors in every
iteration of the PBMGS algorithm. In particular, at the start of every
iteration of PBMGS, we simulate the failures of at most $f$ processors on each
row and each column of the processor grid in a reverse diagonal pattern as
shown in Figure 4. This results in a total of $O(fp)$ processor failures for a
given iteration. Consequently, in every iteration, the coded QR decomposition
framework recovers all the failed processors via checksums before proceeding
with the computation.

(a) Strong scaling of $T_{coding}$ with $p=p_{r}=p_{c}$

(b) Linear relationship between $T_{coding}$ and $f$.

(c) Linear relationship between breakdown overheads (encoding cost $T_{enc}$
and total recovery cost $p\cdot T_{recov}$) and $f$.

Figure 5: Overhead Scaling for Coded QR Decomposition.

###  7.2 Overhead Analysis

We evaluate the total overhead $T_{coding}$ incurred by our coding scheme with
respect to $p$ and $f$, respectively in Figure 5(a) and 5(b). The results
empirically validate that our coding overhead (taken as a proportion of the
total QR decomposition cost $T_{QR}$) scales inversely with $p$ (Figure 5(a))
and linearly with $f$ (Figure 5(b)), thereby confirming Theorem 6.1. Moreover,
the linear trend in Figure 5(b) persists across various processor grid sizes
$p$, providing evidence that the complexity of our fault tolerance mechanism
scales predictably with the degree of redundancy.

Besides requiring less amount of redundancies for fault-tolerance than
checkpointing [6], ABFT is known to support more efficient recovery [12]. We
measure the total recovery overhead of $p\cdot T_{recov}$ across all failures
throughout the $p$ iterations of the QR decomposition, each of which incurs
$T_{recov}$ for recovering in parallel all the failures within such iteration.
We report both the total recovery overhead and the encoding overhead $T_{enc}$
in Figure 5(c) for different degrees of failures $f$ given a fixed number of
processors $p$. The total overhead $p\cdot T_{recov}$ exhibits linear
correlation with $f$ and boundedness independent from $p$, thereby implying
the growth rate of $O(\frac{f}{p})$ with respect to $T_{QR}$ as suggested by
Theorem 6.1. Moreover, we note that as the total overhead matches the encoding
overhead $T_{enc}$, the actual recovery cost per failure of our coding scheme
is minimal.

The empirical findings validate the efficiency of our coding framework in
terms of both the negligible overhead required for maintaining fault-tolerance
and recovery time, all of which result in increasing uptime for large-scale
HPC systems.

##  8 Conclusion and Discussion

In this paper, we extend our prior work [33] and propose a novel coding
strategy for parallel QR decomposition that tolerates multiple-node failures
with negligible overhead under both of the out-of-node and in-node checksum
storage settings. For $R$-factor protection, while the optimal coding schemes
for out-of-node checksum storage have been presented in the literature, we
improve the in-node checksum storage codes in [12], by first proving a
fundamental bound and proposing a strategy that achieves such lower bound and
outperforms the existing strategy [12] by $\sim 2$x. We inherit realistic
assumptions for the HPC computing environment from [12] and bring innovations
from the coding-theoretic perspective. For $Q$-factor protection, which has
only relied on checkpointing-type techniques as the straightforward
application of coding can break the orthogonality of the $Q$-factor, we
propose a coding strategy with low-cost post-orthogonalization that can
recover an orthogonal matrix at the end. Extensive experiments across
different scales of problem sizes, number of processors, and failures validate
our theories and demonstrate the negligible overhead of our coding scheme. We
believe that this work would open up exciting future directions at the
interplay between practical HPC algorithms and coded computing.

Since our scope in this work is limited to full-rank square matrix $A$, coded
QR decomposition for singular matrix or a rectangular matrix $A$ remains as an
open question. Also, designing coding strategies for other commonly used
algorithms for parallel QR decomposition, such as Householder Transformation
[38, 3] and Givens Rotation [17], would be an interesting direction for future
work.

## Acknowledgements

Use was made of computational facilities purchased with funds from the
National Science Foundation (CNS-1725797) and administered by the Center for
Scientific Computing (CSC). The CSC is supported by the California NanoSystems
Institute and the Materials Research Science and Engineering Center (MRSEC;
NSF DMR 2308708) at UC Santa Barbara.

## References

  * [1] M. Aliasgari, J. Kliewer, and O. Simeone, Coded computation against processing delays for virtualized cloud-based channel decoding, IEEE Transactions on Communications, 67 (2019), pp. 28–38.
  * [2] C. Ashcraft and R. G. Grimes, Spooles: An object-oriented sparse matrix library, in PPSC, 1999.
  * [3] G. Ballard, J. Demmel, L. Grigori, M. Jacquelin, H. D. Nguyen, and E. Solomonik, Reconstructing householder vectors from tall-skinny qr, in 2014 IEEE 28th International Parallel and Distributed Processing Symposium, May 2014, https://doi.org/10.1109/IPDPS.2014.120.
  * [4] A. Benoit, A. Cavelan, Y. Robert, and H. Sun, Assessing general-purpose algorithms to cope with fail-stop and silent errors, ACM Transactions on Parallel Computing, 3 (2016), pp. 1–36, https://doi.org/10.1145/2897189.
  * [5] G. Bosilca, R. Delmas, J. Dongarra, and J. Langou, Algorithmic Based Fault Tolerance Applied to High Performance Computing, (2008), https://arxiv.org/abs/0806.3121.
  * [6] G. Bosilca, R. Delmas, J. Dongarra, and J. Langou, Algorithm-based fault tolerance applied to high performance computing, Journal of Parallel and Distributed Computing, 69 (2009), pp. 410–416, https://doi.org/https://doi.org/10.1016/j.jpdc.2008.12.002, https://www.sciencedirect.com/science/article/pii/S0743731508002141.
  * [7] V. Cadambe and P. Grover, Codes for distributed computing: A tutorial, IEEE Information Theory Society Newsletter, 67 (2017), pp. 3–15.
  * [8] E. Chan, M. Heimlich, A. Purkayastha, and R. van de Geijn, Collective communication: Theory, practice, and experience: Research articles, Concurr. Comput. : Pract. Exper., 19 (2007), https://doi.org/10.1002/cpe.v19:13.
  * [9] Z. Chen and J. Dongarra, Algorithm-based checkpoint-free fault tolerance for parallel matrix computations on volatile resources, Proceedings 20th IEEE International Parallel & Distributed Processing Symposium, (2006).
  * [10] T. Davies, C. Karlsson, H. Liu, C. Ding, and Z. Chen, High performance linpack benchmark: A fault tolerant implementation without checkpointing, in Proceedings of the International Conference on Supercomputing, ICS ’11, New York, NY, USA, 2011, ACM, https://doi.org/10.1145/1995896.1995923.
  * [11] J. J. Dongarra, P. Luszczek, and A. Petitet, The linpack benchmark: Past, present, and future, 2002.
  * [12] P. Du, A. Bouteiller, G. Bosilca, T. Herault, and J. Dongarra, Algorithm-based fault tolerance for dense matrix factorizations, in Proceedings of the 17th ACM SIGPLAN Symposium on Principles and Practice of Parallel Programming, PPoPP ’12, New York, NY, USA, 2012, Association for Computing Machinery, p. 225–234, https://doi.org/10.1145/2145816.2145845, https://doi.org/10.1145/2145816.2145845.
  * [13] P. Du, P. Luszczek, S. Tomov, and J. Dongarra, Soft error resilient qr factorization for hybrid system with gpgpu, in Proceedings of the Second Workshop on Scalable Algorithms for Large-scale Systems, ScalA ’11, New York, NY, USA, 2011, ACM, https://doi.org/10.1145/2133173.2133179.
  * [14] S. Dutta, V. Cadambe, and P. Grover, Short-dot: Computing large linear transforms distributedly using coded short dot products, in Advances In Neural Information Processing Systems, 2016.
  * [15] S. Dutta, M. Fahim, F. Haddadpour, H. Jeong, V. Cadambe, and P. Grover, On the optimal recovery threshold of coded matrix multiplication, IEEE Transactions on Information Theory, (2019).
  * [16] S. Dutta, M. Fahim, F. Haddadpour, H. Jeong, V. Cadambe, and P. Grover, On the optimal recovery threshold of coded matrix multiplication, IEEE Transactions on Information Theory, 66 (2020), pp. 278–301.
  * [17] O. Egecioglu, Givens and householder reductions for linear least squares on a cluster of workstations, tech. report, Santa Barbara, CA, USA, 1995.
  * [18] E. Elnozahy, R. Bianchini, T. El-Ghazawi, A. Fox, F. Godfrey, A. Hoisie, K. McKinley, R. Melhem, J. Plank, and P. Ranganathan, System resilience at extreme scale. Defense Advanced Research Project Agency (DARPA), (2008).
  * [19] N. Ferdinand and S. C. Draper, Hierarchical coded computation, in IEEE International Symposium on Information Theory (ISIT), 2018, pp. 1620–1624.
  * [20] N. S. Ferdinand and S. C. Draper, Anytime coding for distributed computation, in Communication, Control, and Computing (Allerton), 2016, pp. 954–960.
  * [21] P. Fitzpatrick and C. C. Murphy, Fault tolerant matrix triangularization and solution of linear systems of equations, in [1992] Proceedings of the International Conference on Application Specific Array Processors, Aug 1992, https://doi.org/10.1109/ASAP.1992.218551.
  * [22] G. H. Golub and C. F. van Loan, Matrix Computations, JHU Press, fourth ed., 2013, http://www.cs.cornell.edu/cv/GVL4/golubandvanloan.htm.
  * [23] T. Herault and Y. Robert, Fault-Tolerance Techniques for High-Performance Computing, 01 2015, https://doi.org/10.1007/978-3-319-20943-2.
  * [24] H. Jeong, T. Low, and P. Grover, Masterless coded computing: A fully-distributed coded fft algorithm, Communication, Control, and Computing (Allerton), (2018).
  * [25] D. D. Joshi, A note on upper bounds for minimum distance codes, (1958).
  * [26] Y. Kim, J. Plank, and J. Dongarra, Fault tolerant matrix operations using checksum and reverse computation, in Proceedings of 6th Symposium on the Frontiers of Massively Parallel Computation (Frontiers ’96), 1996, pp. 70–77, https://doi.org/10.1109/FMPC.1996.558063.
  * [27] K. Lee, M. Lam, R. Pedarsani, D. Papailiopoulos, and K. Ramchandran, Speeding up distributed machine learning using codes, IEEE Transactions on Information Theory, (2017).
  * [28] F. T. Luk and H. Park, An analysis of algorithm-based fault tolerance techniques, Journal of Parallel and Distributed Computing, 5 (1988), https://doi.org/https://doi.org/10.1016/0743-7315(88)90027-5.
  * [29] F. J. MacWilliams and N. J. A. Sloane, The theory of error-correcting codes, Elsevier, 1977.
  * [30] A. Mallick, M. Chaudhari, and G. Joshi, Rateless codes for near-perfect load balancing in distributed matrix-vector multiplication, arXiv preprint arXiv:1804.10331, (2018).
  * [31] O. Maslennikow, J. Kaniewski, and R. Wyrzykowski, Fault tolerant qr-decomposition algorithm and its parallel implementation, in Euro-Par’98 Parallel Processing, D. Pritchard and J. Reeve, eds., Berlin, Heidelberg, 1998, Springer Berlin Heidelberg.
  * [32] Q. M. Nguyen, H. Jeong, and P. Grover, Coded QR decomposition, in 2020 IEEE International Symposium on Information Theory (ISIT), 2020, pp. 191–196, https://doi.org/10.1109/ISIT44484.2020.9173985.
  * [33] Q. M. Nguyen, H. Jeong, and P. Grover, Coded QR Decomposition, in IEEE International Symposium on Information Theory (ISIT), 2020.
  * [34] S. Oliveira, L. Borges, M. Holzrichter, and T. Soma, Analysis of different partitioning schemes for parallel gram-schmidt algorithms, Parallel Algorithms Appl., 14 (2000), https://doi.org/10.1080/10637199808947392.
  * [35] H. Park, On multiple error detection in matrix triangularizations using checksum methods, Journal of Parallel and Distributed Computing, 14 (1992), https://doi.org/https://doi.org/10.1016/0743-7315(92)90100-2.
  * [36] N. Raviv, R. Tandon, A. Dimakis, and I. Tamo, Gradient coding from cyclic mds codes and expander graphs, in International Conference on Machine Learning (ICML), 2018, pp. 4302–4310.
  * [37] A. Reisizadeh, S. Prakash, R. Pedarsani, and S. Avestimehr, Coded computation over heterogeneous clusters, in Information Theory (ISIT), 2017 IEEE International Symposium on, IEEE, 2017.
  * [38] F. Rotella and I. Zambettakis, Block householder transformation for parallel qr factorization, Applied Mathematics Letters, 12 (1999), pp. 29 – 34, https://doi.org/https://doi.org/10.1016/S0893-9659(99)00028-2, http://www.sciencedirect.com/science/article/pii/S0893965999000282.
  * [39] G. Rünger and M. Schwind, Comparison of different parallel modified gram-schmidt algorithms, vol. 3648, 08 2005, https://doi.org/10.1007/11549468_90.
  * [40] G. Rünger and M. Schwind, Comparison of different parallel modified gram-schmidt algorithms, vol. 3648, 08 2005, https://doi.org/10.1007/11549468_90.
  * [41] A. Severinson, A. G. i Amat, and E. Rosnes, Block-Diagonal and LT Codes for Distributed Computing With Straggling Servers, IEEE Transactions on Communications, 67 (2019), pp. 1739–1753.
  * [42] T. Sterling, M. Anderson, and M. Brodowicz, High Performance Computing: Modern Systems and Practices, Morgan Kaufmann Publishers Inc., San Francisco, CA, USA, 1st ed., 2017.
  * [43] K. Swirydowicz, J. Langou, S. Ananthan, U. Yang, and S. Thomas, Low synchronization gmres algorithms, 2018, https://arxiv.org/abs/1809.05805.
  * [44] R. Tandon, Q. Lei, A. G. Dimakis, and N. Karampatziakis, Gradient coding, arXiv preprint arXiv:1612.03301, (2016).
  * [45] R. Tandon, Q. Lei, A. G. Dimakis, and N. Karampatziakis, Gradient Coding: Avoiding Stragglers in Distributed Learning, in International Conference on Machine Learning (ICML), 2017, pp. 3368–3376.
  * [46] R. Thakur, R. Rabenseifner, and W. Gropp, Optimization of collective communication operations in mpich, Int. J. High Perform. Comput. Appl., 19 (2005), https://doi.org/10.1177/1094342005051521.
  * [47] J. L. Träff and A. Ripke, Optimal broadcast for fully connected processor-node networks, J. Parallel Distrib. Comput., 68 (2008), p. 887–901, https://doi.org/10.1016/j.jpdc.2007.12.001, https://doi.org/10.1016/j.jpdc.2007.12.001.
  * [48] L. N. Trefethen and D. Bau III, Numerical linear algebra, vol. 50, Siam, 1997.
  * [49] S. Wang, J. Liu, and N. Shroff, Coded sparse matrix multiplication, in International Conference on Machine Learning (ICML), 2018, pp. 5139–5147.
  * [50] P. Wu and Z. Chen, Ft-scalapack: correcting soft errors on-line for scalapack cholesky, qr, and lu factorization routines, in HPDC, 2014.
  * [51] E. Yao, J. Zhang, M. Chen, G. Tan, and N. Sun, Detection of soft errors in LU decomposition with partial pivoting using algorithm-based fault tolerance, The International Journal of High Performance Computing Applications, 29 (2015), https://doi.org/10.1177/1094342015578487.
  * [52] Q. Yu, M. A. Maddah-Ali, and A. S. Avestimehr, Coded fourier transform, in 2017 55th Annual Allerton Conference on Communication, Control, and Computing (Allerton), IEEE, 2017, pp. 494–501.
  * [53] Q. Yu, M. A. Maddah-Ali, and A. S. Avestimehr, Polynomial Codes: an Optimal Design for High-Dimensional Coded Matrix Multiplication, in Advances In Neural Information Processing Systems (NIPS), 2017.

## APPENDIX

##  Appendix A The PBMGS Algorithm and Its Overhead

Input: $n\times n$ matrix $A$ distributed 2D-block cyclically over
$P=p_{r}\times p_{c}$ processors.

Output: $Q,R$

Result: $A=QR$, where $Q^{T}Q=I$

1 Q=A

2for  _$i=1$ to $n$ Step $b$ _ do

3 $\bar{Q}=Q[:,i:i+b-1]$

4 $ICGS(\bar{Q},R[i:i+b-1,i:i+b-1])$

5 $Q[:,i:i+b-1]=\bar{Q}$

6 if  _$i <n-b$ _ then

7

8 All-reduce: $\bar{R}=\bar{Q}^{T}\cdot Q[:,i+b:n]$

9 $R[i:i+b-1,i:i+b-1]=\bar{R}$

10 $Q[:,i+b:n]-=\bar{Q}\cdot\bar{R}$

11 end if

12

13 end for

Algorithm 2 PBMGS($A,n,n,p_{r},p_{c}$)

We consider the parallel blocked MGS algorithm (PBMGS) in [39], which is
specifically designed for 2D-block cyclic distribution. However, our analysis
can be adapted to other algorithms in [17, 43]. At high level, PBMGS
(Algorithm 2) iterates through $p_{c}=\frac{n}{b}$ block-columns. At the
$i^{th}$ iterations it performs the iterated classical Gram-Schmidt algorithm
(ICGS) on the whole block-column, and updates the remaining part (on the right
side of the block-column) of the matrix.

The ICGS algorithm is essentially the classical Gram-Schmidt (CGS), except
that in each iteration ICGS further reorthogonalizes for one more time,
depending on an easy-to-compute criterion [40]. The cost of ICGS algorithm is
thus lower bounded by the cost of CGS algorithm (Algorithm 3). To formally
analyze the overhead of the PBMGS routine, we denote respectively by:

  * •

$T_{CGS}(A,n,b,p_{r})$ the cost of performing CGS algorithm on the $n\times b$
matrix $A$ distributed 1D-block cyclically over $p_{r}$ processors,

  * •

and $T_{PMGS}(A,n,n,p_{r},p_{c})$ the cost of performing PMGSGS algorithm on
the $n\times n$ matrix $A$ distributed 2D-block cyclically over $P=p_{r}\times
p_{c}$ processors.

In the $i^{th}$ iteration of algorithm $CGS(A,n,b,p_{r})$, the computation of
the $(b-i+1)\times 1$ vector $\bm{v}$ requires the local matrix-vector
multiplication on each of the $p_{r}$ processors incurring $\gamma
2\frac{n}{p_{r}}\cdot(b-i+1)$, and one MPI_Reduce call with overhead
$T_{reduce}(p_{r},b-i+1)$. Summing up over $b$ iterations, the cost of the CGS
algorithm can be lower bounded by:

$\displaystyle T_{CGS}(A,n,b,p_{r})$
$\displaystyle\geq\sum^{b}_{i=1}[T_{reduce}(p_{r},b-i+1)+\gamma
2\frac{n}{p_{r}}(b-i+1)]$ $\displaystyle=\sum^{b}_{i=1}[2\alpha\log
p_{r}+\beta\frac{p_{r}-1}{p_{r}}(b-i+1)$
$\displaystyle+\gamma\frac{p_{r}-1}{p_{r}}(b-i+1)]+\gamma\frac{nb(b+1)}{p_{r}}.$

Hence,

(23) $\displaystyle T_{CGS}(A,n,b,p_{r})$ $\displaystyle\geq 2\alpha b\log
p_{r}+\gamma\frac{nb(b+1)}{p_{r}}.$

In the $i^{th}$ iteration of algorithm $PBMGS(A,n,n,p_{r},p_{c})$, the
computation of the $b\times(n-i-b+1)$ matrix $\bar{R}$ requires the following
main steps (we do not list all steps ) and communication patterns:

  1. 1.

ICGS algorithm to update $\bar{Q}$, incurring $\geq T_{CGS}(A,n,b,p_{r})$
overhead.

  2. 2.

Horizontal broadcast of $\bar{Q}$ via MPI_Broadcast (in preparation for next
step), incurring $T_{broadcast}(p_{c},\frac{n}{p_{r}}b)$ overhead.

  3. 3.

Local matrix-matrix multiplication: each of the $p_{r}$ processors prepares
its local inner-products through a matrix-matrix multiplication of the local
parts of $\bar{Q}$ (transposed) with its local parts of $Q$. This incurs
$\gamma 2b\frac{n}{p_{r}}\cdot\frac{n-i-b+1}{p_{c}}=\gamma
2\frac{n}{P}b(n-i-b+1)$ overhead.

  4. 4.

MPI_Reduce call with overhead $T_{reduce}(p_{r},b(n-i-b+1))$ over $p_{r}$
vertical processors on each block-column. (Note: The local matrix-matrix
multiplication in step 3 results in a $b\times(n-i-b+1)$ submatrix.)

Summing up over all $\frac{n}{b}$ iterations (with step $b$), the cost of the
PBMGS algorithm can be lower bounded by:

$\displaystyle T_{PBMGS}(A,n,n,p_{r},p_{c})$
$\displaystyle\geq\frac{n}{b}T_{CGS}(A,n,b,p_{r})+\frac{n}{b}T_{broadcast}(p_{c},\frac{n}{p_{r}}b)$
$\displaystyle+\sum^{\frac{n}{b}-1}_{\begin{subarray}{c}j=0\\\
i=1+bj\end{subarray}}[T_{reduce}(p_{r},b(n-i-b+1))+\gamma
2\frac{n}{P}b(n-i-b+1)]$ $\displaystyle\geq\frac{n}{b}[2\alpha b\log
p_{r}+\gamma\frac{nb(b+1)}{p_{r}}]+\frac{n}{b}(\alpha\log
p_{c}+\beta\frac{p_{c}-1}{p_{c}}\frac{n}{p_{r}}b)$
$\displaystyle+\sum^{\frac{n}{b}-1}_{\begin{subarray}{c}j=0\\\
i=1+bj\end{subarray}}[\gamma 2\frac{n}{P}b(n-i-b+1)]\text{ (by
(\ref{CGSbound1}))}$ $\displaystyle\geq[2\alpha n\log
p_{r}+\gamma\frac{n^{2}(b+1)}{p_{r}}]+[\beta\frac{(p_{c}-1)n^{2}}{P}]$
$\displaystyle+\gamma(\frac{n^{3}}{P}-\frac{n^{2}b}{P})$
$\displaystyle=2\alpha n\log
p_{r}+\beta\frac{(p_{c}-1)n^{2}}{P}+\gamma(\frac{n^{3}}{P}+\frac{n^{2}(b+1)p_{c}}{P}-\frac{n^{2}b}{P})$
$\displaystyle\geq 2\alpha n\log
p_{r}+\beta\frac{1}{2}\frac{p_{c}n(n+1)}{P}+\gamma\frac{n^{2}(n+1)}{P}.$

Therefore,

(24) $\displaystyle T_{QR}\geq 2\alpha n\log
p_{r}+\beta\frac{1}{2}\frac{p_{c}n(n+1)}{P}+\gamma\frac{n^{2}(n+1)}{P}.$

Input: $n\times b$ matrix $A$ distributed 1D-block cyclically over $p_{r}$
processors.

Output: $Q,R$

Result: $A=QR$, where $Q^{T}Q=I$

1 Q=A

2for  _$i=1$ to $b$ _ do

3 $\bm{q}=Q[:,i]$

4 All-reduce: vector $\bm{v}=Q[:,i:b]^{T}\cdot Q[:,i]$

5 $\bm{q}=\bm{q}/\sqrt{\bm{v}[1]}$

6 $Q[:,i]=\bm{q}$

7 $R[i:,i:b]=\bm{v}^{T}$

8 if  _$i <b$ _ then

9 $Q[:,i+1:b]-=\bm{q}\cdot\bm{v}[i+1:b]^{T}$

10 end if

11

12 end for

Algorithm 3 CGS($A,n,b,p_{r}$)

######  Observation 1.

For the MGS algorithm [39], the overhead $T_{comp}$ for running QR
decomposition on the larger $(m+c)\times n$ encoded matrix $\widetilde{A}$
instead of $m\times n$ matrix $A$ is bounded by:

(25) $\displaystyle T_{comp}$ $\displaystyle\leq\frac{c}{m}\cdot T_{QR}$

Note: We hereby derive the generalized result for $m\times n$ matrix A, and
later set $m=n$ in our overhead analysis.

######  Proof A.1 (Proof of Observation 1).

As any algorithm in [17, 43, 39] iterates from left to right over the width of
the matrix, there are $O(n)$ such iterations. In each iteration, all the
mentioned algorithms use some form of MPI call (in Appendix 6.1) excluding,
MPI_Alltoall. In particular, MPI_Broadcast is required to communicate the data
blocks, or MPI_Reduce/MPI_Allreduce is used to compute inner-products or
matrix-products. As each processor owns $\frac{mn}{P}$ data, the data
transferred in any such MPI call does not exceed $O(\frac{mn}{P})$. Any such
MPI call (excluding MPI_Alltoall) by Appendix 6.1 would incur $O(\alpha\log
P)+O(\beta\frac{mn}{P})$ communication. In terms of computation, all these
efficient algorithms do not exceed $\gamma O(\frac{mn}{P})$ overhead.
Therefore, summing up over $O(n)$ iterations, we can get the asymptotical cost
for $T_{QR}$ as:

(26) $\displaystyle T_{QR}=O(\alpha n\log
P+\beta\frac{mn^{2}}{P}+\gamma\frac{mn^{2}}{P})$

As $m$ increases at most linearly with the asymptotical cost of $T_{QR}$, the
overhead $T_{comp}$ for running QR decomposition on some $(m+c)\times n$
matrix against some $m\times n$ matrix can be bounded by:

$\displaystyle T_{comp}\leq\frac{c}{m}T_{QR}$

######  Remark A.2.

In (26), we also note that $T^{\alpha}_{comp}=O(1)$, as the increase of $m$
does not affect the $\alpha$ term.

##  Appendix B Checksum-preservation

We now proceed to prove the checksum-preservation, equations (3) and (4) in
Lemma 3.2, for sequential MGS (Appendix B.1) and parallel MGS (Appendix B.2)
respectively. Then the proof of Corollary 3.3 is given in Appendix B.3.

###  B.1 Checksum-preservation for sequential MGS

We first consider the columns of the encoded matrix $\widetilde{A}$:

$\displaystyle\widetilde{A}=\begin{bmatrix}\bm{\widetilde{a}_{1}}&\bm{\widetilde{a}_{2}}&\cdots&\bm{\widetilde{a}_{n+d}}\end{bmatrix}\text{,
where $\bm{\widetilde{a}_{i}}=\begin{bmatrix}\bm{a_{i}}\\\
G_{v}\bm{a_{i}}\end{bmatrix}$},$

and the columns of the final $Q$-factor $Q^{(t)}$ at the $t^{th}$ iteration
for $t=1\rightarrow T$:

$\displaystyle
Q^{(t)}=\begin{bmatrix}\bm{q_{1}}^{(t)}&\bm{q_{2}}^{(t)}&\cdots&\bm{q_{n+d}}^{(t)}\end{bmatrix},$
$\displaystyle
Q_{1}^{(t)}=\begin{bmatrix}\bm{q_{11}}^{(t)}&\bm{q_{12}}^{(t)}&\cdots&\bm{q_{1(n+d)}}^{(t)}\end{bmatrix},$
$\displaystyle
Q_{2}^{(t)}=\begin{bmatrix}\bm{q_{21}}^{(t)}&\bm{q_{22}}^{(t)}&\cdots&\bm{q_{2(n+d)}}^{(t)}\end{bmatrix},$
where $\bm{q_{i}}^{(t)}=\begin{bmatrix}\bm{q_{1i}}^{(t)}\\\
\bm{q_{2i}}^{(t)}\end{bmatrix}$.

For completeness, we now present the sequential MGS algorithm on
$\widetilde{A}$ in Algorithm 4.

Input:
$\widetilde{A}=\begin{bmatrix}\bm{\widetilde{a}_{1}}&\bm{\widetilde{a}_{2}}&\cdots&\bm{\widetilde{a}_{n+d}}\end{bmatrix}$
is $(n+c)\times(n+d)$ matrix

Output:
$Q=\begin{bmatrix}\bm{q_{1}}&\bm{q_{2}}&\dots&\bm{q_{n+d}}\end{bmatrix}$ is
$(n+c)\times(n+d)$ matrix , and $R=(r_{ij})$ is $(n+d)\times(n+d)$ matrix

Result: $\widetilde{A}=QR$, where $Q^{T}Q=I$

1 for  _$i=1$ to $n$ _ do

2 $\bm{u_{i}}=\bm{\widetilde{a}_{i}}$

3 end for

4for  _$i=1$ to $n+d$ _ do

5 $r_{ii}=||\bm{u_{i}}||_{2}$

6 $\bm{q_{i}}=\bm{u_{i}}/r_{ii}$

7 for  _$j=i+1$ to $n+d$ _ do

8

9 $r_{ij}=\bm{q_{i}}^{T}\bm{u_{j}}$

10 $\bm{u_{j}}=\bm{u_{j}}-r_{ij}\bm{q_{i}}$

11 end for

12

13 end for

Algorithm 4 Sequential MGS

In the Algorithm 4, we use the temporary $(n+c)\times(n+d)$ matrix
$U=\begin{bmatrix}\bm{u_{1}}&\bm{u_{2}}&\cdots&\bm{u_{n+d}}\end{bmatrix}$. For
analysis, we consider its value $U^{(t)}$ at the end of each iteration
$t=1\rightarrow T$. $U$ is initialized to $U^{(0)}=\widetilde{A}$ by the loop
from line 1 to line 3. We further consider the breakdowns of the intermediate
matrices as follows:

$\displaystyle U^{(t)}=\begin{bmatrix}U^{(t)}_{1}\\\
U^{(t)}_{2}\end{bmatrix}=\begin{bmatrix}\bm{u_{1}}&\bm{u_{2}}&\cdots&\bm{u_{n+d}}\end{bmatrix}\text{,
where $U^{(t)}_{1}:n\times(n+d)$, $U^{(t)}_{2}:c\times(n+d)$,}$ $\displaystyle
U^{(t)}_{1}=\begin{bmatrix}\bm{u_{11}}^{(t)}&\bm{u_{12}}^{(t)}&\cdots&\bm{u_{1(n+d)}}^{(t)}\end{bmatrix},$
$\displaystyle
U^{(t)}_{2}=\begin{bmatrix}\bm{u_{21}}^{(t)}&\bm{u_{22}}^{(t)}&\cdots&\bm{u_{2(n+d)}}^{(t)}\end{bmatrix}\text{,
where $\bm{u_{i}}^{(t)}=\begin{bmatrix}\bm{u_{1i}}^{(t)}\\\
\bm{u_{2i}}^{(t)}\end{bmatrix}$},$ $\displaystyle
R^{(t)}=\begin{bmatrix}R^{(t)}_{1}&R^{(t)}_{2}\end{bmatrix}=(r^{(t)}_{ij})=\begin{bmatrix}\bm{r_{1}}^{(t)}\\\
\vdots\\\ \bm{r_{n+c}}^{(t)}\end{bmatrix},$ $\displaystyle
R^{(t)}_{1}=\begin{bmatrix}\bm{\bar{r}_{11}}^{(t)}\\\ \vdots\\\
\bm{\bar{r}_{(n+c)1}}^{(t)}\end{bmatrix},R^{(t)}_{2}=\begin{bmatrix}\bm{\bar{r}_{12}}^{(t)}\\\
\vdots\\\ \bm{\bar{r}_{(n+c)2}}^{(t)}\end{bmatrix}\text{, where
$\bm{r_{i}}^{(t)}=\begin{bmatrix}\bm{r_{i1}}^{(t)}&\bm{r_{i2}}^{(t)}\end{bmatrix}$}.$

For the sequential MGS (Algorithm 4), there are $T=n+d$ iterations in the main
loop from line 4 to line 11. We proceed to prove by induction that for
$t=0\rightarrow T$, we have $U^{(t)}_{2}=G_{v}U^{(t)}_{1}$,
$Q_{2}^{(t)}=G_{v}Q_{1}^{(t)}$ and $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$.

Base case: For $t=0$:
As $Q$ and $U$ are initialized to
$Q^{(0)}=U^{(0)}=\widetilde{A}=\begin{bmatrix}A&AG_{h}\\\
G_{v}A&G_{v}AG_{h}\end{bmatrix}$, we have $Q_{2}^{(0)}=G_{v}Q_{1}^{(0)}$ and
$U_{2}^{(0)}=G_{v}U_{1}^{(0)}$. As $R$ is initialized to $R^{(0)}=0$, we have
$R_{2}^{(0)}=R_{1}^{(0)}G_{h}=0$.

Inductive step: For the inductive step from $t-1$ to $t$ ($t\geq 1$):
Only the column $\bm{q_{t}}^{(t)}$ of $Q^{(t)}$ is updated on line 6:

$\displaystyle\bm{q_{t}}^{(t)}$
$\displaystyle=\bm{u_{t}}^{(t-1)}/r_{tt}=\frac{1}{r_{tt}}\begin{bmatrix}\bm{u_{1t}}^{(t-1)}\\\
\bm{u_{2t}}^{(t-1)}\end{bmatrix}$
$\displaystyle=\frac{1}{r_{tt}}\begin{bmatrix}\bm{u_{1t}}^{(t-1)}\\\
G_{v}\bm{u_{1t}}^{(t-1)}\end{bmatrix}\text{ (as
$U_{2}^{(t-1)}=G_{v}U_{1}^{(t-1)}$)}$
$\displaystyle=\begin{bmatrix}\bm{u_{1t}}^{(t-1)}/r_{tt}\\\
G_{v}(\bm{u_{1t}}^{(t-1)}/r_{tt})\end{bmatrix}=\begin{bmatrix}\bm{q_{1t}}^{(t)}\\\
\bm{q_{2t}}^{(t)}\end{bmatrix}$

We thus obtain: $\bm{q_{2t}}^{(t)}=G_{v}\bm{q_{1t}}^{(t)}$. Other columns of
$Q$ remain the same: $\bm{q_{j}}^{(t)}=\bm{q_{j}}^{(t-1)}$, so
$\bm{q_{2j}}^{(t)}=G\bm{q_{1j}}^{(t)}$ for all other $j\neq t$. We conclude
that $Q_{2}^{(t)}=G_{v}Q_{1}^{(t)}$.
For the matrix $U$, there are $n-t$ columns $\bm{u_{j}}$ updated (
$j=t+1\rightarrow n$) from line 7 to line 9:

$\displaystyle\bm{u_{j}}^{(t)}$
$\displaystyle=\bm{u_{j}}^{(t-1)}-r_{tj}\bm{q_{j}}^{(t)}=\begin{bmatrix}\bm{u_{1j}}^{(t-1)}\\\
\bm{u_{2j}}^{(t-1)}\end{bmatrix}-r_{tj}\begin{bmatrix}\bm{q_{1j}}^{(t)}\\\
\bm{q_{2j}}^{(t)}\end{bmatrix}$
$\displaystyle=\begin{bmatrix}\bm{u_{1j}}^{(t-1)}\\\
G_{v}\bm{u_{1j}}^{(t-1)}\end{bmatrix}-r_{tj}\begin{bmatrix}\bm{q_{1j}}^{(t)}\\\
G_{v}\bm{q_{1j}}^{(t)}\end{bmatrix}$ (as $U_{2}^{(t-1)}=G_{v}U_{1}^{(t-1)}$
and $Q_{2}^{(t)}=G_{v}Q_{1}^{(t)}$)
$\displaystyle=\begin{bmatrix}\bm{u_{1j}}^{(t-1)}-r_{tj}\bm{q_{1j}}^{(t)}\\\
G_{v}(\bm{u_{1j}}^{(t-1)}-r_{tj}\bm{q_{1j}}^{(t)})\end{bmatrix}=\begin{bmatrix}\bm{u_{1j}}^{(t)}\\\
\bm{u_{2j}}^{(t)}\end{bmatrix}$

We thus obtain: $\bm{u_{2t}}^{(t)}=G_{v}\bm{u_{1t}}^{(t)}$. Other columns of
$U$ remain the same: $\bm{u_{j}}^{(t)}=\bm{u_{j}}^{(t-1)}$, so
$\bm{u_{2j}}^{(t)}=G_{v}\bm{u_{1j}}^{(t)}$ for all other $j\leq t$. We
conclude that $U_{2}^{(t)}=G_{v}U_{1}^{(t)}$.
Now, it is left to show that $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$. We note that only
the row $\bm{r_{t}}^{(t)}$ of $R^{(t)}$ is updated on line 8, and have the
following property.

######  Observation 2.

At the end of iteration $t$, we have:

(27) $\displaystyle\bm{r_{t}}^{(t)}$
$\displaystyle=\bm{q_{t}}^{T}U^{(t)}=\bm{q_{t}}^{T}\widetilde{A}$

Note that here $\bm{q_{t}}$ refers to $\bm{q_{t}}^{(t)}$. However, as the
value of $\bm{q_{t}}$ is last updated in the $t^{th}$ iteration and becomes
the final value afterward, we neglect the superscript $(t)$ for simplicity.

######  Proof B.1.

We first prove that $\bm{r_{t}}^{(t)}=\bm{q_{t}}^{T}U^{(t)}$. By line 5 and 8,
$r^{(t)}_{tj}=\bm{q_{t}}^{T}\bm{u_{j}}^{(t)}$ for $j\geq t$, so it is left to
show that $r^{(t)}_{tj}=\bm{q_{t}}^{T}\bm{u_{j}}^{(t)}$ also for $j<t$. We
note that for $j<t$:

  * •

$r^{(t)}_{tj}=0$ (as $r^{(t)}_{tj}$ is not updated by the algorithm and thus
maintains the initial value $0$).

  * •

$\bm{q_{j}}^{(j)}=\bm{u_{j}}^{(j)}/r_{jj}$, updated on line 8 in the previous
$j^{th}$ iteration.

  * •

$\bm{q_{j}}^{(t)}=\bm{q_{j}}^{(j)}=\bm{q_{j}}$, as $\bm{q_{j}}^{(j)}$ is last
updated in the $j^{th}$ iteration, and becomes the final value of
$\bm{q_{j}}$.

  * •

$\bm{u_{j}}^{(t)}=\bm{u_{j}}^{(j)}=\bm{u_{j}}$, as $\bm{u_{j}}^{(j)}$ is last
updated in the $j^{th}$ iteration, and becomes the final value of
$\bm{u_{j}}$.

For $j<t$, we thus obtain that
$\bm{q_{t}}^{T}\bm{u_{j}}^{(t)}=r_{jj}\bm{q_{t}}^{T}\bm{q_{j}}$. As
$Q=\begin{bmatrix}\bm{q_{1}}&\bm{q_{2}}&\dots&\bm{q_{n+d}}\end{bmatrix}$ is
orthogonal, we have $I=Q^{T}Q=(\bm{q_{i}}^{T}\bm{q_{j}})_{ij}$, implying
$\bm{q_{t}}^{T}\bm{q_{j}}=0$. Therefore,
$\bm{q_{t}}^{T}\bm{u_{j}}^{(t)}=0=r^{(t)}_{tj}$.
To complete the proof of (27), we proceed to prove that
$\bm{q_{t}}^{T}U^{(t)}=\bm{q_{t}}^{T}\widetilde{A}$. This is equivalent to
proving that
$\bm{q_{t}}^{T}\bm{u_{i}}^{(t)}=\bm{q_{t}}^{T}\bm{\widetilde{a}_{i}}$ for any
$i\in[1,n+d]$. Each $\bm{u_{i}}^{(t)}$ is only updated on line 9 in the first
$i$ iterations, and thus formulated as:

$\displaystyle\bm{u_{i}}^{(t)}$
$\displaystyle=\bm{\widetilde{a}_{i}}-\sum_{l=1}^{\min(i,t)-1}r_{lj}\bm{q_{l}}^{(l)}$
$\displaystyle=\bm{\widetilde{a}_{i}}-\sum_{l=1}^{\min(i,t)-1}r_{lj}\bm{q_{l}}$
(as $\bm{q_{l}}^{(l)}$ is last updated in the $l^{th}$ iteration.)

Left multiplying $\bm{q_{t}}^{T}$, we obtain:

$\displaystyle\bm{q_{t}}^{T}\bm{u_{i}}^{(t)}$
$\displaystyle=\bm{q_{t}}^{T}\bm{\widetilde{a}_{i}}-\sum_{l=1}^{\min(i,t)-1}r_{lj}\bm{q_{t}}^{T}\bm{q_{l}}$
$\displaystyle=\bm{q_{t}}^{T}\bm{\widetilde{a}_{i}},$

where the last line holds because $I=Q^{T}Q=(\bm{q_{i}}^{T}\bm{q_{j}})_{ij}$,
implying $\bm{q_{t}}^{T}\bm{q_{j}}=0$ for any $j\neq t$, and $l<\min(i,t)\leq
t$.

Back to our main proof, we now have:

$\displaystyle\bm{r_{t}}^{(t)}$
$\displaystyle=\bm{q_{t}}^{T}\widetilde{A}=\bm{q_{t}}^{T}\begin{bmatrix}A&AG_{h}\\\
G_{v}A&G_{v}AG_{h}\end{bmatrix}$
$\displaystyle=\begin{bmatrix}\bm{q_{t}}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}&\bm{q_{t}}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}G_{h}\end{bmatrix}$
$\displaystyle=\begin{bmatrix}\bm{r_{t1}}^{(t)}&\bm{r_{t2}}^{(t)}\end{bmatrix}$
$\displaystyle\implies$
$\displaystyle\begin{cases}\bm{r_{t1}}^{(t)}=\bm{q_{t}}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}\\\ \bm{r_{t2}}^{(t)}=\bm{q_{t}}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}G_{h}\end{cases}$

We thus obtain: $\bm{r_{t2}}^{(t)}=\bm{r_{t1}}^{(t)}G_{h}$. Other rows of $R$
remain the same: $\bm{r_{j}}^{(t)}=\bm{r_{j}}^{(t-1)}$, so
$\bm{r_{j2}}^{(t)}=\bm{r_{j1}}^{(t)}G_{h}$ for all other $j\neq t$. We
conclude that $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$.

###  B.2 Checksum-preservation for parallel MGS

We consider the columns of the encoded matrix $\widetilde{A}$:

$\displaystyle\widetilde{A}=\begin{bmatrix}\bm{\widetilde{a}_{1}}&\bm{\widetilde{a}_{2}}&\cdots&\bm{\widetilde{a}_{n+d}}\end{bmatrix}\text{,
where $\bm{\widetilde{a}_{i}}=\begin{bmatrix}\bm{a_{i}}\\\
G_{v}\bm{a_{i}}\end{bmatrix}$}$

and the columns of the final $Q$-factor $Q^{(t)}$ at the $t^{th}$ iteration
for $t=1\rightarrow T$:

$\displaystyle
Q^{(t)}=\begin{bmatrix}\bm{q_{1}}^{(t)}&\bm{q_{2}}^{(t)}&\cdots&\bm{q_{n+d}}^{(t)}\end{bmatrix}$
$\displaystyle
Q_{1}^{(t)}=\begin{bmatrix}\bm{q_{11}}^{(t)}&\bm{q_{12}}^{(t)}&\cdots&\bm{q_{1(n+d)}}^{(t)}\end{bmatrix}$
$\displaystyle
Q_{2}^{(t)}=\begin{bmatrix}\bm{q_{21}}^{(t)}&\bm{q_{22}}^{(t)}&\cdots&\bm{q_{2(n+d)}}^{(t)}\end{bmatrix}\text{,
where $\bm{q_{i}}^{(t)}=\begin{bmatrix}\bm{q_{1i}}^{(t)}\\\
\bm{q_{2i}}^{(t)}\end{bmatrix}$.}$

Next, we present the BMGS algorithm on $\widetilde{A}$:

Input:
$\widetilde{A}=\begin{bmatrix}\widetilde{a}_{1}&\widetilde{a}_{2}&\cdots&\widetilde{a}_{n}\end{bmatrix}$
is $(m+c)\times n$ matrix

Output: $Q=\begin{bmatrix}q_{1}&q_{2}&\dots&q_{n}\end{bmatrix}$ is
$(m+c)\times n$ matrix , and $R=(r_{ij})$ is $n\times n$ matrix

Result: $\widetilde{A}=QR$, where $Q^{T}Q=I$

1 $Q=\widetilde{A}$

2 for  _$i=1$ to $n$ Step $b$_ do

3 $\bar{Q}=Q[:,i:i+b-1]$

4 MGS($\bar{Q},R[i:i+b-1,i:i+b-1]$)

5 $Q[:,i:i+b-1]=\bar{Q}$

6 if  _$i <n-b$_ then

7 $\bar{R}=\bar{Q}^{T}Q[:,i+b:n]$

8 $R[i:i+b-1,i+b:n]=\bar{R}$

9 $Q[:,i+b:n]-=\bar{Q}\bar{R}$

10 end if

11

12 end for

Algorithm 5 BMGS

The BMGS is essentially the reformulation of MGS beneficial for
parallelization. The parallel BMGS, called PBMGS, is obtained by parallelizing
certain steps in BMGS via MPI operations to adapt to the 2D cyclic
distribution, where step $b$ is set to the block size. As the computation and
thus numerical output at each iteration of PBMGS remain the same as that of
BMGS through the parallelization, it suffices to prove the checksum-
preservation for BMGS.

In the Algorithm 5, we use the temporary $(n+c)\times b$ matrix $\bar{Q}$. For
analysis, we consider its value $\bar{Q}^{(t)}$ at the end of each iteration
$t=1\rightarrow T$. $\bar{Q}$ is initialized to $\bar{Q}^{(0)}=\widetilde{A}$
on line 1 of Algorithm 5.

Next, we further consider the breakdowns of the intermediate matrices as
follows:

$\displaystyle\bar{Q}^{(t)}=\begin{bmatrix}\bar{Q}^{(t)}_{1}\\\
\bar{Q}^{(t)}_{2}\end{bmatrix}\text{, where $\bar{Q}_{1}:n\times b$ and
$\bar{Q}_{2}:c\times b$},$ $\displaystyle
R^{(t)}=\begin{bmatrix}R^{(t)}_{1}&R^{(t)}_{2}\end{bmatrix}=(r^{(t)}_{ij})=\begin{bmatrix}\bm{r_{1}}^{(t)}\\\
\vdots\\\ \bm{r_{n+c}}^{(t)}\end{bmatrix},$ $\displaystyle
R^{(t)}_{1}=\begin{bmatrix}\bm{\bar{r}_{11}}^{(t)}\\\ \vdots\\\
\bm{\bar{r}_{(n+c)1}}^{(t)}\end{bmatrix},R^{(t)}_{2}=\begin{bmatrix}\bm{\bar{r}_{12}}^{(t)}\\\
\vdots\\\ \bm{\bar{r}_{(n+c)2}}^{(t)}\end{bmatrix}\text{, where
$\bm{r_{i}}^{(t)}=\begin{bmatrix}\bm{r_{i1}}^{(t)}&\bm{r_{i2}}^{(t)}\end{bmatrix}$}.$

We proceed to prove by induction that for $t=0\rightarrow T$, we have
$Q_{2}^{(t)}=G_{v}Q_{1}^{(t)}$ and $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$.

Base case: For $t=0$:
As $Q$ is initialized to $Q^{(0)}=\widetilde{A}=\begin{bmatrix}A&AG_{h}\\\
G_{v}A&G_{v}AG_{h}\end{bmatrix}$, we have $Q_{2}^{(0)}=G_{v}Q_{1}^{(0)}$.

Inductive step: For the inductive step from $t-1$ to $t$ ($t\geq 1$):
At the $t^{th}$ iteration, we have $i=1+tb$ in the for-loop of BMGS. The
$Q$-factor is updated on line 5 for the $b$ columns $Q[:,i:i+b-1]$ and on line
9 for the $n-i-b+1$ columns $Q[:,i+b:n]$. The matrix $\bar{Q}$ is initialized
to $\bar{Q}=Q^{(t-1)}[:,i:i+b-1]$ on line 3 and then serves as the data input
for the MGS algorithm. As $Q_{2}^{(t-1)}=G_{v}Q_{1}^{(t-1)}$, we obtain that
$\bar{Q}_{2}=G_{v}\bar{Q}_{1}$ on line 3. By the checksum preservation
properties of MGS proven in Appendix B.1, the checksum relation
$\bar{Q}_{2}=G_{v}\bar{Q}_{1}$ of $\bar{Q}$ is maintained throughout and at
the end of the computation of MGS on line 4. Then, for the $b$ columns
$Q[:,i:i+b-1]$ updated on line 5, we have:

$\displaystyle Q^{(t)}[:,i:i+b-1]$
$\displaystyle=\bar{Q}=\begin{bmatrix}\bar{Q}_{1}\\\
\bar{Q}_{2}\end{bmatrix}=\begin{bmatrix}\bar{Q}_{1}\\\
G_{v}\bar{Q}_{1}\end{bmatrix}$
$\displaystyle=\begin{bmatrix}Q^{(t)}_{1}[:,i:i+b-1]\\\
Q^{(t)}_{2}[:,i:i+b-1]\end{bmatrix}.$

We thus obtain: $Q^{(t)}_{2}[:,i:i+b-1]=G_{v}Q^{(t)}_{1}[:,i:i+b-1]$. Now, for
the $n-i-b+1$ columns $Q[:,i+b:n]$ updated on line 9:

$\displaystyle Q^{(t)}[:,i+b:n]$
$\displaystyle=Q^{(t-1)}[:,i+b:n]-\bar{Q}\bar{R}$
$\displaystyle=\begin{bmatrix}Q^{(t-1)}_{1}[:,i+b:n]\\\
Q^{(t-1)}_{2}[:,i+b:n]\end{bmatrix}-\begin{bmatrix}\bar{Q}_{1}\\\
\bar{Q}_{2}\end{bmatrix}\bar{R}$
$\displaystyle=\begin{bmatrix}Q^{(t-1)}_{1}[:,i+b:n]\\\
GQ^{(t-1)}_{1}[:,i+b:n]\end{bmatrix}-\begin{bmatrix}\bar{Q}_{1}\\\
G\bar{Q}_{1}\end{bmatrix}\bar{R}$ (as $Q_{2}^{(t-1)}=GQ_{1}^{(t-1)}$ and
$\bar{Q}_{2}=G\bar{Q}_{1}$)
$\displaystyle=\begin{bmatrix}Q^{(t-1)}_{1}[:,i+b:n]-\bar{Q}_{1}\bar{R}\\\
G(Q^{(t-1)}_{1}[:,i+b:n]-\bar{Q}_{1}\bar{R})\end{bmatrix}$
$\displaystyle=\begin{bmatrix}Q^{(t)}_{1}[:,i+b:n]\\\
Q^{(t)}_{2}[:,i+b:n]\end{bmatrix}.$

We thus obtain: $Q^{(t)}_{2}[:,i+b:n]=GQ^{(t)}_{1}[:,i+b:n]$.
All in all, we have just proved that $q^{(t)}_{2j}=Gq^{(t)}_{1j}$ for
$j\in[i,n]$, where each of these columns are updated on either line 5 or line
9 of the BMGS algorithm. Other columns of $Q$ remain the same:
$q^{(t)}_{j}=q^{(t-1)}_{j}$, so $q^{(t)}_{2j}=Gq^{(t)}_{1j}$ for all other
$j<i=1+tb$. Consequently, we conclude that $Q_{2}^{(t)}=GQ_{1}^{(t)}$.

Now, it is left to show that $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$. We note that the
$R$-factor is updated on line 8 for the $b$ rows $R[i:i+b-1,:]$, and have the
following property.

######  Observation 3.

At the end of iteration $t$, we have:

(28) $\displaystyle R^{(t)}[i:i+b-1,:]=\bar{Q}^{T}\widetilde{A}$

######  Proof B.2.

The same reasoning as Observation 2’s proof can be adapted. Here we would like
to present an intuitive proof to show the validity of (28). From the
algorithm, we know that $R^{(t)}[i:i+b-1,:]$, and $Q^{(t)}[:,i:i+b-1]=\bar{Q}$
are last updated in this $t^{th}$ iteration and thus also the final values:

$\displaystyle R[i:i+b-1,:]=R^{(t)}[i:i+b-1,:],$ $\displaystyle
Q[:,i:i+b-1]=Q^{(t)}[:,i:i+b-1]=\bar{Q}.$

If the algorithm correctly computes the QR decomposition of $\widetilde{A}$,
at the end we would have $\widetilde{A}=QR$, implying:

$\displaystyle R=Q^{T}\widetilde{A},$ $\displaystyle
R[i:i+b-1,:]=Q[:,i:i+b-1]^{T}\widetilde{A},$ $\displaystyle
R^{(t)}[i:i+b-1,:]=Q^{(t)}[:,i:i+b-1]^{T}\widetilde{A},$ $\displaystyle
R^{(t)}[i:i+b-1,:]=\bar{Q}^{T}\widetilde{A}.$

Back to the main proof, we have:

$\displaystyle R^{(t)}[i:i+b-1,:]$
$\displaystyle=\bar{Q}^{T}\widetilde{A}=\bar{Q}^{T}\begin{bmatrix}A&AG_{h}\\\
G_{v}A&G_{v}AG_{h}\end{bmatrix}$
$\displaystyle=\begin{bmatrix}\bar{Q}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}&\bar{Q}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}G_{h}\end{bmatrix}$
$\displaystyle=\begin{bmatrix}R_{1}^{(t)}[i:i+b-1,:]&R_{2}^{(t)}[i:i+b-1,:]\end{bmatrix}$
$\displaystyle\implies$
$\displaystyle\begin{cases}R_{1}^{(t)}[i:i+b-1,:]=\bar{Q}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}\\\ R_{2}^{(t)}[i:i+b-1,:]=\bar{Q}^{T}\begin{bmatrix}A\\\
G_{v}A\end{bmatrix}G_{h}\end{cases}.$

We thus obtain: $\bm{r_{j2}}^{(t)}=\bm{r_{j1}}^{(t)}G_{h}$ for
$j\in[i,i+b-1]$. Other rows of $R$ remain the same:
$\bm{r_{j}}^{(t)}=\bm{r_{j}}^{(t-1)}$, so
$\bm{r_{j2}}^{(t)}=\bm{r_{j1}}^{(t)}G_{h}$ for all other $j\notin[i,i+b-1]$.
We conclude that $R_{2}^{(t)}=R_{1}^{(t)}G_{h}$.

###  B.3 Corollary 3.3: final form of the QR decomposition on $\widetilde{A}$

By Lemma 3.2, we know that at the end of the algorithm,
$Q_{2}^{(T)}=G_{v}Q_{1}^{(T)}$ and $R_{2}^{(T)}=R_{1}^{(T)}G_{h}$. As the QR
decomposition is retrieved as $\widetilde{A}=QR$, where
$Q=Q^{(T)}=\begin{bmatrix}Q_{1}^{(T)}\\\
Q_{2}^{(T)}\end{bmatrix}=\begin{bmatrix}Q_{1}\\\ Q_{2}\end{bmatrix}$ and
$R=R^{(T)}=\begin{bmatrix}R_{1}^{(T)}&R_{2}^{(T)}\end{bmatrix}=\begin{bmatrix}R_{1}&R_{2}\end{bmatrix}$,
equation (5) directly holds:

$\displaystyle\widetilde{A}=\begin{bmatrix}.Q_{1}\\\
Q_{2}\end{bmatrix}\begin{bmatrix}R_{1}&R_{2}\end{bmatrix}=\begin{bmatrix}Q_{1}\\\
G_{v}Q_{1}\end{bmatrix}\begin{bmatrix}R_{1}&R_{1}G_{h}\end{bmatrix}.$

##  Appendix C Example for Code Construction 2 and Theorem 5.4

######  Example 1 (Code construction for $f=1$ and $n=p_{r}=p_{c}=2$).

Consider the $2\times 2$ matrix $A=\begin{bmatrix}A_{11}&A_{12}\\\
A_{21}&A_{22}\end{bmatrix}$, each entry of which is held by a processor.

We now describe in details the steps of encoding, failure recovery and
decoding.

Encoding: For $R$-factor protection, we use
$G_{h}=\begin{bmatrix}1&1\end{bmatrix}$, which is MDS. For $Q$-factor
protection, we follow Construction 2 to initiate $\widetilde{V}=[v]$ and
$\widetilde{G}_{1}=[g_{1}]$ whereby $v\sim Unif(0,1)$ and
$g_{1}=-\frac{1}{2}v^{2}$ (since
$\widetilde{G}_{1}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}$). Given
$\widetilde{G}=\begin{bmatrix}g_{1}&v\end{bmatrix}$, we construct the vertical
checksum matrix in view of (18) as $G_{v}=\widetilde{G}\otimes
I_{n/p_{r}}=\begin{bmatrix}g_{1}&v\end{bmatrix}$. Now, we encode the input
matrix $A$ as:

(29) $A\xrightarrow{\text{Encode}}\widetilde{A}=\begin{bmatrix}A&AG_{h}\\\
G_{v}A&G_{v}AG_{h}\end{bmatrix}=\begin{bmatrix}A_{11}&A_{12}&\widetilde{A}_{13}\\\
A_{21}&A_{22}&\widetilde{A}_{23}\\\
\widetilde{A}_{31}&\widetilde{A}_{32}&\widetilde{A}_{33}\end{bmatrix}$

,
where, from our checksum design, the input matrix $A$ is protected via:

  * •

Horizontal checksums: $\widetilde{A}_{13}=A_{11}+A_{12}$ and
$\widetilde{A}_{23}=A_{21}+A_{22}$

  * •

Vertical checksums: $\widetilde{A}_{31}=g_{1}A_{11}+vA_{21}$ and
$\widetilde{A}_{32}=g_{1}A_{12}+vA_{22}$

Failure Recovery: At the end of the first iteration of the MGS algorithm, we
have the following view of the $Q$-factor and $R$-factor:

(30) $Q^{(1)}=\begin{bmatrix}Q_{11}&A_{12}&\widetilde{A}_{13}\\\
Q_{21}&A_{22}&\widetilde{A}_{23}\\\
Q_{31}&\widetilde{A}_{32}&\widetilde{A}_{33}\end{bmatrix},\quad
R^{(1)}=\begin{bmatrix}R_{11}&R_{12}&R_{13}\\\ 0&0&0\\\ 0&0&0\end{bmatrix}$

.
Now, from Lemma 3.2, we have $Q_{31}=g_{1}Q_{11}+vQ_{12}$ and
$R_{13}=R_{11}+R_{12}$, which, if complemented with the horizontal checksums
of $A$ as in the encoding phase, would provide the full protection for
$Q^{(1)}$ and $R^{(1)}$; e.g. assume that the processor $\Pi(1,1)$ holding
$Q_{11}$ and $R_{11}$ fails. Then we can recover
$Q_{11}=\frac{1}{g_{1}}(Q_{31}-vQ_{12})$ and $R_{11}=R_{13}-R_{12}$.

Decoding: At the end of the MGS algorithm, we have the following view of the
$Q$-factor and $R$-factor:

(31) $Q^{(3)}=\begin{bmatrix}Q_{11}&Q_{12}&Q_{13}\\\ Q_{21}&Q_{22}&Q_{23}\\\
Q_{31}&Q_{32}&Q_{33}\end{bmatrix},\quad
R^{(3)}=\begin{bmatrix}R_{11}&R_{12}&R_{13}\\\ 0&R_{22}&R_{23}\\\
0&0&0\end{bmatrix}$

from which we retrieve $A=Q_{1}R$ with $Q_{1}=\begin{bmatrix}Q_{11}&Q_{12}\\\
Q_{21}&Q_{22}\end{bmatrix}$ and $R=\begin{bmatrix}R_{11}&R_{12}\\\
0&R_{22}\end{bmatrix}$. Next, we construct $G_{0}=\begin{bmatrix}1+g_{1}&v\\\
v&-1\end{bmatrix}$ as in (10). Finally, to solve the full-rank square system
$A\bm{x}=\bm{b}$, we compute the orthogonal matrix $G_{0}Q_{1}$ and the
transformed vector $G_{0}\bm{b}$ in order to solve for $\bm{x}$ from (16).

##  Appendix D Omitted Proofs

###  D.1 Proof of Theorem 4.2

######  Proof D.1.

We first prove a lemma that provides a condition on $K$.

######  Lemma D.2.

Under the load-balanced in-node checksum distribution, to recover from any $r$
failed nodes out of $P$ nodes, the number of checksum blocks $K$ must satisfy:

(32) $a\lceil\frac{K}{P}\rceil+b\lfloor\frac{K}{P}\rfloor\geq\;f\frac{L}{P}$

for any $a+b=P-f$ and $a\leq\overline{a}$ and $b\leq\overline{b}$.

######  Proof D.3 (Proof of Lemma D.2).

Assume there exists $a^{\prime}$ and $b^{\prime}$ such that
$a^{\prime}+b^{\prime}=P-f$, $a^{\prime}\leq\overline{a}$,
$b^{\prime}\leq\overline{b}$ that satisfy
$a^{\prime}\lceil\frac{K}{P}\rceil+b^{\prime}\lfloor\frac{K}{P}\rfloor<f\frac{L}{P}$.
If $a^{\prime}$ nodes with $\lceil\frac{K}{P}\rceil$ checksum blocks and
$b^{\prime}$ nodes with $\lfloor\frac{K}{P}\rfloor$ checksum blocks are the
only surviving nodes, then we have fewer than $f\frac{L}{P}$ from the $P-f$
surviving nodes. However, we have $f\frac{L}{P}$ data blocks missing from the
$f$ failed nodes. Since we have fewer equations than the number of variables,
we cannot recover all the data blocks (Contradiction).

By Lemma 1, we have:

(33) $a\lceil\frac{K}{P}\rceil+b\lfloor\frac{K}{P}\rfloor\geq\;f\frac{L}{P}$

for any $a+b=P-f$ and $a\leq\overline{a}$ and $b\leq\overline{b}$. As this has
to be satisfied for the worst case $a$ and $b$, let us consider:

(34) $\min_{\begin{subarray}{c}a,b\geq 0,\;a+b=f,\\\
a\leq\overline{a},\;b\leq\overline{b}\end{subarray}}a\lceil\frac{K}{P}\rceil+b\lfloor\frac{K}{P}\rfloor=\begin{cases}(P-f)\lfloor\frac{K}{P}\rfloor&\text{
if }\overline{b}\geq P-f\\\
(P-f-\overline{b})\lceil\frac{K}{P}\rceil+\overline{b}\lfloor\frac{K}{P}\rfloor&\text{
otherwise. }\end{cases}$

Note that for any $K$,
$(P-f)\lceil\frac{K}{P}\rceil\geq(P-f)\lfloor\frac{K}{P}\rfloor$ and
$(P-f)\lceil\frac{K}{P}\rceil\geq(P-f-\overline{b})\lceil\frac{K}{P}\rceil+\overline{b}\lfloor\frac{K}{P}\rfloor$.
Hence,

(35) $(P-f)\lceil\frac{K}{P}\rceil\geq\min_{\begin{subarray}{c}a,b\geq
0,\;a+b=f,\\\
a\leq\overline{a},\;b\leq\overline{b}\end{subarray}}a\lceil\frac{K}{P}\rceil+b\lfloor\frac{K}{P}\rfloor\geq
f\frac{L}{P}.$

Because $\lceil\frac{K}{P}\rceil$ is an integer,

(36) $\lceil\frac{K}{P}\rceil\geq\lceil\frac{fL}{(P-f)P}\rceil.$

Now, let us find $K$ such that
$\lceil\frac{K}{P}\rceil=\lceil\frac{fL}{(P-f)P}\rceil$ that satisfy (33).
Then, $K$ be expressed as $K=P\lceil\frac{fL}{(P-f)P}\rceil-m$ where
$m=0,\ldots,P-1$. We will try to find the largest $m$, and we divide these
into three cases: i) $m\geq P-f$, ii) $0<m<P-f$, and iii) $m=0$.

Case i) $m\geq P-f$: Notice that $\overline{b}=m\geq P-f$, and thus we have:

(37) $\displaystyle\min_{\begin{subarray}{c}a,b\geq 0,\;a+b=f,\\\
a\leq\overline{a},\;b\leq\overline{b}\end{subarray}}a\lceil\frac{K}{P}\rceil+b\lfloor\frac{K}{P}\rfloor$
$\displaystyle=(P-f)\lfloor\frac{K}{P}\rfloor=(P-f)\left(\lceil\frac{fL}{(P-f)P}\rceil-1\right).$

When $\frac{fL}{(P-f)P}$ is an integer,
$\text{\eqref{eq:thm_lb_prf1}}=(P-f)\left(\frac{fL}{(P-f)P}-1\right)<\frac{fL}{P}$.
When $\frac{fL}{(P-f)P}$ is not an integer,
$\text{\eqref{eq:thm_lb_prf1}}=(P-f)\lfloor\frac{fL}{(P-f)P}\rfloor<\frac{fL}{P}$.
Thus, there does not exist $m\geq P-f$ that satisfies (33).

Case ii) $0<m<P-f$: As $\overline{b}=m<P-f$, we have:

$\displaystyle\min_{\begin{subarray}{c}a,b\geq 0,\;a+b=f,\\\
a\leq\overline{a},\;b\leq\overline{b}\end{subarray}}a\lceil\frac{K}{P}\rceil+b\lfloor\frac{K}{P}\rfloor$
$\displaystyle=(P-f-\overline{b})\lceil\frac{K}{P}\rceil+\overline{b}\lfloor\frac{K}{P}\rfloor=K-f\left(\lfloor\frac{K}{P}\rfloor+1\right)$
$\displaystyle=K-f\lceil\frac{K}{P}\rceil=K-f\lceil\frac{fL}{(P-f)P}\rceil\geq\frac{fL}{P}$

Hence, $K=f\lceil\frac{fL}{(P-f)P}\rceil+\frac{fL}{P}$ satisfy (33). We can
also check that:

$\displaystyle m$
$\displaystyle=P\lceil\frac{fL}{(P-f)P}\rceil-K=P\lceil\frac{fL}{(P-f)P}\rceil-f\lceil\frac{fL}{(P-f)P}\rceil-\frac{fL}{P}$
$\displaystyle=(P-f)\lceil\frac{fL}{(P-f)P}\rceil-\frac{fL}{P}<P-f.$

As we already found $m$ in Case ii), we do not have to consider Case iii).

###  D.2 Proof of Theorem 4.3

######  Proof D.4.

From any $P-f$ surviving nodes, we have $(P-f)\frac{L}{P}$ data blocks and
$f\frac{L}{P}$ checksum blocks as $K$ satisfies the condition in (32). Let
$p_{1},\ldots,p_{r}$ be the index of failed nodes. After subtracting the data
blocks from the survived nodes, the checksum equations given in (7) reduces to
equations in $f\frac{L}{P}$ data blocks as follows:

(38)
$C_{i}=\sum_{j\in\cup_{k=1}^{r}\mathcal{I}_{p_{k}}}\widetilde{G}_{i,j}D_{j}.$

We have $f\frac{L}{P}$ such equations and since any square submatrix of
$\widetilde{G}$ is full-rank, we can solve the set of linear equations to
retrieve $D_{j}$’s for $j\in\cup_{k=1}^{r}\mathcal{I}_{p_{k}}$.

###  D.3 Proof of Theorem 5.1

######  Proof D.5.

Recall that given $G=\begin{bmatrix}G_{1}&V\end{bmatrix}$, the $n\times n$
matrix $G_{0}$ is computed by:

$\displaystyle G_{0}=\begin{bmatrix}I_{c}+G_{1}&V\\\
V^{T}&-I_{n-c}\end{bmatrix}$

The following Observation 4 is used in our proof later:

######  Observation 4.

Under the restriction (11) that $G_{1}=-\frac{1}{2}VV^{T}$, the following
equation holds: $G_{0}^{T}G_{0}=I+G^{T}G$.

######  Proof D.6 (Proof of Observation 4).

From (11), we have $G_{1}^{T}=G_{1}=\frac{-1}{2}VV^{T}$, and obtain that:

$\displaystyle G_{0}^{T}G_{0}$ $\displaystyle=\begin{bmatrix}I+G_{1}&V\\\
V^{T}&-I\end{bmatrix}^{T}\begin{bmatrix}I+G_{1}&V\\\
V^{T}&-I\end{bmatrix}=\begin{bmatrix}(I+G_{1})^{2}+VV^{T}&(I+G_{1})V-V\\\
V^{T}(I+G_{1})-V^{T}&V^{T}V+I\end{bmatrix}$
$\displaystyle=\begin{bmatrix}I+G_{1}^{2}&G_{1}V\\\
V^{T}G_{1}&V^{T}V+I\end{bmatrix}=I+\begin{bmatrix}G_{1}^{2}&G_{1}V\\\
V^{T}G_{1}&V^{T}V\end{bmatrix}=I+G^{T}G.$

At the end of QR decomposition of encoded $\widetilde{A}$ in equation (8), the
left factor $\begin{bmatrix}Q_{1}\\\ GQ_{1}\end{bmatrix}$ is orthogonal :

$\displaystyle I$ $\displaystyle=\begin{bmatrix}Q_{1}\\\
GQ_{1}\end{bmatrix}^{T}\begin{bmatrix}Q_{1}\\\
GQ_{1}\end{bmatrix}=Q_{1}^{T}Q_{1}+Q_{1}^{T}G^{T}GQ_{1}$
$\displaystyle=Q_{1}^{T}(I+G^{T}G)Q_{1}=Q_{1}^{T}G_{0}^{T}G_{0}Q_{1}\text{ (by
Observation \ref{G0_relation})}$ $\displaystyle=(G_{0}Q_{1})^{T}(G_{0}Q_{1}).$

Therefore,$G_{0}Q_{1}$ is orthogonal. The following Observation 5 will be
useful for computing the determinant of $G_{0}$.

######  Observation 5.

Consider matrices $A,B,C,D$ of size $n\times n,n\times m,m\times n,m\times m$
respectively. If D is invertible then:

$det\begin{bmatrix}A&B\\\ C&D\end{bmatrix}=det(D)det(A-BD^{-1}C)$

We compute the determinant of $G_{0}$:

$\displaystyle det(G_{0})$ $\displaystyle=det\begin{bmatrix}I_{c}+G_{1}&V\\\
V^{T}&-I_{m-c}\end{bmatrix}$
$\displaystyle=det(-I_{m-c})det((I+G_{1})-V(-I)^{-1}V^{T})$
$\displaystyle=(-1)^{m-c}det(I+G_{1}+VV^{T})$

where the second line follows Observation 5. Under our setting
$G_{1}=\frac{-1}{2}VV^{T}$, the determinant of $G_{0}$ becomes:

$\displaystyle det(G_{0})=(-1)^{m-c}det(I+\frac{1}{2}VV^{T}).$

Now taking any non-zero column vector $\bm{z}$, we have
$\bm{z}^{T}(I+\frac{1}{2}VV^{T})\bm{z}=\bm{z}^{T}\bm{z}+\frac{1}{2}(V^{T}\bm{z})^{T}V^{T}\bm{z}\geq\bm{z}^{T}\bm{z}>0$.
This means that $I+\frac{1}{2}VV^{T}$ is positive definite, and
$det(I+\frac{1}{2}VV^{T})\neq 0$. Therefore, $det(G_{0})\neq 0$ and $G_{0}$ is
invertible.

###  D.4 Proof of Theorem 5.4

######  Proof D.7.

Let $n_{s}$ denote the dimension of $S$. We can divide this into two cases:
(i) When $S$ is a submatrix of $\widetilde{V}$, and (ii) when $t$ columns of
$S$ are from $\widetilde{G}_{1}$ and $(n_{s}-t)$ columns are from
$\widetilde{V}$ ($1\leq n_{s}\leq f$, $0\leq t\leq n_{s}$). For the first
case, any square submatrix of $\widetilde{V}$ is a random matrix where each
entry was chosen iid. Hence, is is full rank with probability 1. Thus, we
focus on the second case.

Let $\mathrm{Perm}(n_{s})$ be all permutations of $\\{1,\ldots,n_{s}\\}$, and
for permutation $\sigma$, $\sigma_{i}$ is the $i$-th element in the sequence.
For example when $n_{s}=3$, if $\sigma=[2\;\;3\;\;1]$, then
$\sigma_{1}=2,\sigma_{2}=3,\sigma_{3}=1$. The Leibniz formula for matrix
determinant is:

(39)
$\mathrm{det}(S)=\sum_{\sigma\in\mathrm{Perm}(n_{s})}\mathrm{sgn}(\sigma)\prod_{i=1}^{n_{s}}s_{i,\sigma_{i}}=\sum_{\sigma\in\mathrm{Perm}(n_{s})}\xi(\sigma),\qquad\xi(\sigma)\triangleq\mathrm{sgn}(\sigma)\prod_{i=1}^{n_{s}}s_{i,\sigma_{i}}.$

$\mathrm{sgn}(\sigma)$ is $+1$ or $-1$, but since this term is not important
for our argument, we refer to [48] for more explanation.

We now prove that $\mathrm{det}(S)$ is a non-zero polynomial in $v_{i,j}$’s.
Let us consider the term:

(40) $\xi(\sigma=[1\;\;\ldots\;\;n_{s}])=\prod_{i=1}^{n_{s}}s_{i,i}.$

For simplicity, let us assume that the $t$ columns from $\widetilde{G}_{1}$
are $\mathbf{g}_{1},\ldots\mathbf{g}_{t}$ and the $n_{s}-t$ columns of
$\widetilde{V}$ are $\mathbf{v}_{1},\ldots\mathbf{v}_{n_{s}-t}$. The argument
naturally extends to a general case. Now, $\prod_{i=1}^{n_{s}}s_{i,i}$ can be
expanded as:

(41) $\displaystyle\prod_{i=1}^{n_{s}}s_{i,i}$
$\displaystyle=\prod_{i=1}^{t}\widetilde{g}_{i,i}\prod_{i=1}^{n_{s}-t}v_{i,t+i}$
(42)
$\displaystyle=\left(-\frac{1}{2}\right)^{t}\cdot\left(\prod_{i=1}^{t}\sum_{j=1}^{p_{r}}v_{i,j}^{2}\right)\cdot\prod_{i=1}^{n_{s}-t}v_{i,t+i}.$

In (42), let us focus on one term:
$(v_{1,p_{r}})^{2}(v_{2,p_{r}-1})^{2}\cdots(v_{t,p_{r}-t+1})^{2}\cdot\prod_{l=1}^{n_{s}-t}v_{l,t+l}\triangleq\psi$.
We show that this term does not appear anywhere else in $\xi(\sigma^{\prime})$
when $\sigma^{\prime}\neq[1\;\;\ldots\;\;n_{s}]$, and thus this term cannot be
canceled out. We prove this by contradiction.

Assume there exists a $\sigma^{\prime}$ such that
$\psi\in\xi(\sigma^{\prime})$. First, let us assume that there exists $1\leq
i\leq t$ such that $\sigma^{\prime}_{i}=j\neq i$. We get the $v_{i,p_{r}-i+1}$
term only from $\widetilde{g}_{i,j}$ and $\widetilde{g}_{k,i}$ ($k\neq i$) as:

(43) $\displaystyle\widetilde{g}_{i,j}$
$\displaystyle=v_{i,1}v_{j,1}+\cdots+v_{i,p_{r}-i+1}v_{j,p_{r}-i+1}+\cdots
v_{i,p_{r}}v_{j,p_{r}}$ (44) $\displaystyle\widetilde{g}_{k,i}$
$\displaystyle=v_{k,1}v_{i,1}+\cdots+v_{k,p_{r}-i+1}v_{i,p_{r}-i+1}+\cdots
v_{k,p_{r}}v_{i,p_{r}}.$

Then, the only term that has $(v_{i,p_{r}-i+1})^{2}$ is the product of the
following terms $v_{i,p_{r}-i+1}v_{j,p_{r}-i+1}$ and
$v_{k,p_{r}-i+1}v_{i,p_{r}-i+1}$, i.e.

(45) $v_{i,p_{r}-i+1}^{2}v_{j,p_{r}-i+1}v_{k,p_{r}-i+1}.$

Notice that this term has $v_{j,p_{r}-i+1}v_{k,p_{r}-i+1}$. We next proceed to
prove that $v_{j,p_{r}-i+1}$ and similarly $v_{k,p_{r}-i+1}$ are not present
in $\psi$. First, since $j\neq i$, the term $v_{j,p_{r}-i+1}$ does not appear
in $\prod_{l=1}^{t}(v_{l,p_{r}-l+1})^{2}$ of $\psi$. It is left to show that
the term $v_{j,p_{r}-i+1}$ also does not appear in
$\prod_{l=1}^{n_{s}-t}v_{l,t+l}$. For the sake of contradiction, suppose that
$v_{j,p_{r}-i+1}$ is one of those $v_{l,t+l}$’s for $l=1\to n_{s}-t$. Then,
this is equivalent to having $j\in[1,n_{s}-t]$ and $p_{r}-i+1=t+j$, which
would imply:

$\displaystyle n_{s}-t$ $\displaystyle\geq j=p_{r}-i+1-t$ $\displaystyle
n_{s}$ $\displaystyle\geq p_{r}-i+1\geq p_{r}-n_{s}+1$
$\displaystyle\therefore 2n_{s}$ $\displaystyle\geq p_{r}+1>p_{r}.$

However, the last line is a contradiction since $2n_{s}\leq 2f\leq p_{r}$.
Therefore, we have proven that $v_{j,p_{r}-i+1}$ is not present in $\psi$.
Similarly, $v_{k,p_{r}-i+1}$ and thus (45) are not present in $\psi$. Hence,
$\psi\notin\xi(\sigma)$. Now, assume that there exists $t+1\leq i\leq n_{s}$
such that $\sigma^{\prime}_{i}=j\neq i$. Then, $\xi(\sigma^{\prime})$ always
has the term $v_{t+i,j}$ ($j\neq i$), which was not present in $\psi$. Hence,
$\psi\notin\xi(\sigma)$. This proves that $\psi$ appears only in
$\xi(\sigma=[1\;\;\cdots\;\;n_{s}])$ and hence this term has a a non-zero
coefficient, $(-1/2)^{t}$. As
$\mathrm{Pr}(\\{v_{i,j}\\}_{\begin{subarray}{c}i=1,\ldots,r,\\\
j=1,\ldots,p_{r}\end{subarray}}:p(\\{v_{i,j}\\})=0)=0$ for a non-zero
polynomial $p$, this completes the proof.

###  D.5 Proof of Theorem 6.1

####  D.5.1 Supplementary Lemmas

First, we present the following Lemmas that characterize the analytical
overhead of coded QR decomposition.

######  Lemma D.8.

Under the out-of-node checksum storage, the coding strategy in Construction 2
can achieve the following coding overheads:

(46) $\displaystyle T_{\text{enc}}=f\cdot
T_{\text{reduce}}(p_{r}+1,\frac{n^{2}}{P})+\gamma\frac{fn^{2}}{P}+\gamma
f(f+1)(p_{r}-f),$ $\displaystyle T_{\text{post}}=f\cdot
T_{\text{broadcast}}(p_{r},\frac{n(n+p_{c})}{P})+f\cdot
T_{\text{reduce}}(p_{r}-f+1,\frac{n(n+p_{c})}{P})$ (47)
$\displaystyle\quad\quad\quad+\gamma(2f-1)\frac{n(n+p_{c})}{P},$ (48)
$\displaystyle T_{\text{comp}}\leq\frac{f}{p_{r}}T_{\text{QR}}.$

######  Proof D.9.

Now we analyze the overhead of encoding, recovery (decoding), and post-
processing.

1) Encoding overhead: For encoding, we instantiate
$\widetilde{G}=(\widetilde{g}_{i,j})$ following Construction 2, and compute
$C_{i,j}=\sum_{t=0}^{p_{r}-1}\widetilde{g}_{i,t}A_{t,j}$ to be sent to the
processor $\Pi_{C}(i,j)$. The former takes $\gamma f(p_{r}-f)$ computational
cost for generating the $f\times(p_{r}-f)$ random matrix
$\widetilde{V}=[\widetilde{v}_{1},...,\widetilde{v}_{p_{r}-f}]$ (where
$\widetilde{v}_{i}\in\mathbb{R}^{f\times 1}$) on all processors555In this
step, due to small computational burden, we can let all processors to generate
in parallel the same random matrix $\widetilde{V}$ via the same random seed.
and $\gamma f^{2}(p_{r}-f)$ cost to compute
$\widetilde{G}_{1}=-\frac{1}{2}\widetilde{V}\widetilde{V}^{T}$. The later step
computing $C_{i,j}$ can be done by sequentially performing linear-combination
reduce in $f$ iterations. Hence, the total encoding overhead is:

$\displaystyle T_{\text{enc}}$ $\displaystyle=\gamma f(f+1)(p_{r}-f)+f\cdot
T_{\text{lin-comb}}(p_{r},\frac{n^{2}}{P})$ $\displaystyle=\gamma
f(f+1)(p_{r}-f)+f\cdot\left[T_{\text{reduce}}(p_{r}+1,\frac{n^{2}}{P})+\gamma\frac{n^{2}}{P}\right]$
$\displaystyle=f\cdot
T_{\text{reduce}}(p_{r}+1,\frac{n^{2}}{P})+\gamma\frac{fn^{2}}{P}+\gamma
f(f+1)(p_{r}-f).$

2) Post-processing overhead: The overhead for post-processing is attributed to
the cost for performing matrix multiplications $(G_{0}Q_{1})$ and
$(G_{0}\bm{b})$ (refer to Section 5.1). We hereby present an efficient way for
such computation. Since
$G_{0}Q_{1}=[G_{0}Q_{1}[:,0],Q_{1}[:,1],...,Q_{1}[:,p_{r}-1]]$, we compute all
the $G_{0}Q_{1}[:,j]$’s in parallel. Furthermore, this can be combined with
the computation of $(G_{0}\bm{b})$ for efficiency Specifically, to formally
describe our computational pipeline, we first consider the
$n\times(\frac{n}{p_{c}}+1)$ matrix
$\overline{Z}^{j}=\begin{bmatrix}Q_{1}[:,j]&\bm{b}\end{bmatrix}=\begin{bmatrix}\overline{Z}^{j}_{0}\\\
\overline{Z}^{j}_{1}\\\ \vdots\\\ \overline{Z}^{j}_{p_{r}-1}\end{bmatrix}$,
where each sub-block $\overline{Z}^{j}_{i}$ is of dimension
$\frac{n}{p_{r}}\times(\frac{n}{p_{c}}+1)$ . Then, we aim to compute
$G_{0}\overline{Z}^{0},G_{0}\overline{Z}^{1},...,G_{0}\overline{Z}^{p_{r}-1}$
in parallel across the $p_{r}$ block-rows. Furthermore, thanks to
parallelisim, the post-processing overhead is the cost for performing matrix
multiplication $G_{0}\overline{Z}^{j}$, which is further expressed as the
followings to utilize the sparsity of $G_{0}$:

$\displaystyle G_{0}\overline{Z}^{j}$
$\displaystyle=\begin{bmatrix}I_{\frac{fn}{p_{r}}}+G_{1}&V\\\
V^{T}&-I_{n-\frac{fn}{p_{r}}}\end{bmatrix}\overline{Z}^{j}$
$\displaystyle=\begin{bmatrix}H\otimes
I_{\frac{n}{p_{r}}}&\widetilde{V}\otimes I_{\frac{n}{p_{r}}}\\\
\widetilde{V}^{T}\otimes I_{\frac{n}{p_{r}}}&-I_{p_{r}-f}\otimes
I_{\frac{n}{p_{r}}}\end{bmatrix}\begin{bmatrix}\overline{Z}^{j}_{0}\\\
\overline{Z}^{j}_{1}\\\ \vdots\\\
\overline{Z}^{j}_{p_{r}-1}\end{bmatrix}\text{ where
$H=(I_{f}+\widetilde{G}_{1})$}$
$\displaystyle=\begin{bmatrix}\sum_{i=0}^{f-1}h_{0,i}\cdot\overline{Z}^{j}_{i}+\sum_{t=f}^{p_{r}-1}\widetilde{v}_{0,t}\cdot\overline{Z}_{t}\\\
\sum_{i=0}^{f-1}h_{1,i}\cdot\overline{Z}^{j}_{i}+\sum_{t=f}^{p_{r}-1}\widetilde{v}_{1,t}\cdot\overline{Z}^{j}_{t}\\\
\vdots\\\
\sum_{i=0}^{f-1}h_{f-1,i}\cdot\overline{Z}^{j}_{i}+\sum_{t=f}^{p_{r}-1}\widetilde{v}_{f-1,t}\cdot\overline{Z}^{j}_{t}\\\
\sum_{i=0}^{f-1}\widetilde{v}_{i,0}\cdot\overline{Z}^{j}_{i}-\overline{Z}^{j}_{f}\\\
\sum_{i=0}^{f-1}\widetilde{v}_{i,1}\cdot\overline{Z}^{j}_{i}-\overline{Z}^{j}_{f+1}\\\
\vdots\\\
\sum_{i=0}^{f-1}\widetilde{v}_{i,p_{r}-f-1}\cdot\overline{Z}^{j}_{i}-\overline{Z}^{j}_{p_{r}-1}\\\
\end{bmatrix}\triangleq\begin{bmatrix}W_{0}\\\ W_{1}\\\ \vdots\\\
W_{p_{r}-1}\end{bmatrix}.$

Notice that all $W_{i}$’s for $i=0,\ldots,p_{r}-1$ have some linear
combination of $\overline{Z}^{j}_{t}$ ($t=0,\ldots f-1$). To compute this, we
first broadcast $\overline{Z}^{j}_{t}$’s ($t=0,\ldots f-1$) to all the nodes.
This requires performing $f$ operations of MPI_broadcast($p_{r}$,
$\frac{n}{p_{r}}\times(\frac{n}{p_{c}}+1)$). Next, we describe the steps of
computing these $W_{i}$’s for $i=0,\ldots,p_{r}-1$ in details for algorithmic
implementation in practice. Before we begin the series of broadcast
operations, we first initialize $W_{i}$ at each node as:

$\displaystyle W_{i}=h_{i,i}\overline{Z}^{j}_{i}\;\;\text{ for
}i<f,\;\;W_{i}=-\overline{Z}^{j}_{i}\;\;\text{ for }i\geq f.$

After the $l$-th broadcast of the matrix $\overline{Z}^{j}_{l}$
($l=0,\ldots,f-1$), nodes will update its $W_{i}$ as follows:

$\displaystyle W_{i}=W_{i}+h_{i,l}\overline{Z}^{j}_{l}\;\;\text{ for
}i<f(i\neq
l),\;\;W_{i}=W_{i}+\widetilde{v}_{l,i}\overline{Z}^{j}_{l}\;\;\text{ for
}i\geq f.$

The initialization and updating $W_{i}$ after each iteration requires a total
of $(2f-1)\frac{n}{p_{r}}\times(\frac{n}{p_{c}}+1)=(2f-1)\frac{n(n+p_{c})}{P}$
flops. After this, the computation for $W_{f},\ldots,W_{p_{r}-1}$ is complete,
and we only need more computation for $W_{0},\ldots W_{f-1}$ to obtain
$\sum_{t=f}^{p_{r}-1}\widetilde{v}_{i,t}\cdot\overline{Z}^{j}_{t}$ for
$i=0,\ldots,f-1$. This can be achieved by $f$ operations of linear-combination
reduce from the last $p_{r}-f$ nodes to one of the first $f$ nodes, which
incurs $fT_{\text{lin-comb}}(p_{r}-f,\frac{n(n+p_{c})}{P})$. Hence, the total
cost of post-processing is given by:

$\displaystyle T_{\text{post}}=$ $\displaystyle f\cdot
T_{\text{broadcast}}(p_{r},\frac{n(n+p_{c})}{P})+f\cdot
T_{\text{reduce}}(p_{r}-f+1,\frac{n(n+p_{c})}{P})$
$\displaystyle+\gamma(2f-1)\frac{n(n+p_{c})}{P}.$

3) Increased computation cost for QR decomposition of the encoded matrix:
Using Observation 1, we have $T_{\text{comp}}\leq\frac{c}{n}T_{\text{QR}}$
where $T_{\text{comp}}$ is the overhead for running QR decomposition on the
larger $(m+c)\times n$ encoded matrix $\widetilde{A}$ instead of $m\times n$
matrix $A$. In our code construction, we use $c=\frac{fn}{p_{r}}$ checksums
and thus obtain that:

(49) $\displaystyle T_{\text{comp}}$
$\displaystyle\leq\frac{fn/p_{r}}{n}T_{\text{QR}}=\frac{f}{p_{r}}T_{\text{QR}}.$

######  Lemma D.10.

Under the out-of-node checksum storage, the coding strategy in Construction 2
can achieve the following recovery overhead:

(50) $\displaystyle T_{\text{recov}}=f\cdot
T_{reduce}(p_{r}+1,\frac{n^{2}}{P})+\gamma(f_{1}^{2}+\frac{2}{3}f_{1}^{3}+\frac{fn^{2}}{P}).$

######  Proof D.11.

Before we begin our analysis we define a commonly-occurring operation,
“linear-combination reduce”. The goal of this operation is to compute the
linear combination of the data blocks at $p$ nodes and then send it to a
separate destination node. I.e., we compute
$\alpha_{1}D_{1}+\cdots\alpha_{p}D_{p}$ for some
$\alpha_{i}\in\mathbb{R}\;\;(i=1,\ldots,p)$ distributedly, and send it to the
$(p+1)$-th node. This can be achieved through MPI_reduce($p+1$, $w$) operation
by setting the data block at the destination node to be $\bm{0}$. The overhead
of this operation when the size of the data blocks is $w$ is: $T_{\text{lin-
com}}(p,w)=T_{\text{reduce}}(p+1,w)+\gamma w$ where $\gamma w$ comes from
computing $\alpha_{i}\cdot D_{i}$ for $i=1,\ldots,p$ in parallel.

Recovery requires first decoding the lost data blocks from the failed
systematic nodes and then recomputing the lost checksum blocks from the failed
checksum nodes. Let $f_{1}$ and $f_{2}$ be the number of failed systematic
nodes and failed checksum nodes. In the proof we only consider the case when
$f=f_{1}+f_{2}$ because the overhead when the number of failures is smaller
than $f$ is upper bounded by the worst-case recovery cost of having $f$
failures.

We partition the index set
$\\{0,1,...,p_{r}-1\\}=\mathcal{A}_{\text{fail}}\cup\mathcal{A}_{\text{succ}}$
, where $\mathcal{A}_{\text{fail}}$ denotes the indices of failed nodes among
the systematic nodes of the $j$-th column, i.e.,

$\mathcal{A}_{\text{fail}}=\\{t|\Pi(t,j)\text{ is a failed systematic
processor}\\}$

and $\mathcal{A}_{\text{succ}}$ denotes the indices of successful nodes among
the systematic nodes of the $j$-th column. Similarly, we partition
$\\{0,1,...,f-1\\}=\mathcal{C}_{\text{fail}}\cup\mathcal{C}_{\text{succ}}$ to
denote the failed and successful node indices within the checksum nodes on the
$j$-th column. Note that $|\mathcal{A}_{\text{fail}}|=f_{1}$ and
$|\mathcal{C}_{\text{fail}}|=f_{2}$.

First, (19) can be rewritten as:

(51) $\displaystyle
C_{i,j}-\sum_{t\in\mathcal{A}_{succ}}\widetilde{g}_{i,t}A_{t,j}$
$\displaystyle=\sum_{t\in\mathcal{A}_{\text{fail}}}\widetilde{g}_{i,t}A_{t,j}\forall
i\in\mathcal{C}_{\text{succ}}$

Let
$\widehat{G}=[\widetilde{g}_{i,t}]_{i\in\mathcal{C}_{\text{succ}},t\in\mathcal{A}_{\text{fail}}}$
bethe $f_{1}\times f_{1}$ matrix. To decode $A_{t,j}$ for
$t\in\mathcal{A}_{\text{fail}}$, we need to invert $\widehat{G}$. Note that
$\widehat{G}$ is always invertible as any square submatrix of $\widetilde{G}$
is invertible. Denoting the $(i,j)$-th entry of $\widehat{G}^{-1}$ as
$\widehat{g}^{-1}_{i,j}$, $A_{t,j}$ ($t\in\mathcal{A}_{\text{fail}}$) can be
recovered from:

$\displaystyle A_{t,j}$
$\displaystyle=\sum_{i\in\mathcal{C}_{\text{succ}}}\widehat{g}^{-1}_{t,i}(C_{i,j}-\sum_{l\in\mathcal{A}_{succ}}\widetilde{g}_{t,l}A_{l,j})$
(52)
$\displaystyle=\sum_{i\in\mathcal{C}_{\text{succ}}}\widehat{g}^{-1}_{t,i}\cdot
C_{i,j}-\sum_{l\in\mathcal{A}_{succ}}\left(\sum_{i\in\mathcal{C}_{\text{succ}}}\widehat{g}^{-1}_{t,i}\cdot\widetilde{g}_{t,l}\right)\cdot
A_{l,j},$

which is a linear combination of $C_{i,j}$ ($i\in\mathcal{C}_{\text{succ}}$)
and $A_{l,j}$ ($l\in\mathcal{A}_{succ}$). Computing the coefficients in (52)
requires inverting $\widehat{G}$, which costs $\gamma\frac{2}{3}f_{1}^{3}$ and
computing
$\sum_{i\in\mathcal{C}_{\text{succ}}}\widehat{g}^{-1}_{t,i}\cdot\widetilde{g}_{t,l}$,
which costs $\gamma f_{1}^{2}$. Once we have the coefficients, we can perform
$f_{1}$ iterations of linear-combination reduce. The total cost of recovering
the systematic data blocks is:

(53) $f_{1}T_{\text{lin-
comb}}(p_{r},\frac{n^{2}}{P})+\gamma(\frac{2}{3}f_{1}^{3}+f_{1}^{2})=f_{1}\cdot
T_{\text{reduce}}(p_{r}+1,\frac{n^{2}}{P})+\gamma(\frac{2}{3}f_{1}^{3}+f_{1}^{2}+\frac{f_{1}n^{2}}{P}).$

It is left to recover $f_{2}$ failed checksum blocks, so that the QR
decomposition is still fault-tolerant in the later iterations. This
computation is the same as encoding the checksums (with the exception that we
only have to compute a subset of checksums now). We thus defer readers to the
analysis encoding overhead for more details, and only state the total cost of
this step, which is This requires: $f_{2}T_{\text{lin-
comb}}(p_{r},\frac{n^{2}}{P})=f_{2}\cdot
T_{reduce}(p_{r}+1,\frac{n^{2}}{P})+\gamma\frac{f_{2}n^{2}}{P}$ .

In total, the recovery overhead is thus:

$\displaystyle T_{\text{recov}}=f\cdot
T_{reduce}(p_{r}+1,\frac{n^{2}}{P})+\gamma(f_{1}^{2}+\frac{2}{3}f_{1}^{3}+\frac{fn^{2}}{P}).$

####  D.5.2 Proof of the main Theorem

######  Proof D.12.

We have the following bound of $T_{QR}$ from (24):

(54) $\displaystyle T_{QR}\geq 2\alpha n\log
p_{r}+\beta\frac{1}{2}\frac{p_{c}n(n+1)}{P}+\gamma\frac{n^{2}(n+1)}{P}.$

From Lemma D.8, we have:

$\displaystyle T_{\text{enc}}$ $\displaystyle=f\cdot
T_{\text{reduce}}(p_{r}+1,\frac{n^{2}}{P})+\gamma\frac{fn^{2}}{P}+\gamma
f(f+1)(p_{r}-f)$ (55) $\displaystyle\leq 2\alpha
f\log(p_{r}+1)+2\beta\frac{fn(n+1)}{P}+\gamma\frac{2fn^{2}}{P}+\gamma
f(f+1)(p_{r}-f)$ $\displaystyle=(\frac{4f}{p_{c}}\cdot 2\alpha n\log
p_{r})\frac{p_{c}}{n}\frac{\log(p_{r}+1)}{\log(p_{r}^{4})}+(\frac{4f}{p_{c}}\cdot\beta\frac{1}{2}\frac{p_{c}n(n+1)}{P})$
$\displaystyle\quad+(\frac{4f}{p_{c}}\cdot\gamma\frac{n^{2}(n+1)}{P})\frac{p_{c}}{2(n+1)}+(\frac{f}{p_{c}}\cdot\gamma\frac{n^{2}(n+1)}{P})\frac{(f+1)p_{c}(p_{r}-f)}{n^{2}(n+1)}$
(56) $\displaystyle\leq\frac{5f}{p_{c}}T_{QR},$

where for the last inequality we use $f\leq\frac{1}{2}\min\\{p_{r},p_{c}\\}$
and $P=p_{r}p_{c}\leq n$. From Lemma D.8, we have:

$\displaystyle T_{\text{post}}$ $\displaystyle\leq f\cdot
T_{\text{broadcast}}(p_{r},\frac{2n^{2}}{P})+f\cdot
T_{\text{reduce}}(p_{r},\frac{2n^{2}}{P})+\gamma\frac{4fn^{2}}{P}$
$\displaystyle\leq
f(\sqrt{\alpha\log(p_{r})}+\sqrt{\beta\frac{2n^{2}}{P}})^{2}+2\alpha
f\log(p_{r})+2\beta\cdot\frac{2fn^{2}}{P}+\gamma\frac{6fn^{2}}{P}$ (57)
$\displaystyle\leq 4\alpha
f\log(p_{r})+8\beta\cdot\frac{fn^{2}}{P}+\gamma\frac{6fn^{2}}{P}$
$\displaystyle\leq(\frac{16f}{p_{c}}\cdot 2\alpha
n\log(p_{r}))\cdot\frac{p_{c}}{8n}+\frac{16f}{p_{c}}\cdot\beta\frac{1}{2}\frac{p_{c}n(n+1)}{P}$
$\displaystyle\quad+(\frac{16f}{p_{c}}\cdot\gamma\frac{n^{2}(n+1)}{P})\frac{3p_{c}}{8(n+1)})$
(58) $\displaystyle\leq\frac{16f}{p_{c}}T_{QR}\text{ (using $p_{c}\leq n$)}.$

Using (A.2), we have:

(59) $\displaystyle T^{\alpha}_{\text{comp}}$ $\displaystyle=O(1).$

Combining the bounds (56), (58) and (49), we obtain that:

$\displaystyle T_{\text{coding}}$
$\displaystyle=T_{\text{enc}}+T_{\text{post}}+T_{\text{comp}}\leq(\frac{6f}{p_{c}}+\frac{16f}{p_{c}}+\frac{f}{p_{r}})T_{QR}=(\frac{22f}{p_{c}}+\frac{f}{p_{r}})T_{QR},$
$\displaystyle\therefore T_{\text{coding}}$
$\displaystyle=O\left(\frac{f}{p_{r}}+\frac{f}{p_{c}}\right)T_{\text{QR}}.$

From Lemma D.8, we have:

$\displaystyle T_{\text{recov}}$ $\displaystyle\leq
T_{\text{reduce}}(p_{r},\frac{fn^{2}}{P})+\gamma[\frac{fn^{2}}{P}+\frac{2}{3}\frac{fn^{2}}{P}+\frac{fn^{2}}{P}]$
(using $f_{1}^{2}\leq f_{1}^{3}\leq f^{3}\leq
fp_{r}p_{c}=fP\leq\frac{fn^{2}}{P}$) $\displaystyle\leq
3T_{\text{enc}}\leq\frac{12f}{p_{c}}T_{\text{QR}},$ $\displaystyle\therefore
T_{\text{recov}}$ $\displaystyle=O\left(\frac{f}{p_{c}}\right)T_{\text{QR}}$

From (55), (57) and (59) , we have:

$\displaystyle T^{\alpha}_{\text{enc}}$ $\displaystyle\leq 2\alpha
f\log(p_{r}+1)+2\alpha\log(p_{r}+f)=O(\alpha f\log(p_{r})),$ $\displaystyle
T^{\alpha}_{\text{post}}$ $\displaystyle\leq 4\alpha f\log(p_{r})=O(\alpha
f\log(p_{r})),$ $\displaystyle T^{\alpha}_{\text{coding}}$
$\displaystyle=T^{\alpha}_{\text{enc}}+T^{\alpha}_{\text{post}}+T^{\alpha}_{\text{comp}}=O(\alpha
f\log(p_{r}))$ $\displaystyle=O\left(\frac{f}{n}\right)2\alpha
n\log(p_{r})\leq O\left(\frac{f}{n}\right)T^{\alpha}_{\text{QR}}\text{ (by
\eqref{QR_bound1})},$ $\displaystyle T^{\alpha}_{\text{coding}}$
$\displaystyle=O\left(\frac{f}{n}\right)T^{\alpha}_{\text{QR}}.$