Electronics Guide

Rotation Mode

Rotation mode drives the angle accumulator to zero. The algorithm chooses d_i from the sign of the residual angle z_i, so the sequence of micro-rotations performs a binary search over the elementary angles until the requested rotation has been consumed.

The step is not a true rotation. A true rotation by α_i multiplies the coordinates by cos α_i as well as shearing them, and CORDIC drops that factor because it is not a power of two. What remains is a pseudo-rotation that turns the vector through the correct angle but lengthens it by 1/cos α_i = √(1 + 2⁻²ⁱ) at every step. Starting from (1, 0) and rotating through θ, circular rotation mode therefore produces

(K cos θ, K sin θ),

where K is the accumulated gain described below, and not (cos θ, sin θ). An implementation must either pre-scale the starting vector or compensate afterward. Subject to that compensation, the mode delivers sine and cosine together from one pass through the pipeline, which is its principal advantage over methods that compute the two functions separately.

Circular rotation mode converges only for angles it can reach. The reachable set is bounded by the sum of all the elementary angles, roughly 1.7433 radians or about 99.88 degrees. That exceeds a right angle but falls well short of a half turn, so a general argument still needs pre-reduction: the implementation folds the argument into the first quadrant by exchanging and negating coordinates, an exact operation that costs no accuracy, and records the quadrant for reconstruction.

Vectoring Mode

Vectoring mode performs the inverse operation. The direction decision now comes from the sign of y_i, so the iteration rotates the input vector onto the positive x axis while the accumulator collects the total angle turned. Given coordinates (x, y), the algorithm drives y toward zero and leaves

x_n ≈ K √(x² + y²), z_n ≈ z₀ + arctan(y/x).

Only the magnitude carries the gain. The angle accumulator is a running sum of tabulated constants and is never scaled, so an arctangent or a phase measurement taken from vectoring mode needs no gain compensation at all. The asymmetry is useful in practice: a phase detector built on CORDIC is accurate to its word length, while a magnitude taken from the same hardware is too large by the factor K unless the designer corrects for it.

Vectoring mode underpins rectangular-to-polar conversion in signal processing. Phase discriminators, automatic gain control, frequency estimators, and vector normalizers all reduce to a vectoring CORDIC, and the absence of a multiplier lets the block run at the sample rate in fabric that has no multiplier to spare.

Hyperbolic CORDIC

Setting m = −1 and replacing the elementary angles with the inverse hyperbolic tangents of 2⁻ⁱ yields hyperbolic CORDIC, which computes the hyperbolic sine, cosine, and tangent and their inverses. The exponential follows from the sum of the hyperbolic sine and cosine, the natural logarithm from twice the inverse hyperbolic tangent of (x − 1)/(x + 1), and the square root from a related vectoring initialization.

The hyperbolic elementary angles shrink fast enough that the naive sequence fails to converge: the angles that remain can become too small to cancel the residual. The standard remedy repeats selected iterations, at indices 4, 13, and 40, and thereafter at each index 3k + 1 derived from the previously repeated index k. With those repetitions, an iteration starting at i = 1 covers residual arguments up to roughly 1.1182 and accumulates a gain of about 0.82816, whose reciprocal is about 1.20750. Note the direction: the hyperbolic gain shrinks the vector, whereas the circular gain grows it.

Linear CORDIC

Setting m = 0 makes the elementary angles the plain powers of two and turns the recurrence into shift-and-add multiplication in rotation mode or division in vectoring mode. The linear mode has a gain of exactly one, because the pseudo-rotation no longer changes the vector length, and its convergence range follows from the sum of the shift weights in the chosen sequence. Dedicated multipliers have made linear CORDIC uncommon in new designs, but it still matters in two places: it lets a unified CORDIC unit supply multiplication and division without extra hardware, and it makes the structure of the algorithm easy to see.

CORDIC Scaling Factor

The pseudo-rotations of circular CORDIC multiply the vector length at every step, and the gain accumulated over an unbounded number of iterations is the infinite product

K = ∏ √(1 + 2⁻²ⁱ) ≈ 1.646760,

taken over i = 0, 1, 2, and onward. Its reciprocal, roughly 0.607253, is the constant an implementation usually stores. The product converges quickly: twelve iterations reach K to about seven decimal places, and twenty-four reach the limit to the last bit of double precision. Descriptions that quote 0.6073 as the scaling factor are quoting the reciprocal, not the gain.

