数理解析研究所講究録
第 666 巻 1988 年 115-124

Area-Time Efficient Evaluation of Elementary Functions
Yasuo OKABE and Shuzo YAJIMA
Faculty of Engineering, Kyoto University

ABSTRACTThis paper describes area-time efficient VLSI algorithms for evaluating elementaryfunctions. For square rooting, a VLSI implementation of A -optimal, i.e., A $T^{2}=0(n^{2})$ ,is the computation time in theis the chip area andcircuits are presented, where. For logarithms and exponentials, anrange $[\Omega((\log n)^{1+c}, O(\sqrt{n})]$ for arbitrary
upper bound $AT^{2}=0(n^{2}(\log n)^{2})$ is given for any $T\in[\Omega((\log n)^{2+\epsilon}, 0(\sqrt{n}\log n)]$. VLSI
circuits for these functions witlt minimum computation time $0(1ogn)$ and almost-preformance, $i.e.,$ A $T^{2}=O(n^{2+c})$ , are also exhibited. These are achieved by
optimalusing an extension of the Beame-Cook-Hoover Method which we have proposed before.
$T^{2}$$T$$A$$\epsilon>0$$A$$T^{2}$

1. Introduction

Much research has long been conducted on arithmetic operations, such as addition,
multiplication, division, square rooting and evaluation of elementary functions. Recent

advances in large-scale integration technology of circuits have especially motivated theresearch on hardware algorithms for arithmetic operations, suitable for VLSI implemen-tations.
In VLSI circuits, the chip area is a more reasonable cost measure than the numberof gates. Considering this, VLSI models are proposed as more practical computation
models $[1][2]$ , and much attention has been paid about area-time tradeoffs for many

basic operations such as arithmetic operations.
For multiplication,$A$$T^{2}$-optimal, $i.e.,$ A $T^{2}=O(n^{2})$ , multipliers are known for allcomputation times in the range
Preparata exhibited an

range

$A$

$T^{2}$[3].

$[\Omega(\log n), O(\sqrt{n})]$For division, Mehlhorn and-optimal design of dividers with computation times in the$[\Omega((\log n)^{1+\epsilon}),O(\sqrt{n})]$for arbitrary constant

Hoover division technique. Mehlhorn also exhibited$e>0$$A$ $T^{2}$[4], utilizing Beame-Cook--optimal square rooting circuits

for all computation time in the range $[\Omega(.(\log n)^{3}),O(\sqrt{n})][5]$ .
In this paper, we present several fast and area-time efficient algorithms for evaluat-

ing elementary functions, and give some new upper bounds of thethese functions. First we describetion times in the range$A$$T^{2}$$A$$T^{2}$-complexity of-optimal square rooting circuits for all computa-$[\Omega((\log n)^{1+\epsilon}),O(\sqrt{n})]$-1-

for arbitrary $e>0$ . Next we exhibit

$O(n^{2+\epsilon})$for arbitraryvalue for$x$$\epsilon>0$.

in the above domains. We adopt the VLSI model developed by Brent andKung [2] as our computation model, and consider time and area as our cost measures.

2. Area-Time Optinial Square Rooting Circuits for$T=\Omega((\log n)^{1+\epsilon})$In this section, we present an area-time optimal square rooting circuit for$T=\Omega((\log n)^{1+\epsilon})$.

Our square rooting algorithm is a modification of Mehlhorn-

Preparata’s area-time optimal division [4]. In our algorithm,$1/\sqrt{x}$the Newton iteration as a root of the equation $u^{-2}-x$ , andis calculated by the equa-$\sqrt{x}$is frst evaluated bytion $\sqrt{x}=x’(1/\sqrt{x})$ . The iteration rule is$u;+1= \frac{1}{2}u_{j}(3-xu_{i^{2}})$.

Since the Newton iteration has a self-correcting property, it is sufficient to computewith 2 $i+1$ -bit precision at the i-th stage of the iteration. Using this rule, it is immediately shown that there exists an$A$$T^{2}$-optimal, i.e., A $T^{2}=O(n^{2})$ , n-bit square rooting cir-

cuit for any $T\in[\Omega((\log n)^{2}), O(\sqrt{n})]$ . (See [5].)We now describe an area-time optimal square rooting circuit with computation

time

$T$

in the range $[\Omega((\log n)^{1+\epsilon}),O((\log n)^{2})]$ for any$\epsilon>0$, which is a modification of

Mehlhorn-Preparata’s area-time optimal divider with computation time in the samerange [4]. It is easily shown that almost all techniques for division are also applicable tosquare rooting, except for the polynomial approximation of inverses; the approximation,
however, essentially depends on the peculiarity of the Taylor expansion of $1/x$ .

Instead of the polynomial approximation, we present a new technique for approximation of inverse of square roots, based on the Newton iteration. Consider the algo-rithm described below:
$-2-$

$|$

117

AIgorithm INS $QRT2(x)$Given: an integer $s\in[2,l$) and an integer)$/4,1$ )Input: an l-bit number
Output: an $l+3$-bit number $v\in(1,2$ ] s.t.

