  We discuss, study, and compare experimentally three methods for solving the
system of algebraic equations $\mathbb{A}^\alpha \bf{u}=\bf{f}$, $0< \alpha
<1$, where $\mathbb{A}$ is a symmetric and positive definite matrix obtained
from finite difference or finite element approximations of second order
elliptic problems in $\mathbb{R}^d$, $d=1,2,3$. The first method, introduced by
Harizanov et.al, based on the best uniform rational approximation (BURA)
$r_\alpha(t)$ of $t^{1-\alpha}$ for $0 \le t \le 1$, is used to get the
rational approximation of $t^{-\alpha}$ in the form $t^{-1}r_\alpha(t)$. Here
we develop another method, denoted by R-BURA, that is based on the best
rational approximation $r_{1-\alpha}(t)$ of $t^\alpha$ on the interval $[0,1]$
and approximates $t^{-\alpha}$ via $r^{-1}_{1-\alpha}(t)$. The third method,
introduced and studied by Bonito and Pasciak, is based on an exponentially
convergent quadrature scheme for the Dundord-Taylor integral representation of
the fractional powers of elliptic operators. All three methods reduce the
solution of the system $\mathbb{A}^\alpha \bf{u}=\bf{f}$ to solving a number of
equations of the type $(\mathbb{A} +c\mathbb{I})\bf{u}= \bf{f}$, $c \ge 0$.
Comprehensive numerical experiments on model problems with $\mathbb A$ obtained
by approximation of elliptic equations in one and two spatial dimensions are
used to compare the efficiency of these three algorithms depending on the
fractional power $\alpha$. The presented results prove the concept of the new
R-BURA method, which performs well for $\alpha$ close to $1$ in contrast to
BURA, which performs well for $\alpha $ close to $0$. As a result, we show
theoretically and experimentally, that they have mutually complementary
advantages.