Three compensation strategies are in use. Pre-scaling loads the initial vector as (0.607253, 0) rather than (1, 0), which costs nothing at all when the iteration count is fixed, because K then depends only on that count and is a constant known at design time rather than a runtime quantity. Post-scaling multiplies the two outputs by the same constant, at the price of the multiplier CORDIC was chosen to avoid. Absorption folds the constant into whatever follows, such as the amplitude coefficient of a numerically controlled oscillator or the gain of the next filter stage, and is generally free.

Scaling-free variants attack the problem differently, adding correction terms or repeating iterations so that the product collapses to a power of two, which a shift undoes exactly. Whether the extra iterations cost less than one multiplication depends entirely on the target technology.

Convergence Rate

CORDIC resolves roughly one bit of the result per iteration. Convergence is therefore linear, and an n-bit result costs about n iterations plus a few guard iterations and guard bits to absorb the rounding accumulated in the coordinate updates. This is the algorithm's defining trade-off, and it stands in sharp contrast to the Newton-Raphson and Goldschmidt methods discussed below, whose convergence is quadratic and which roughly double the number of correct bits per step. CORDIC wins when a multiplier is unavailable or expensive and the per-iteration cost of an adder and a shifter is low, so that many cheap steps beat a few expensive ones. It loses in any datapath that already contains a fast fused multiply-add.

Unified CORDIC

Unified CORDIC architectures carry the coordinate-system parameter m as an input, so a single datapath serves the circular, linear, and hyperbolic modes. The elementary-angle table grows to hold all three sequences, and the control logic must handle the repeated hyperbolic iterations, but the adders, shifters, and registers are shared. The approach maximizes hardware reuse in systems that need several function types, accepting some control complexity and a slightly longer critical path in exchange for reduced area.

Taylor Series

A Taylor expansion represents the function as an infinite sum of terms built from successive derivatives at one point. The series for sine about zero is

sin x = x − x³/3! + x⁵/5! − x⁷/7! + …

and it converges for every argument. Convergence of the infinite series is not the property an implementer needs, however. A hardware approximation truncates at a fixed degree, and the truncated Taylor polynomial concentrates its accuracy at the expansion point, where it matches the function to high order, while its error grows steeply toward the ends of the interval.

A minimax polynomial of the same degree spreads the error evenly and does markedly better. Take the odd polynomial of degree seven for sine on the interval from 0 to π/4, which is the natural interval after trigonometric reduction. Truncating the Taylor series leaves a worst-case absolute error of about 3.1 × 10⁻⁷, roughly 21.6 bits. The minimax polynomial with the same four coefficients, and therefore the same three multiplications and three additions to evaluate, holds its worst-case error to about 1.2 × 10⁻⁹, roughly 29.6 bits. The improvement is a factor of about 260, close to eight bits, at identical runtime cost.

Use the Taylor series to understand the function, to bound a remainder analytically, or to supply a starting point for a coefficient search. Do not ship it as the approximation unless the interval is small enough that the difference no longer matters.

Chebyshev Approximations

Chebyshev polynomials of the first kind are orthogonal on the interval from −1 to 1 with respect to the weight 1/√(1 − x²). Expanding a function in that basis and truncating gives an approximation whose error, unlike the Taylor remainder, oscillates across the whole interval instead of piling up at the ends.

Chebyshev truncation is near-minimax rather than minimax. Its worst-case error exceeds the true optimum by a factor that grows only logarithmically with the degree, which at the degrees used in function evaluation means a modest penalty, and many designs stop there. Coefficients follow from a discrete cosine transform of the function sampled at the Chebyshev nodes, a calculation cheap enough to run inside a design loop, after which the polynomial converts to the power basis for evaluation. That conversion becomes ill-conditioned at high degree, so tools that need many coefficients evaluate in the Chebyshev basis directly using the Clenshaw recurrence.

Minimax Polynomials

The Remez exchange algorithm computes the true minimax polynomial, the one that minimizes the maximum error over the interval. Its solution is characterized by equioscillation: for a polynomial of degree n the error reaches its maximum magnitude, with alternating sign, at no fewer than n + 2 points. The algorithm iterates by guessing those points, solving a small linear system for the coefficients and the error level, then moving the points to the actual error extrema and repeating.