$k=\lceil\log_{4}s\rceil$$x\in[$$v$

gives the leading $(l+3)$ bits of$1/\sqrt{x}$In line (3) of Algorithm INSQRT2,$f_{x}^{k}(u)=f_{x}(f_{x}$

(

,

. .

$f_{x}(u)$

. .

,

)

$)$

is evaluated in one step, instead of k-time evaluation of$f_{x}$.

The following theorem tells us the area-time preformance of computing
lines$(fl)-(f2)$$f_{x}^{k}(u)$inThis theorem is proved by using an extension of Beame-Cook-[7].

Hoover’s method [8].
PROPOSITION 1.

[Okabe, Takagi and Yajima, 1987] There exists a circuit whichin an n-bit number,computes the n-bit approximate value of a polynomial of degree$k$in time $O(\log(nk))$ and with area $O((nk)^{2})$ if $k\geq O(\log n)$ or

klog n) if $k<O(\log n)$ .

Since$f_{x}^{k}(u)$is a polynomial of degree$4^{k}(\approx s)$

in

$u$

and

$O$

(

$x,$

$n^{2}$$f_{x}^{k}(u)$can be computedin time $O(\log l)$ and with area $O(l^{2}s(s+\log l))$ . The number of iterations required inlines (2) $-(4)$ is$\lceil\log_{2}l\rceil/k$, and thus the total computation time of the square rooting cir-

cuit is $O(\log^{2}l/\log s)$ and the chip area is $O(l^{2}s(s+\log l))$ . We now obtain the following

lemma:LEMMA 1.

For any $s\in[2,l]$ , there exists a circuit which computes the l-bit inverse of

the square root of an l-bit number in time $O((\log l)^{2}/\log s)$ and has area $O((ls)^{2})$ if$s\geq O(\log l)$

or

$O$

( slog l) if $s<O(\log l)$ .
$l^{2}$We are now ready to describe our square rooting algorithm with optimal

performance for the computation time $T=O(\log^{1+\epsilon}n)$ for any $e>0$ .

-3-$AT^{2_{-}}$

technique [4] for square rooting can be written as:

Th

AIgorithm $INSQRT3(x)$Given: an integer sequence)
Input: an l-bit number$esuccessiver_{2}e_{i}finement(4)z.=vt$
Output: an l-bit number $v\in(1,2$] s.t.

$|$$l_{1}<l_{2}<\ldots<l_{J}=l$$x\in[)_{4}’,1$$vx=1+\epsilon,$$\epsilon<2^{-l}$(1) begin $v:=1$ ;(2) for $i:=1$ to do
begin $t.:=leftmost(l_{1}+1)$ bits of ;
(3)
$J$

$x$

(5)
(6)

bits of ;$leftmost(l_{1}+1)$-bit overinverse of square root of ); {i.e.,($v.:=$$x_{i}:=leftmost(l_{l}+1)$$z_{i}$$\lceil 2^{l_{l}+1}/\sqrt{x_{j}}\rceil\cdot 2^{-\langle l_{\iota}+1)}$

$x:$

}

as a

is

subroutine. To proceed the same discussion as in the case of division, we need the next
LEMMA 2.

If an l-bit number$x\in[l_{4}^{/},1$) has $(l’-1)$ zeros immediately to the right of

the leading 1, the l-bit inverse of the square root of$x$

