1*f99ccf65SEugene Zhulenev //===- PolynomialApproximation.cpp - Approximate math operations ----------===//
2*f99ccf65SEugene Zhulenev //
3*f99ccf65SEugene Zhulenev // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4*f99ccf65SEugene Zhulenev // See https://llvm.org/LICENSE.txt for license information.
5*f99ccf65SEugene Zhulenev // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6*f99ccf65SEugene Zhulenev //
7*f99ccf65SEugene Zhulenev //===----------------------------------------------------------------------===//
8*f99ccf65SEugene Zhulenev //
9*f99ccf65SEugene Zhulenev // This file implements expansion of math operations to fast approximations
10*f99ccf65SEugene Zhulenev // that do not rely on any of the library functions.
11*f99ccf65SEugene Zhulenev //
12*f99ccf65SEugene Zhulenev //===----------------------------------------------------------------------===//
13*f99ccf65SEugene Zhulenev 
14*f99ccf65SEugene Zhulenev #include "mlir/Dialect/Math/IR/Math.h"
15*f99ccf65SEugene Zhulenev #include "mlir/Dialect/Math/Transforms/Passes.h"
16*f99ccf65SEugene Zhulenev #include "mlir/Dialect/Vector/VectorOps.h"
17*f99ccf65SEugene Zhulenev #include "mlir/IR/Builders.h"
18*f99ccf65SEugene Zhulenev #include "mlir/Transforms/DialectConversion.h"
19*f99ccf65SEugene Zhulenev #include "mlir/Transforms/GreedyPatternRewriteDriver.h"
20*f99ccf65SEugene Zhulenev 
21*f99ccf65SEugene Zhulenev using namespace mlir;
22*f99ccf65SEugene Zhulenev using namespace mlir::vector;
23*f99ccf65SEugene Zhulenev 
24*f99ccf65SEugene Zhulenev static bool isValidFloatType(Type type) {
25*f99ccf65SEugene Zhulenev   if (auto vectorType = type.dyn_cast<VectorType>())
26*f99ccf65SEugene Zhulenev     return vectorType.getElementType().isa<FloatType>();
27*f99ccf65SEugene Zhulenev   return type.isa<FloatType>();
28*f99ccf65SEugene Zhulenev }
29*f99ccf65SEugene Zhulenev 
30*f99ccf65SEugene Zhulenev //----------------------------------------------------------------------------//
31*f99ccf65SEugene Zhulenev // A PatternRewriter wrapper that provides concise API for building expansions
32*f99ccf65SEugene Zhulenev // for operations on float scalars or vectors.
33*f99ccf65SEugene Zhulenev //----------------------------------------------------------------------------//
34*f99ccf65SEugene Zhulenev 
35*f99ccf65SEugene Zhulenev namespace {
36*f99ccf65SEugene Zhulenev class FloatApproximationBuilder {
37*f99ccf65SEugene Zhulenev public:
38*f99ccf65SEugene Zhulenev   FloatApproximationBuilder(Location loc, Type type, PatternRewriter &rewriter);
39*f99ccf65SEugene Zhulenev 
40*f99ccf65SEugene Zhulenev   Value constant(double value) const;
41*f99ccf65SEugene Zhulenev 
42*f99ccf65SEugene Zhulenev   Value abs(Value a) const;
43*f99ccf65SEugene Zhulenev   Value min(Value a, Value b) const;
44*f99ccf65SEugene Zhulenev   Value max(Value a, Value b) const;
45*f99ccf65SEugene Zhulenev   Value mul(Value a, Value b) const;
46*f99ccf65SEugene Zhulenev   Value div(Value a, Value b) const;
47*f99ccf65SEugene Zhulenev 
48*f99ccf65SEugene Zhulenev   // Fused multiple-add operation: a * b + c.
49*f99ccf65SEugene Zhulenev   Value madd(Value a, Value b, Value c) const;
50*f99ccf65SEugene Zhulenev 
51*f99ccf65SEugene Zhulenev   // Compares values `a` and `b` with the given `predicate`.
52*f99ccf65SEugene Zhulenev   Value cmp(CmpFPredicate predicate, Value a, Value b) const;
53*f99ccf65SEugene Zhulenev 
54*f99ccf65SEugene Zhulenev   // Selects values from `a` or `b` based on the `predicate`.
55*f99ccf65SEugene Zhulenev   Value select(Value predicate, Value a, Value b) const;
56*f99ccf65SEugene Zhulenev 
57*f99ccf65SEugene Zhulenev private:
58*f99ccf65SEugene Zhulenev   Location loc;
59*f99ccf65SEugene Zhulenev   PatternRewriter &rewriter;
60*f99ccf65SEugene Zhulenev   VectorType vectorType; // can be null for scalar type
61*f99ccf65SEugene Zhulenev   FloatType elementType;
62*f99ccf65SEugene Zhulenev };
63*f99ccf65SEugene Zhulenev } // namespace
64*f99ccf65SEugene Zhulenev 
65*f99ccf65SEugene Zhulenev FloatApproximationBuilder::FloatApproximationBuilder(Location loc, Type type,
66*f99ccf65SEugene Zhulenev                                                      PatternRewriter &rewriter)
67*f99ccf65SEugene Zhulenev     : loc(loc), rewriter(rewriter) {
68*f99ccf65SEugene Zhulenev   vectorType = type.dyn_cast<VectorType>();
69*f99ccf65SEugene Zhulenev 
70*f99ccf65SEugene Zhulenev   if (vectorType)
71*f99ccf65SEugene Zhulenev     elementType = vectorType.getElementType().cast<FloatType>();
72*f99ccf65SEugene Zhulenev   else
73*f99ccf65SEugene Zhulenev     elementType = type.cast<FloatType>();
74*f99ccf65SEugene Zhulenev }
75*f99ccf65SEugene Zhulenev 
76*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::constant(double value) const {
77*f99ccf65SEugene Zhulenev   auto attr = rewriter.getFloatAttr(elementType, value);
78*f99ccf65SEugene Zhulenev   Value scalar = rewriter.create<ConstantOp>(loc, attr);
79*f99ccf65SEugene Zhulenev 
80*f99ccf65SEugene Zhulenev   if (vectorType)
81*f99ccf65SEugene Zhulenev     return rewriter.create<BroadcastOp>(loc, vectorType, scalar);
82*f99ccf65SEugene Zhulenev   return scalar;
83*f99ccf65SEugene Zhulenev }
84*f99ccf65SEugene Zhulenev 
85*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::abs(Value a) const {
86*f99ccf65SEugene Zhulenev   return rewriter.create<AbsFOp>(loc, a);
87*f99ccf65SEugene Zhulenev }
88*f99ccf65SEugene Zhulenev 
89*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::min(Value a, Value b) const {
90*f99ccf65SEugene Zhulenev   return select(cmp(CmpFPredicate::OLT, a, b), a, b);
91*f99ccf65SEugene Zhulenev }
92*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::max(Value a, Value b) const {
93*f99ccf65SEugene Zhulenev   return select(cmp(CmpFPredicate::OGT, a, b), a, b);
94*f99ccf65SEugene Zhulenev }
95*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::mul(Value a, Value b) const {
96*f99ccf65SEugene Zhulenev   return rewriter.create<MulFOp>(loc, a, b);
97*f99ccf65SEugene Zhulenev }
98*f99ccf65SEugene Zhulenev 
99*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::div(Value a, Value b) const {
100*f99ccf65SEugene Zhulenev   return rewriter.create<DivFOp>(loc, a, b);
101*f99ccf65SEugene Zhulenev }
102*f99ccf65SEugene Zhulenev 
103*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::madd(Value a, Value b, Value c) const {
104*f99ccf65SEugene Zhulenev   return rewriter.create<FmaFOp>(loc, a, b, c);
105*f99ccf65SEugene Zhulenev }
106*f99ccf65SEugene Zhulenev 
107*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::cmp(CmpFPredicate predicate, Value a,
108*f99ccf65SEugene Zhulenev                                      Value b) const {
109*f99ccf65SEugene Zhulenev   return rewriter.create<CmpFOp>(loc, predicate, a, b);
110*f99ccf65SEugene Zhulenev }
111*f99ccf65SEugene Zhulenev 
112*f99ccf65SEugene Zhulenev Value FloatApproximationBuilder::select(Value predicate, Value a,
113*f99ccf65SEugene Zhulenev                                         Value b) const {
114*f99ccf65SEugene Zhulenev   return rewriter.create<SelectOp>(loc, predicate, a, b);
115*f99ccf65SEugene Zhulenev }
116*f99ccf65SEugene Zhulenev 
117*f99ccf65SEugene Zhulenev //----------------------------------------------------------------------------//
118*f99ccf65SEugene Zhulenev // TanhOp approximation.
119*f99ccf65SEugene Zhulenev //----------------------------------------------------------------------------//
120*f99ccf65SEugene Zhulenev 
121*f99ccf65SEugene Zhulenev namespace {
122*f99ccf65SEugene Zhulenev struct TanhApproximation : public OpRewritePattern<math::TanhOp> {
123*f99ccf65SEugene Zhulenev public:
124*f99ccf65SEugene Zhulenev   using OpRewritePattern::OpRewritePattern;
125*f99ccf65SEugene Zhulenev 
126*f99ccf65SEugene Zhulenev   LogicalResult matchAndRewrite(math::TanhOp op,
127*f99ccf65SEugene Zhulenev                                 PatternRewriter &rewriter) const final;
128*f99ccf65SEugene Zhulenev };
129*f99ccf65SEugene Zhulenev } // namespace
130*f99ccf65SEugene Zhulenev 
131*f99ccf65SEugene Zhulenev LogicalResult
132*f99ccf65SEugene Zhulenev TanhApproximation::matchAndRewrite(math::TanhOp op,
133*f99ccf65SEugene Zhulenev                                    PatternRewriter &rewriter) const {
134*f99ccf65SEugene Zhulenev   if (!isValidFloatType(op.operand().getType()))
135*f99ccf65SEugene Zhulenev     return rewriter.notifyMatchFailure(op, "unsupported operand type");
136*f99ccf65SEugene Zhulenev 
137*f99ccf65SEugene Zhulenev   Value operand = op.operand();
138*f99ccf65SEugene Zhulenev   FloatApproximationBuilder builder(op->getLoc(), operand.getType(), rewriter);
139*f99ccf65SEugene Zhulenev 
140*f99ccf65SEugene Zhulenev   // Clamp operand into [plusClamp, minusClamp] range.
141*f99ccf65SEugene Zhulenev   Value plusClamp = builder.constant(7.90531110763549805);
142*f99ccf65SEugene Zhulenev   Value minusClamp = builder.constant(-7.9053111076354980);
143*f99ccf65SEugene Zhulenev   Value x = builder.max(builder.min(operand, plusClamp), minusClamp);
144*f99ccf65SEugene Zhulenev 
145*f99ccf65SEugene Zhulenev   // Mask for tiny values that are approximated with `operand`.
146*f99ccf65SEugene Zhulenev   Value tiny = builder.constant(0.0004f);
147*f99ccf65SEugene Zhulenev   Value tinyMask = builder.cmp(CmpFPredicate::OLT, builder.abs(operand), tiny);
148*f99ccf65SEugene Zhulenev 
149*f99ccf65SEugene Zhulenev   // The monomial coefficients of the numerator polynomial (odd).
150*f99ccf65SEugene Zhulenev   Value alpha1 = builder.constant(4.89352455891786e-03);
151*f99ccf65SEugene Zhulenev   Value alpha3 = builder.constant(6.37261928875436e-04);
152*f99ccf65SEugene Zhulenev   Value alpha5 = builder.constant(1.48572235717979e-05);
153*f99ccf65SEugene Zhulenev   Value alpha7 = builder.constant(5.12229709037114e-08);
154*f99ccf65SEugene Zhulenev   Value alpha9 = builder.constant(-8.60467152213735e-11);
155*f99ccf65SEugene Zhulenev   Value alpha11 = builder.constant(2.00018790482477e-13);
156*f99ccf65SEugene Zhulenev   Value alpha13 = builder.constant(-2.76076847742355e-16);
157*f99ccf65SEugene Zhulenev 
158*f99ccf65SEugene Zhulenev   // The monomial coefficients of the denominator polynomial (even).
159*f99ccf65SEugene Zhulenev   Value beta0 = builder.constant(4.89352518554385e-03);
160*f99ccf65SEugene Zhulenev   Value beta2 = builder.constant(2.26843463243900e-03);
161*f99ccf65SEugene Zhulenev   Value beta4 = builder.constant(1.18534705686654e-04);
162*f99ccf65SEugene Zhulenev   Value beta6 = builder.constant(1.19825839466702e-06);
163*f99ccf65SEugene Zhulenev 
164*f99ccf65SEugene Zhulenev   // Since the polynomials are odd/even, we need x^2.
165*f99ccf65SEugene Zhulenev   Value x2 = builder.mul(x, x);
166*f99ccf65SEugene Zhulenev 
167*f99ccf65SEugene Zhulenev   // Evaluate the numerator polynomial p.
168*f99ccf65SEugene Zhulenev   Value p = builder.madd(x2, alpha13, alpha11);
169*f99ccf65SEugene Zhulenev   p = builder.madd(x2, p, alpha9);
170*f99ccf65SEugene Zhulenev   p = builder.madd(x2, p, alpha7);
171*f99ccf65SEugene Zhulenev   p = builder.madd(x2, p, alpha5);
172*f99ccf65SEugene Zhulenev   p = builder.madd(x2, p, alpha3);
173*f99ccf65SEugene Zhulenev   p = builder.madd(x2, p, alpha1);
174*f99ccf65SEugene Zhulenev   p = builder.mul(x, p);
175*f99ccf65SEugene Zhulenev 
176*f99ccf65SEugene Zhulenev   // Evaluate the denominator polynomial q.
177*f99ccf65SEugene Zhulenev   Value q = builder.madd(x2, beta6, beta4);
178*f99ccf65SEugene Zhulenev   q = builder.madd(x2, q, beta2);
179*f99ccf65SEugene Zhulenev   q = builder.madd(x2, q, beta0);
180*f99ccf65SEugene Zhulenev 
181*f99ccf65SEugene Zhulenev   // Divide the numerator by the denominator.
182*f99ccf65SEugene Zhulenev   Value res = builder.select(tinyMask, x, builder.div(p, q));
183*f99ccf65SEugene Zhulenev 
184*f99ccf65SEugene Zhulenev   rewriter.replaceOp(op, res);
185*f99ccf65SEugene Zhulenev 
186*f99ccf65SEugene Zhulenev   return success();
187*f99ccf65SEugene Zhulenev }
188*f99ccf65SEugene Zhulenev 
189*f99ccf65SEugene Zhulenev //----------------------------------------------------------------------------//
190*f99ccf65SEugene Zhulenev 
191*f99ccf65SEugene Zhulenev void mlir::populateMathPolynomialApproximationPatterns(
192*f99ccf65SEugene Zhulenev     OwningRewritePatternList &patterns, MLIRContext *ctx) {
193*f99ccf65SEugene Zhulenev   patterns.insert<TanhApproximation>(ctx);
194*f99ccf65SEugene Zhulenev }
195