Two refinements matter in practice. First, minimizing relative rather than absolute error is usually what the specification demands, and Remez accommodates that through a weight function. Second, the coefficients must eventually be rounded to the target number format, and rounding an optimal real coefficient does not give the optimal machine-representable one. Design tools therefore re-optimize over the representable coefficient set, recovering a fraction of a bit that naive rounding throws away. The search runs once, offline; only the evaluation runs on the hardware.

Horner's Method

Horner's method evaluates a polynomial by repeated factoring:

a_nxⁿ + a_n−1xⁿ⁻¹ + … + a₁x + a₀ = (…((a_nx + a_n−1)x + a_n−2)x + …)x + a₀

This costs n multiplications and n additions for degree n, and it maps directly onto a fused multiply-add unit that computes ab + c with a single rounding. The single rounding is worth having, since it removes one error source per term from the accuracy budget.

Horner's weakness is latency. Every step depends on the previous one, so the evaluation takes n dependent multiply-add latencies end to end. Where throughput is what matters, a pipeline hides this by overlapping independent evaluations. Where a single result must arrive quickly, Estrin's scheme rebalances the same polynomial into a tree of independent products, shortening the dependency chain to roughly the logarithm of the degree at the cost of a few extra multiplications. Even and odd functions offer a further saving: writing sine on the reduced interval as x times a polynomial in x² halves the number of Horner steps.

Rational Approximations

A rational approximation is the quotient of two polynomials, and for a given total number of coefficients it often beats any single polynomial, especially for functions with poles or strong asymptotic behavior, such as the tangent near a right angle or the inverse hyperbolic tangent near its endpoints.

The Pade approximant is the rational function whose Taylor coefficients match the target's to the highest possible order. It inherits the Taylor weakness exactly: excellent at the expansion point, poor at the ends of the interval. The useful form is the minimax rational approximation, computed by a rational variant of the Remez algorithm, which equioscillates in the same way as the polynomial case.

The cost is a division, which in most datapaths carries several times the latency of a multiplication and is not pipelined as deeply. A rational form pays for itself when the division hardware is already present and idle, or when the polynomial degree it saves is large enough to offset the division. Otherwise a higher-degree polynomial is usually cheaper.

Direct Table Lookup

The simplest table stores the function at uniformly spaced arguments and returns the nearest entry. Accuracy follows from the spacing alone: indexing with k bits of the argument and storing the best constant on each subinterval gives roughly k bits of result, so a 256-entry reciprocal table over the interval from 1 to 2 delivers about eight to nine bits. Reaching the twenty-four bits of single precision by direct lookup alone would take about sixteen million entries, which is why direct lookup survives only where the argument is already short, as in a fixed-point sine table for a numerically controlled oscillator, or where a coarse estimate is all that is wanted.

Linear Interpolation

Interpolating linearly between adjacent entries changes the error from first order in the spacing to second order. The worst-case error of linear interpolation is bounded by one eighth of the maximum second derivative times the square of the spacing, so halving the spacing gains two bits rather than one. For the reciprocal over the interval from 1 to 2, whose second derivative is at most 2, thirty-two segments give about twelve bits and sixty-four give about fourteen. The hardware reads two entries, or one entry and a stored slope, forms the difference, multiplies by the fractional part of the argument, and adds. Storing the slope beside the value avoids the second read and shortens the critical path, at the cost of doubling the table width.

Bipartite Tables

Bipartite tables exploit a first-order expansion to replace one large table with two small ones. Split the n-bit argument into three fields of about n/3 bits each. One table, addressed by the upper two fields, holds the function value at the truncated argument; the second, addressed by the outer and the lowest fields, holds the first-order correction. Adding the two outputs reconstructs the result, and the neglected second-order term is of order 2⁻ⁿ, so the method is accurate to roughly n bits and can be built to round faithfully.

Each table is addressed by about 2n/3 bits, so the pair holds on the order of 2^2n/3 entries instead of 2ⁿ. At sixteen bits of accuracy that is a few thousand entries against sixty-five thousand for direct lookup, a compelling saving. The scaling is still exponential: pushing the same structure to the twenty-four bits of single precision raises the pair to the order of a hundred thousand entries. Single-precision reciprocal and reciprocal-square-root units therefore add interpolation, raise the polynomial degree, or split the argument further rather than rely on a plain bipartite table.