can be computed in time$i^{Called}$$T=O(\log(l/l’)\cdot\log l/\log s))$and with area $A=O((ls)^{2})$ , for any $s\in[O(\log l),l/l’]$ .

$1mm_{p^{a:}lbit^{v}s^{;}.\square }^{In1i_{isi\circ n.of^{v_{A1gor^{;}ithm_{it}}}}^{)_{en}:=v}}t^{e_{oar^{(}e^{7}c_{l=^{d_{\frac{n}{(\log n)^{1+\epsilon_{S}};}}},=(\frac{n}{(\log n)^{1+\epsilon}})^{1}(i=1J)}^{n_{;}e(6).’Proof.
least

INS_{i}QRT2presented_{/}i_{s}ntheprevioussubs,ection}}}e.nd.;-(\log n)^{\epsilon(i-1)}$

In Algorithm INSQRT2, the first approximation $u_{0}=1$ has a precision of atbits. This implies that $0(\log(l/l’))$ iterations are sufficient to compute$l’$

For $i=2$

dition of

$x$

$J$

, it is obviously verified thatin the above lemma for$l=t_{i}$$x_{j}$$1/\sqrt{x}$$|$in Algorithm INSQRT3 satisfies the con-and $l’=l_{i-1}$ .

Only a straightforward discussion remains. Choose

where

$J$

is a largest value of$i$for which $s_{i}>2$ . (Indeed$J=\theta(1/\epsilon).$) Then it is provedthat the successive refinement modules based on Algorithm INSQRT3 can be implemented as a VLSI circuit with $T=O((\log n)^{1+\epsilon}),$ $A=O(n^{2}/((\log n)^{1+\epsilon})^{2})$ , though wewill omit the proof.

The Newton iteration techniques utilized here are completely equal to those of$|$

Mehlhorn-Preparata’s for their division algorithms [4]. Thus we have:
THEOREM 1.

For any fixed $e>0$ , the n-bit square root of an n-bit number can be$|$evaluated with optimal$A$ $T^{2}$-performance $o(n^{2})$ for any $T\in[\Omega((\log n)^{1+\epsilon}),O(\sqrt{n})]$ .
$- 4-$$\ovalbox{\tt\small REJECT}$

119

Note that similar results can also be derived for the computation of$k\sqrt{x}$(the k-th

root of x) for any fixed integer .
$k$

On the other hand, by choosing

$s$

as

$s=l^{\epsilon}(1\geq\epsilon>0)$

in Algorithm INSQRT2, theresulting circuit achieves $T=O((1/e)\log l)$ and $A=O(l^{2(1+e)})$ . Thus we get the following
result:THEOREM 2.

For any $e>0$ , there exists a circuit which computes the n-bit squareroot of an n-bit number in time $O(\log n)$ and with area $A=O(n^{2+\epsilon})$ .

3. Area-Time Efficient Evaluation of Exponentials and Logarithms
In this section, we will consider the area-time efficient implementation of circuits

for exponentials and logarithms.
First
$A$

we

$T^{2}=O(n^{2}\log^{2}n)$present

algorithm

anforevaluating

n-bit

logarithms

withbased on the arithmetic-geometric mean iteration of Gauss [6]. The

algorithms is as follows:
AIgorithm LOGI $(x)$Input: an n-bit number $x\in[1,2$ )Output: an n-bit number $y\in[0,\ln 2$ ) s.t. $|y-\ln x|<2^{-}"$(1) begin $a_{0}:=1;b_{0}:=2^{2-n}\cdot x^{-1};i:=0$ ;(2) while $a_{i}-b_{i}>2^{-n}$ do
begin$a_{t}+1:=(a_{t}+b_{j})/2$ ;
(3)
;
(4)
$i:=i+1$ ;(5)
end;(6) $y:=\pi/(2a_{i})$ -nln 2end.
$b_{i+1}:=\sqrt{a_{i}b_{i}}$Using the result on area-time optimal square rooting circuits, the following theorem
follows immediately.
THEOREM 3.

culatedwith

For any fixed$\epsilon>0$, the n-bit logarithm of an n-bit number can be calfor$AT^{2}=O(n^{2}(\log n)^{2})$$T\in[\Omega((\log n)^{2+\epsilon}), O(\sqrt{n}\log n)]$any

computation

time

.

Once an effcient evaluation of logarithms is established, the exponential function$y=\exp x$ can be evaluated as the root of the equation,$f(y)=\ln y-x=0$

for any , by applying the Newton method. The algorithm can be written as:
$x$

$-5-$

$s\sim v$

AIgorithm EXPI $(x)$Input: an n-bit number $x\in[0,h2$ )Output: an n-bit number $y\in[1,2$ ) s.t. I $y-\exp x|<2^{-n}$begin( -bit approximation of $\exp x$ ) $;i:=0$ ;while $|hy-x|>2^{-n}$ do

(1)
(2)

$l$

$y_{0}:=$

begin$y:+1:=y_{i}(1+x-\ln y)$ ;

(3)
(4)

$i:=i+1$ ;

end;

(5)

$y:=y_{i}$end.

We now consider the efficient implementation of this algorithm for each computa-tion time$T$in the range $T\in[\Omega((\log n)^{2+\epsilon}), O(\sqrt{n}\log n)]$ .
For any fixed $e>0$ , the n-bit exponential of an n-bit number can be

THEOREM 4.

evaluated

with$A$ $T^{2}=O(n^{2}(\log$T\in[\Omega((\log n)^{2+\epsilon}), O(\sqrt{n}\log n)]$forn)^{2})$

any

time

computation.

The theorem is easily proved for$T$in the range $[\Omega((\log n)^{4}),O(\sqrt{n}\log n)]$ , bychoosing $y_{0}:=1$ as an initial approximate value in line (1), utilizing the efficient Newton
iteration technique in [5] and evaluating logarithms in line (3) by Algorithm LOGI.

Thus

we

implementationtheconsider$T\in[\Omega((\log n)^{2+\epsilon}), O((\log n)^{4})]$

Algorithm

of

EXPI

for

.

First we propose a new algorithm for fast evaluation of a good initial l-bit approxi-

mate value of $\exp x$ in line(l). Let$x_{0}$

be the leftmost bit of , and$x$$x=x_{0}+x_{1}+x_{2}+\cdots+x_{k}$$x_{j}$$x$

be an l-bit binary number, and supposebe the $(s^{i-1}+1)$ -th to$s^{i}$$l=s^{k}$. Let-th bits of $x(i=1,\ldots,k)$ . Thenand therefore$\exp x=\exp(x_{0}+x_{1}+x_{2}+’\cdot\cdot+x_{k})=(\exp x_{0})(\exp x_{1})(\exp x_{2})\cdots(\exp x_{k})$.

This leads to the algorithm shown below:$\rangle\backslash$$- 6-$$\ovalbox{\tt\small REJECT}^{r}I$121

AIgorith $mEXP2(x)$Given: integers $s,k$ s.t.

$s\leq l$

and$k=\lceil\log_{s}l\rceil$

Input: an l-bit number $x\in[0,h2$ ),Output: an l-bit number $y\in[1,2$ ), s.t. I $y-\exp x|<2^{-l}$(1) begin(2) for $i:=1$ to

begin

$k$

{in parallel} do-th to -th bits of ;$y_{i}:=2$ l-bit approximation of $\exp x_{i}$ ;$end$ ;

(3)(4)$s^{i}$$x_{i}:=tl\iota e(s^{-1}+1)$

(5)

$x$

$y:= \prod_{1=1}y_{i}$end.

Consider the Taylor expansion of$\exp x_{i}$,$\exp x_{j}=1+x_{i}+(1/2!)x;^{2}+\cdot\cdot’+(1/n!)x_{i}^{n}+\cdots$

where$x_{j}<2^{-s^{i-1}}$

This means that the sum of the first $2l/s^{i-1}$ terms is an approximate value ofwith the precision of 21 bits, since$\exp x_{j}$hne (4) is a polynomial ofSince$x_{i}$$x_{i}$$x_{i^{2l/s^{i-1}}}<2^{-2l}$. Thus we may assume that$y_{i}$

inof degree 2 $l/s^{i-1}$ .

is an $(s^{i}-s^{i-1})$ -bit number, Step (4) is carried out in time $O(\log l)$ andwith area $O(l^{2}s(s+\log l))$ from Proposition 1. (Note that $(s^{j}-s^{i-1})’(2l/s^{i-1})=O(ls).$ )All$x_{i}’ s$can be calculated in parallel (lines (2) $-(4)$ ), and multiplied up pairwise in abinary-tree-form (line (5)). Thus computation time required in line (2) $-(4)$ is $O(\log l)$
and area $O(k\cdot l^{2}s(s+\log l))$ . Using $(2k-1)$ multipliers with time $O(\log l)$ and area

in line (5), their product can be obtained in time$O(k\cdot l^{2})$.

$O(k\prime l^{2}s(s+\log l))$LEMMA 3.

tial

of

the

Hence

an

computation

time$O((\log k)\cdot(\log l))$$O((\log k)\cdot(\log l))$and$O(l^{2})$and with area
the

chip

area

in line (5). Thus we have the following lemma:

For any $s\in[2,l]$ , there exists a circuit which computes the l-bit exponen-l-bit$O((ls)^{2}\log l/\log s)$

numberintime$O((\log l)\log(\log l/\log s))$

and

has

area

if $s\geq O(\log l)$ or $0(l^{2}s\log^{2}l/\log s)$ if $s<O(\log l)$ .

Let us continue the proof of the theorem for $T\in[\Omega(\log^{2+\epsilon}n), O(\log^{4}n)]$ . Let $l=n/T$,

and consider the Algorithm EXPI for this . From the above lemma for $s=2$ , it follows$l$that there is a circuit, say$F_{E}$, which computes the initial l-bit approximation of $\exp x$

in line (1) in time
$T_{B}=O((\log l)\log\log l)=O((\log n)\log\log n)$

and has area
$-7-$

122
$*\ovalbox{\tt\small REJECT}_{;}^{y}$$A_{E}=O(l^{2}\log^{2}l)=O(n^{2}\log^{2}n/T^{2})$.

feedback, and the later $m-i_{0}$ iterations on their own circuits

former$i_{0}$$F_{i_{0}+1},\ldots,F_{m}$respectively.

iterations in time $T_{F}=O(\log^{2+c}n)$ and has area

$A_{F}=O((l\prime 2\eta 2,(i_{0}+\log_{2}l)^{2}/(\log^{2+c}n)^{2})=O(n^{2}\log^{2}n/T^{2})$The circuit.

$F$choose the computation time of $F_{m-j}$ as$T_{m-j}=T/2^{j/2}$

for each $j=0,\ldots,m-\iota_{0}-1$ . Note that$\tau_{;_{0+1}}=T/l^{(m-i_{0}+1)/2}=O((T\log^{2+c}n)^{l4})>O(\log^{2+c}n)$From the result of Theorem 3, there exists a circuit$F_{j}$.

which has area

$A_{m-j}=O((2^{m-j})^{2}(m-j)^{2}/T_{j}^{2})=O(((2^{m})^{2}m^{2}/T^{2})2^{-j})=O((n^{2}\log^{2}n/T^{2})2^{-j})$for any $j=0,\ldots,m-i_{0}-1$ .
These imply that our circuit has area$|$

$A_{r}=A_{E}+A_{P}+ \sum_{j}A_{j}=O(n^{2}\log^{2}n/T^{2})+\sum_{j}O((n^{2}\log^{2}n/T^{2})2^{-j})=O(n^{2}\log^{2}n/T^{2})$$1$and computes the n-bit logarithm in time
$T_{r}=T_{E}+T_{F} \cdot i_{0}+\sum_{j}T_{j}=O((\log n)\log\log n)+O(\log^{1+c}n)\cdot O(\log n)+\sum_{j}T/2^{-j/2}=O(T)$,since $T\geq O(\log^{2+c}n)$ . Thus the theorem follows.
Next, let us consider the time-optimal, i.e.
and logarithms, with suboptimal$A$$T^{2}$$T=O(\log n)$ ,circuits for exponentials-performance.

Consider Algorithm EXP2 in the previous subsection, and choose$1\geq\epsilon>0$.

The Circuit in Lemma 3 operates in time

$A=O(l^{2+(2/3)\epsilon})$. Thus we have:-8-$T=O(\log l)$$s=l^{\epsilon/3}$

for any

and has area

123

For any $e>0$ , there exists a circuit which computes the n-bit exponen-

THEOREM 5.

tial of an n-bit number in time $O(\log n)$ and with area$O(n^{2+\epsilon})$We will show the corresponding result for logarithms.

number and suppose
where$\eta_{j}$is an$s^{i}$$l=s^{k}$. It is easily verified that-bit number which has at least$x$$s^{i-1}$

.

Let

$x$

be an l-bit binary

can be written as$x\simeq\eta_{1}\cdot\eta_{2}\cdot\cdot’$$\eta_{k}$,consecutive $0’ s$ immediately to theright of the leading 1. Since$\ln x=\ln(\eta_{1}\cdot\eta_{2}’\cdot’\cdot\eta_{k})=\ln\eta_{1}+\ln\eta_{2}+\cdots+\ln\eta_{k}$$\ln x$,is obtained as:
AIgorithm $LOG2(x)$Given: integers
s.t.
$k$

$s,$

and

$s\leq l$$k=\lceil\log_{s}l\rceil$

Input: an l-bit number $x\in[1,2$ )Output: an l-bit number $y\in[0,h2$ ), s.t.

(1) begin $x_{0}:=x;y:=0$ ;(2) for $i:=1$ to do
begin$\eta_{i}:=leftmosts^{i}$ bits of;(3)$x_{i}:=x/\eta;$ ;(4)-bit approximation of(5)$y:=y+y_{i}$ ;
(6)
end
end.

$|y-hx|<2^{-l}$$k$$x_{i-1}$$y_{i}:=l$$h\eta_{i}$Consider the Taylor expansion of;$\ln\eta_{j}$,$\ln\eta_{i}=\ln(1+h_{i})=h_{i}-(1/2)h_{j}^{2}+\cdot’\cdot+(-1)^{n}(1/n)h_{i}^{n}+\cdots$

As is obvious in line (3) $-(4),$an approximate value of ln
is a polynomial of

imation of$\ln\eta$$h_{i}$$h_{i}<2^{s^{i-1}}$$\eta_{i}$This means that sum of the frst 2 $l/s^{i-1}$ terms iswith the precision of $2l$ bits. Thus we may assume thatof degree 2 $l/s^{i-1}$ . Since; can be computed in time$h_{i}$$y_{j}$

is an $(s^{i}-s^{i-1})$ -bit number, l-bit approx-$O(\log(s^{i}\cdot(2l/s^{;-1})))=O(\log(ls))$

and witharea $O((s^{i}\cdot(2l/s^{i-1}))^{2})=O((ls)^{2})$ (Proposition 1).
Choosing$s=l^{\epsilon/3}$

for any $1\geq e>0$ , this satisfies $k=\theta(1/e)$ . Step (4) can be carriedout in time $O(\log l)$ using a divider with area$y_{j}$$O(l^{2+(2/3)\epsilon})$(see [4]). For each $i=1,\ldots,k$ ,can be evaluated (line (5)) and added up (line (6)) in time $O(\log(ls))=O(\log l)$ andwith area$O((ls)^{2})=O(l^{2})$Thus the total computation time of our circuit is$O(k\log l)=O(\log l)$ and the chip area is

THEOREM 6.

$O(l^{2+(2/3)\epsilon})$. We get:For any $e>0,$ there exists a circuit which computes the n-bit loga-

rithms of an n-bit number in time $O(\log n)$ and with area

-9-$O(n^{2+\epsilon})$.

$x,\dot{}..\ovalbox{\tt\small REJECT}$

4. Considerations$\ovalbox{\tt\small REJECT}$

Trigonometric functions such as sines, cosines and arctangents can be evaluatedfrom exponentials and logarithms with complex arguments, since$\log(v+iw)=\log|v+iw|+i’\arctan(w/v)$ ,$\exp i\theta=\cos\theta+i’\sin\theta$

Acknowledgementscomments and criticism on parallel and VLSI algorithms.

References
Sci., Carnegie-Mellon Univ., Pittsburgh, Pa $(J_{an}$ . 1979).

521-534 (1981).
Computation$Time^{\prime t}$, Inform. and Control, 58, 137-156 (1983).

Comput. 72, 270-282 (1987).

163-167 (1984).
Press, New York, 151-176 (1976).

Number System“, 1987 LA Symposium in Summer, (July 1987), in Japanese.
lems“, SIAM J. Comput., 15-4, 994-1003 (1986).