numerical approximation of the inverse of normal distribution function original algorithm: https://github.com/scipy/scipy/blob/master/scipy/special/cephes/ndtri.c case 1: 0 < x < exp(-2) z = sqrt(-2 * log(x)) t = 1 / z res = log(z) / z - z + t * P(t) / Q(t) where coefficients of P and Q are different for z < 8 and for z >= 8 case2: exp(-2) <= x <= 1 - exp(-2) w = x - 0.5 res = sqrt(2pi) * (w + w^
| 107 | // ndtri(x) = -ndtri(1 - x) |
| 108 | // fallback to case 1 |
| 109 | mlir::Value ndtri_approx(ValueBuilderHelper& helper, mlir::Value x) { |
| 110 | // polynomial P |
| 111 | auto P = [&](mlir::Value i, mlir::Value cond) { |
| 112 | std::vector<mlir::Value> coeff0 = { |
| 113 | helper.const_f32(4.05544892305962419923E0), |
| 114 | helper.const_f32(3.15251094599893866154E1), |
| 115 | helper.const_f32(5.71628192246421288162E1), |
| 116 | helper.const_f32(4.40805073893200834700E1), |
| 117 | helper.const_f32(1.46849561928858024014E1), |
| 118 | helper.const_f32(2.18663306850790267539E0), |
| 119 | helper.const_f32(-1.40256079171354495875E-1), |
| 120 | helper.const_f32(-3.50424626827848203418E-2), |
| 121 | helper.const_f32(-8.57456785154685413611E-4)}; |
| 122 | std::vector<mlir::Value> coeff1 = { |
| 123 | helper.const_f32(3.23774891776946035970E0), |
| 124 | helper.const_f32(6.91522889068984211695E0), |
| 125 | helper.const_f32(3.93881025292474443415E0), |
| 126 | helper.const_f32(1.33303460815807542389E0), |
| 127 | helper.const_f32(2.01485389549179081538E-1), |
| 128 | helper.const_f32(1.23716634817820021358E-2), |
| 129 | helper.const_f32(3.01581553508235416007E-4), |
| 130 | helper.const_f32(2.65806974686737550832E-6), |
| 131 | helper.const_f32(6.23974539184983293730E-9)}; |
| 132 | return helper.select( |
| 133 | cond, polynomial(helper, i, coeff0), polynomial(helper, i, coeff1)); |
| 134 | }; |
| 135 | |
| 136 | // polynomial Q |
| 137 | auto Q = [&](mlir::Value i, mlir::Value cond) { |
| 138 | std::vector<mlir::Value> coeff0 = { |
| 139 | helper.const_f32(1.f), |
| 140 | helper.const_f32(1.57799883256466749731E1), |
| 141 | helper.const_f32(4.53907635128879210584E1), |
| 142 | helper.const_f32(4.13172038254672030440E1), |
| 143 | helper.const_f32(1.50425385692907503408E1), |
| 144 | helper.const_f32(2.50464946208309415979E0), |
| 145 | helper.const_f32(-1.42182922854787788574E-1), |
| 146 | helper.const_f32(-3.80806407691578277194E-2), |
| 147 | helper.const_f32(-9.33259480895457427372E-4)}; |
| 148 | std::vector<mlir::Value> coeff1 = { |
| 149 | helper.const_f32(1.f), |
| 150 | helper.const_f32(6.02427039364742014255E0), |
| 151 | helper.const_f32(3.67983563856160859403E0), |
| 152 | helper.const_f32(1.37702099489081330271E0), |
| 153 | helper.const_f32(2.16236993594496635890E-1), |
| 154 | helper.const_f32(1.34204006088543189037E-2), |
| 155 | helper.const_f32(3.28014464682127739104E-4), |
| 156 | helper.const_f32(2.89247864745380683936E-6), |
| 157 | helper.const_f32(6.79019408009981274425E-9)}; |
| 158 | return helper.select( |
| 159 | cond, polynomial(helper, i, coeff0), polynomial(helper, i, coeff1)); |
| 160 | }; |
| 161 | |
| 162 | // polynomial R |
| 163 | auto R = [&](mlir::Value i) { |
| 164 | std::vector<mlir::Value> coeff = { |
| 165 | helper.const_f32(-5.99633501014107895267E1), |
| 166 | helper.const_f32(9.80010754185999661536E1), |