Multipartite Tables

Multipartite methods generalize the split to several correction tables, each addressed by a different combination of argument fields and each contributing a finer term. The design problem becomes a discrete optimization over the number of tables, the field boundaries, and the word width of each table, subject to a target error and a faithful-rounding requirement. Automated generators search that space and emit the table contents together with the adder network. The payoff is real but bounded: past roughly twenty-four bits the table sizes grow faster than a polynomial evaluation of equivalent accuracy, and the polynomial wins.

Symmetric Table Methods

Symmetry shrinks tables without costing accuracy. Sine and cosine need only a quarter turn of storage, because the remaining quadrants follow from sign changes, and using the identity that relates the sine of the complement to the cosine reduces the stored range to a single octant. The odd symmetry of the hyperbolic tangent halves its table in the same way, and the correction table of a symmetric bipartite scheme folds about its midpoint to halve again.

The exponential and the logarithm do not share a table, despite being inverses, because a table uniform in the argument of one is not uniform in the argument of the other. What they share is a decomposition. Writing the exponential of a sum as the product of exponentials lets two small tables, indexed by disjoint fields of the reduced argument, replace one large table, and the logarithm admits the mirror-image decomposition through the logarithm of a product.

Additive Range Reduction

Additive reduction subtracts an integer multiple of a constant. For the trigonometric functions the constant is π/2 or π, and the reduction computes r = x − N C, where N is the integer nearest x/C and the low bits of N select the octant for reconstruction.

The difficulty is that C cannot be represented exactly. If the stored constant differs from the true π/2 by δ, the reduced argument inherits an absolute error of about Nδ, and N grows in proportion to the argument, while r itself stays small. The result is catastrophic cancellation: the absolute error holds roughly constant as the quantity it contaminates shrinks, so the relative error of the reduced argument grows without bound as the argument grows. Obtaining p correct bits in r for an argument of magnitude 2^E requires the reduction constant to be carried to roughly E + p bits, plus a margin for the worst case in which x lands very close to a multiple of π/2.

Cody-Waite reduction meets this by splitting the constant into two or three parts, each with enough trailing zero bits that the product N C₁ is exact in the working format. Subtracting the parts in decreasing order keeps the intermediate results exact until the last and smallest term. The technique is cheap and handles arguments up to a moderate size, after which the exactness condition fails and the accuracy collapses; the threshold depends on how many bits the split reserves for N.

Multiplicative Range Reduction

Multiplicative reduction divides out a scale factor rather than subtracting an offset, and for the exponential it is exact. Writing the exponential of x as the exponential of a reduced argument r multiplied by 2ⁿ, where n is the integer nearest x/ln 2, leaves r in a short interval around zero. The final scaling by a power of two is exact in binary floating point and introduces no error of its own. Only the overflow and underflow boundaries need care, where the scaling must be split into two steps to avoid a spurious intermediate subnormal.

The logarithm inverts the process. A binary floating-point input already carries its own multiplicative split, so the exponent contributes the exponent times ln 2 and the significand contributes the logarithm of a value in the interval from 1 to 2. Implementations usually shift that interval to be centered on 1, because the logarithm has a zero there and a symmetric interval keeps the relative error bounded near it.

Payne-Hanek Reduction

Payne-Hanek reduction removes the argument-size limit. Instead of subtracting a rounded constant, it multiplies the argument by a stored high-precision representation of 2/π and keeps only the bits of the product that survive into the fractional part, discarding both the leading bits, which contribute whole multiples of a half turn, and the trailing bits, which fall below the target precision. Because the extracted window depends on the argument's exponent, the stored constant must span the whole exponent range: for IEEE 754 binary64, whose exponents exceed 1000, the table runs to more than a thousand bits.

The method returns a reduced argument accurate to full precision for every finite input, including the notorious arguments whose product with 2/π falls within a hair of an integer. The cost is a multiword multiplication and the storage for the constant, so implementations run the Cody-Waite path for small arguments and fall back to Payne-Hanek only when the exponent demands it.

Trigonometric Identities

The identities sin(x + π/2) = cos x, sin(−x) = −sin x, and cos(−x) = cos x map any argument onto the interval from 0 to π/4, where both sine and cosine are smooth and well conditioned and a short odd or even polynomial suffices. The low three bits of the reduction multiple N give the octant, which selects whether to return the sine or the cosine polynomial and with which sign.

Reducing to π/4 rather than π/2 matters. On the wider interval the cosine approaches zero, where relative accuracy becomes expensive, while the sine approaches 1, where its vanishing derivative makes a fixed absolute error harmless. Splitting at π/4 keeps whichever function is evaluated away from its zero. The same reasoning explains why the tangent is usually built from separately reduced sine and cosine values rather than from a single polynomial that must chase a pole.

Exponential and Logarithmic Identities

The exponential of a sum equals the product of the exponentials, which separates the reduced argument from the power-of-two scaling and, in table-assisted designs, separates it again into a table index and a very small residue. A common arrangement splits the reduction multiple into a whole part and a fraction over 2^k, looks the exponential of the fraction up in a table of 2^k entries, and evaluates a polynomial of degree two or three on the residue that remains. The residue interval shrinks by a factor of 2^k, so a few table bits buy several polynomial degrees.

The logarithm of a product equals the sum of the logarithms, which is what lets the exponent and the significand be handled separately. Near an argument of 1 the naive form loses accuracy to cancellation, because the logarithm is close to zero while the argument is close to 1; libraries therefore expose a separate routine for the logarithm of one plus a small quantity, and use it internally whenever the significand lands near unity.

Hyperbolic Identities

The hyperbolic sine and cosine follow from the exponential as (exp x − exp(−x))/2 and (exp x + exp(−x))/2. For large arguments the negative exponential is negligible and the formulas are well conditioned. For small arguments the hyperbolic sine formula subtracts two nearly equal quantities and loses most of its significant bits, so implementations switch to a direct odd polynomial below a threshold of the order of the square root of the machine epsilon. The hyperbolic tangent needs the same treatment at both ends: a polynomial near zero, and saturation to ±1 once the argument grows large enough that the correction falls below the last bit.

Newton-Raphson Iteration

Applying Newton's method to a function whose root is the reciprocal of d gives a recurrence free of division:

x_n+1 = x_n(2 − d x_n)

Its behavior is transparent in terms of the residual e_n = 1 − d x_n. Substituting gives

1 − d x_n+1 = (1 − d x_n)², that is e_n+1 = e_n²,

which is the whole story. Convergence is quadratic, so each step squares the residual and approximately doubles the number of correct bits. Convergence also requires |e₀| < 1, which for a positive divisor means the initial estimate must satisfy 0 < x₀ < 2/d. Outside that interval the residual grows and the iteration diverges, so the seed table is not merely an optimization but a correctness requirement.

The doubling is approximate rather than exact. Each step rounds two products, so a step that should produce 2n correct bits typically produces a bit or so fewer, and implementations either carry guard bits or budget one extra iteration. Starting from an eight-bit table estimate, three iterations pass through roughly sixteen, thirty-two, and sixty bits, which covers the fifty-three-bit significand of double precision with margin. The step needs two dependent multiplications and one subtraction, all of which a fused multiply-add executes cleanly.

The quotient follows by multiplying the dividend by the reciprocal, but that product is not necessarily the correctly rounded quotient: an error of a fraction of a unit in the last place in the reciprocal can push the product to the wrong side of a rounding boundary. Implementations that must meet the IEEE 754 requirement for division therefore finish with a residual computation, forming the remainder with a fused multiply-add and using its sign to choose between the candidate quotient and its neighbor.

Goldschmidt Division

Goldschmidt division scales dividend and divisor by the same factor at every step, preserving their ratio while driving the divisor toward unity. When the divisor reaches one, the dividend has become the quotient.

The principal advantage over Newton-Raphson is parallelism. The two multiplications within a Goldschmidt iteration are independent and can execute simultaneously, whereas the two multiplications of a Newton-Raphson step are sequentially dependent. This independence keeps a pipelined multiplier fully utilized and lets several divisions overlap in the pipeline, raising throughput. The trade-off is accuracy: Goldschmidt is not self-correcting, so rounding errors accumulate across iterations, and implementations retain extra guard bits to compensate. Newton-Raphson, by contrast, recomputes the residual from the original divisor at every step, so an error introduced in one iteration is squared away by the next.

Series Expansion Methods

Writing y = 1 − d expands the reciprocal as a geometric series,

1/d = 1/(1 − y) = 1 + y + y² + y³ + …, valid when |1 − d| < 1,

a condition that holds once the divisor has been scaled into a neighborhood of one. Summing term by term converges linearly and is of little practical use. The series factors, however:

1 + y + y² + y³ + … = (1 + y)(1 + y²)(1 + y⁴)(1 + y⁸) …

Each new factor costs one squaring and doubles the number of series terms accumulated, which recovers quadratic convergence. This factored form is the Goldschmidt iteration written another way, which is why the two methods share both their parallelism and their sensitivity to accumulated rounding.

Table-Based Initial Approximation

A table indexed by the leading significand bits of the divisor gives the fastest seed. Eight index bits and 256 entries yield roughly eight bits of accuracy, comfortably inside the convergence interval and enough that three Newton-Raphson iterations reach double precision. The design trade is direct: each extra index bit doubles the table and buys one bit, whereas each iteration costs two multiplication latencies and buys a factor of two in bits. Tables win at low precision and iterations at high.

Linear Approximation

Storing a slope beside each entry and interpolating makes the error quadratic in the segment width, so a table of thirty-two to sixty-four segments across the significand range reaches about twelve to fourteen bits. That is the accuracy a direct table of several thousand entries would need, obtained instead from a few dozen entries and one multiply-add, and it removes one full iteration from a double-precision sequence.

Quadratic Approximation

A second-order fit stores three coefficients per segment and evaluates ax² + bx + c with two multiply-adds. The error now falls with the cube of the segment width: for the reciprocal across the significand range, thirty-two segments reach about twenty bits and sixty-four about twenty-three. Twenty bits of seed accuracy removes another iteration from a double-precision divider, and it is enough on its own for a single-precision result of graphics quality.

Hardware Reciprocal Instructions

Most general-purpose instruction sets expose the seed generator directly. The x86 reciprocal and reciprocal-square-root estimate instructions return a result whose relative error is below roughly 2⁻¹², and the later AVX-512 forms tighten that to roughly 2⁻¹⁴. Arm NEON pairs an eight-bit estimate instruction with a dedicated step instruction that performs one Newton-Raphson refinement, so the software sequence is an estimate followed by as many steps as the required precision demands. Compilers emit these for fast reciprocal and inverse-square-root idioms when the language permits relaxed accuracy, and fall back to the fully rounded divide instruction when it does not.

Graphics Processing Units

GPUs include special function units that compute the reciprocal, reciprocal square root, sine, cosine, base-two exponential, and base-two logarithm. These serve the transcendental demand of rasterization and shading, where coordinate transforms, lighting, and texture filtering call for billions of such operations per second and the result often feeds an eight-bit or ten-bit display channel.

The units implement quadratic interpolation over a segmented table with minimax-fitted coefficients, returning a result in a few clock cycles. They deliver reduced precision rather than the correctly rounded single-precision value: the hardware-accelerated intrinsics carry an error of a few units in the last place, leaving roughly twenty-two of the twenty-four significand bits trustworthy, and the trigonometric intrinsics additionally restrict the argument range they accept, leaving the caller to reduce large arguments. For rendering, the trade is invisible in the output. Compute workloads that need the full single-precision guarantee call the library routine instead, which performs proper reduction and refinement on top of the same hardware estimate.

Digital Signal Processors

DSP architectures often include specialized support for the trigonometric functions used in communication signal processing. Phase rotators, coordinate converters, and numerically controlled oscillators implement optimized algorithms for their specific functions, typically a CORDIC or a quarter-wave sine table with interpolation. Fixed-point arithmetic predominates, with word widths matched to the spurious-free dynamic range the application requires. The phase accumulator is usually far wider than the table index, and truncating the phase to the table depth is itself a designed source of spurious tones.

Scientific Computing Accelerators

High-performance computing systems may include elementary-function support that maintains full double-precision accuracy while achieving high throughput. These implementations combine table-assisted polynomial evaluation with complete range reduction, and they are vectorized, because the cost of a scalar transcendental call is dominated by the branches that handle special cases. Vector math libraries meet that by evaluating the common path for every lane and patching the rare cases afterward, which keeps the pipeline full at the cost of some wasted work.

Neural Network Accelerators

Machine learning workloads require activation functions including the logistic sigmoid, the hyperbolic tangent, the Gaussian error linear unit, and the exponential inside a softmax. Neural network accelerators implement these with reduced-precision approximations matched to the noise tolerance of the models, usually piecewise linear segments or small lookup tables in the same narrow format as the weights. The softmax is the exception that needs care: its exponentials must be shifted by the maximum of the inputs before evaluation, or they overflow, and that maximum is a serial dependency across the whole vector.

FPGA Implementations

Field-programmable gate arrays allow the implementation to be chosen per application rather than per chip. CORDIC maps well onto fabric, since a fully unrolled pipeline is nothing but adders and hard-wired shifts, and it delivers one result per clock at the cost of one adder stage per bit of precision. Where DSP blocks are available, a table-assisted polynomial is usually smaller and shorter in latency for the same accuracy, and where block RAM is available, a bipartite or multipartite table can beat both.

The freedom to choose extends to the number format. An FPGA implementation need not use a standard floating-point width; it can carry exactly the integer and fraction bits the error budget demands, which often halves the resource count compared with a single-precision implementation of the same specification.

Accuracy Requirements

Required accuracy varies by orders of magnitude across applications. A shader tolerates several units in the last place. A control loop needs monotonicity more than accuracy, because a non-monotone approximation can make a feedback system chatter. A numerical library aims at a fraction of a unit in the last place. The total error is the sum of three contributions, and all three must be counted: the reduction error, which the earlier discussion showed can dominate; the approximation error of the polynomial or table; and the rounding error of the evaluation itself. A design that budgets only the second will miss its target.

Correct Rounding and the Table-Maker's Dilemma

IEEE 754 requires correct rounding for addition, subtraction, multiplication, division, and the square root. It does not require it for the elementary functions, which the standard lists among the recommended operations, and most historical libraries do not provide it. The reason is the table-maker's dilemma. The exact value of a transcendental function at a representable argument can fall arbitrarily close to the boundary between two representable results, and no general bound on that closeness is known in advance. Deciding which way to round may therefore require the intermediate computation to carry far more than the target precision, with no way to know how much before trying.

The practical answer is Ziv's strategy: evaluate with a modest number of extra bits, test whether the result is far enough from a rounding boundary to be decided, and fall back to a slower and wider evaluation in the rare cases that are not. Exhaustive search has established the hardest cases for single precision and for many double-precision functions, which allows a fixed working precision to be proven sufficient, and recent library work has produced correctly rounded implementations of the common functions at speeds competitive with faithful ones. A claim that a routine is accurate to the last bit should be read as a claim about which of these regimes it implements.

Performance Requirements

Latency-critical paths favor short dependency chains: table lookup with low-degree interpolation, Estrin evaluation rather than Horner, and seed accuracy bought with memory rather than iterations. Throughput-critical paths favor deep pipelines and methods whose steps are independent, which is why Goldschmidt appears in dividers that must sustain a high issue rate. Energy-constrained designs count total operations and memory accesses, and often reach the opposite conclusion from the latency-driven one, because a large table costs more energy per access than several arithmetic operations.

Resource Constraints

Available resources narrow the field quickly. Without a multiplier, CORDIC or a digit-recurrence method is the only option. With a fast fused multiply-add and little memory, polynomial evaluation wins outright. With abundant on-chip memory and a latency target, table methods win. The interesting designs sit in the middle and combine all three: a table for the seed, a short polynomial for the correction, and an exact power-of-two scaling to put the result back in range.

Special Cases

Elementary functions must handle zero, the infinities, quiet and signaling NaNs, signed zeros, and subnormal inputs and outputs. IEEE 754 specifies the required result for many of these, including the cases where a function of a NaN must return a quiet NaN and where a limit gives a finite value at an infinite argument. Boundary conditions inside the implementation deserve equal attention: the argument where range reduction changes branch, the segment boundaries of a piecewise table, and the threshold between the polynomial and exponential forms of a hyperbolic function are all places where a careless implementation produces a visible discontinuity or loses monotonicity.