Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 1 | // Copyright (c) 1999, Silicon Graphics, Inc. -- ALL RIGHTS RESERVED |
| 2 | // |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 3 | // Permission is granted free of charge to copy, modify, use and distribute |
| 4 | // this software provided you include the entirety of this notice in all |
| 5 | // copies made. |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 6 | // |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 7 | // THIS SOFTWARE IS PROVIDED ON AN AS IS BASIS, WITHOUT WARRANTY OF ANY |
| 8 | // KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, WITHOUT LIMITATION, |
| 9 | // WARRANTIES THAT THE SUBJECT SOFTWARE IS FREE OF DEFECTS, MERCHANTABLE, FIT |
| 10 | // FOR A PARTICULAR PURPOSE OR NON-INFRINGING. SGI ASSUMES NO RISK AS TO THE |
| 11 | // QUALITY AND PERFORMANCE OF THE SOFTWARE. SHOULD THE SOFTWARE PROVE |
| 12 | // DEFECTIVE IN ANY RESPECT, SGI ASSUMES NO COST OR LIABILITY FOR ANY |
| 13 | // SERVICING, REPAIR OR CORRECTION. THIS DISCLAIMER OF WARRANTY CONSTITUTES |
| 14 | // AN ESSENTIAL PART OF THIS LICENSE. NO USE OF ANY SUBJECT SOFTWARE IS |
| 15 | // AUTHORIZED HEREUNDER EXCEPT UNDER THIS DISCLAIMER. |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 16 | // |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 17 | // UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, WHETHER TORT (INCLUDING, |
| 18 | // WITHOUT LIMITATION, NEGLIGENCE OR STRICT LIABILITY), CONTRACT, OR |
| 19 | // OTHERWISE, SHALL SGI BE LIABLE FOR ANY DIRECT, INDIRECT, SPECIAL, |
| 20 | // INCIDENTAL, OR CONSEQUENTIAL DAMAGES OF ANY CHARACTER WITH RESPECT TO THE |
| 21 | // SOFTWARE INCLUDING, WITHOUT LIMITATION, DAMAGES FOR LOSS OF GOODWILL, WORK |
| 22 | // STOPPAGE, LOSS OF DATA, COMPUTER FAILURE OR MALFUNCTION, OR ANY AND ALL |
| 23 | // OTHER COMMERCIAL DAMAGES OR LOSSES, EVEN IF SGI SHALL HAVE BEEN INFORMED OF |
| 24 | // THE POSSIBILITY OF SUCH DAMAGES. THIS LIMITATION OF LIABILITY SHALL NOT |
| 25 | // APPLY TO LIABILITY RESULTING FROM SGI's NEGLIGENCE TO THE EXTENT APPLICABLE |
| 26 | // LAW PROHIBITS SUCH LIMITATION. SOME JURISDICTIONS DO NOT ALLOW THE |
| 27 | // EXCLUSION OR LIMITATION OF INCIDENTAL OR CONSEQUENTIAL DAMAGES, SO THAT |
| 28 | // EXCLUSION AND LIMITATION MAY NOT APPLY TO YOU. |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 29 | // |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 30 | // These license terms shall be governed by and construed in accordance with |
| 31 | // the laws of the United States and the State of California as applied to |
| 32 | // agreements entered into and to be performed entirely within California |
| 33 | // between California residents. Any litigation relating to these license |
| 34 | // terms shall be subject to the exclusive jurisdiction of the Federal Courts |
| 35 | // of the Northern District of California (or, absent subject matter |
| 36 | // jurisdiction in such courts, the courts of the State of California), with |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 37 | // venue lying exclusively in Santa Clara County, California. |
Hans Boehm | 7b2383f | 2014-11-26 18:07:25 -0800 | [diff] [blame] | 38 | // |
Hans Boehm | b849a8e | 2015-05-20 17:53:27 -0700 | [diff] [blame] | 39 | // 5/2014 Added Strings to ArithmeticExceptions. |
| 40 | // 5/2015 Added support for direct asin() implementation in CR. |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 41 | |
Hans Boehm | 97767fa | 2014-11-26 17:16:40 -0800 | [diff] [blame] | 42 | package com.hp.creals; |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 43 | // import android.util.Log; |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 44 | |
| 45 | import java.math.BigInteger; |
| 46 | |
| 47 | /** |
| 48 | * Unary functions on constructive reals implemented as objects. |
| 49 | * The <TT>execute</tt> member computes the function result. |
| 50 | * Unary function objects on constructive reals inherit from |
| 51 | * <TT>UnaryCRFunction</tt>. |
| 52 | */ |
| 53 | // Naming vaguely follows ObjectSpace JGL convention. |
| 54 | public abstract class UnaryCRFunction { |
| 55 | abstract public CR execute(CR x); |
| 56 | |
| 57 | /** |
| 58 | * The function object corresponding to the identity function. |
| 59 | */ |
| 60 | public static final UnaryCRFunction identityFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 61 | new identity_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 62 | |
| 63 | /** |
| 64 | * The function object corresponding to the <TT>negate</tt> method of CR. |
| 65 | */ |
| 66 | public static final UnaryCRFunction negateFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 67 | new negate_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 68 | |
| 69 | /** |
| 70 | * The function object corresponding to the <TT>inverse</tt> method of CR. |
| 71 | */ |
| 72 | public static final UnaryCRFunction inverseFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 73 | new inverse_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 74 | |
| 75 | /** |
| 76 | * The function object corresponding to the <TT>abs</tt> method of CR. |
| 77 | */ |
| 78 | public static final UnaryCRFunction absFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 79 | new abs_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 80 | |
| 81 | /** |
| 82 | * The function object corresponding to the <TT>exp</tt> method of CR. |
| 83 | */ |
| 84 | public static final UnaryCRFunction expFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 85 | new exp_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 86 | |
| 87 | /** |
| 88 | * The function object corresponding to the <TT>cos</tt> method of CR. |
| 89 | */ |
| 90 | public static final UnaryCRFunction cosFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 91 | new cos_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 92 | |
| 93 | /** |
| 94 | * The function object corresponding to the <TT>sin</tt> method of CR. |
| 95 | */ |
| 96 | public static final UnaryCRFunction sinFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 97 | new sin_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 98 | |
| 99 | /** |
| 100 | * The function object corresponding to the tangent function. |
| 101 | */ |
| 102 | public static final UnaryCRFunction tanFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 103 | new tan_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 104 | |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 105 | /** |
| 106 | * The function object corresponding to the inverse sine (arcsine) function. |
| 107 | * The argument must be between -1 and 1 inclusive. The result is between |
| 108 | * -PI/2 and PI/2. |
| 109 | */ |
| 110 | public static final UnaryCRFunction asinFunction = |
Hans Boehm | b849a8e | 2015-05-20 17:53:27 -0700 | [diff] [blame] | 111 | new asin_UnaryCRFunction(); |
| 112 | // The following also works, but is slower: |
| 113 | // CR half_pi = CR.PI.divide(CR.valueOf(2)); |
| 114 | // UnaryCRFunction.sinFunction.inverseMonotone(half_pi.negate(), |
| 115 | // half_pi); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 116 | |
| 117 | /** |
| 118 | * The function object corresponding to the inverse cosine (arccosine) function. |
| 119 | * The argument must be between -1 and 1 inclusive. The result is between |
| 120 | * 0 and PI. |
| 121 | */ |
| 122 | public static final UnaryCRFunction acosFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 123 | new acos_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 124 | |
| 125 | /** |
| 126 | * The function object corresponding to the inverse cosine (arctangent) function. |
| 127 | * The result is between -PI/2 and PI/2. |
| 128 | */ |
| 129 | public static final UnaryCRFunction atanFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 130 | new atan_UnaryCRFunction(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 131 | |
| 132 | /** |
| 133 | * The function object corresponding to the <TT>ln</tt> method of CR. |
| 134 | */ |
| 135 | public static final UnaryCRFunction lnFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 136 | new ln_UnaryCRFunction(); |
| 137 | |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 138 | /** |
| 139 | * The function object corresponding to the <TT>sqrt</tt> method of CR. |
| 140 | */ |
| 141 | public static final UnaryCRFunction sqrtFunction = |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 142 | new sqrt_UnaryCRFunction(); |
| 143 | |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 144 | /** |
| 145 | * Compose this function with <TT>f2</tt>. |
| 146 | */ |
| 147 | public UnaryCRFunction compose(UnaryCRFunction f2) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 148 | return new compose_UnaryCRFunction(this, f2); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 149 | } |
| 150 | |
| 151 | /** |
| 152 | * Compute the inverse of this function, which must be defined |
| 153 | * and strictly monotone on the interval [<TT>low</tt>, <TT>high</tt>]. |
| 154 | * The resulting function is defined only on the image of |
| 155 | * [<TT>low</tt>, <TT>high</tt>]. |
| 156 | * The original function may be either increasing or decreasing. |
| 157 | */ |
| 158 | public UnaryCRFunction inverseMonotone(CR low, CR high) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 159 | return new inverseMonotone_UnaryCRFunction(this, low, high); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 160 | } |
| 161 | |
| 162 | /** |
| 163 | * Compute the derivative of a function. |
| 164 | * The function must be defined on the interval [<TT>low</tt>, <TT>high</tt>], |
| 165 | * and the derivative must exist, and must be continuous and |
| 166 | * monotone in the open interval [<TT>low</tt>, <TT>high</tt>]. |
| 167 | * The result is defined only in the open interval. |
| 168 | */ |
| 169 | public UnaryCRFunction monotoneDerivative(CR low, CR high) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 170 | return new monotoneDerivative_UnaryCRFunction(this, low, high); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 171 | } |
| 172 | |
| 173 | } |
| 174 | |
| 175 | // Subclasses of UnaryCRFunction for various built-in functions. |
| 176 | class sin_UnaryCRFunction extends UnaryCRFunction { |
| 177 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 178 | return x.sin(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 179 | } |
| 180 | } |
| 181 | |
| 182 | class cos_UnaryCRFunction extends UnaryCRFunction { |
| 183 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 184 | return x.cos(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 185 | } |
| 186 | } |
| 187 | |
| 188 | class tan_UnaryCRFunction extends UnaryCRFunction { |
| 189 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 190 | return x.sin().divide(x.cos()); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 191 | } |
| 192 | } |
| 193 | |
Hans Boehm | b849a8e | 2015-05-20 17:53:27 -0700 | [diff] [blame] | 194 | class asin_UnaryCRFunction extends UnaryCRFunction { |
| 195 | public CR execute(CR x) { |
| 196 | return x.asin(); |
| 197 | } |
| 198 | } |
| 199 | |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 200 | class acos_UnaryCRFunction extends UnaryCRFunction { |
| 201 | public CR execute(CR x) { |
Hans Boehm | b849a8e | 2015-05-20 17:53:27 -0700 | [diff] [blame] | 202 | return x.acos(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 203 | } |
| 204 | } |
| 205 | |
| 206 | // This uses the identity (sin x)^2 = (tan x)^2/(1 + (tan x)^2) |
| 207 | // Since we know the tangent of the result, we can get its sine, |
| 208 | // and then use the asin function. Note that we don't always |
| 209 | // want the positive square root when computing the sine. |
| 210 | class atan_UnaryCRFunction extends UnaryCRFunction { |
| 211 | CR one = CR.valueOf(1); |
| 212 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 213 | CR x2 = x.multiply(x); |
| 214 | CR abs_sin_atan = x2.divide(one.add(x2)).sqrt(); |
| 215 | CR sin_atan = x.select(abs_sin_atan.negate(), abs_sin_atan); |
Hans Boehm | b849a8e | 2015-05-20 17:53:27 -0700 | [diff] [blame] | 216 | return sin_atan.asin(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 217 | } |
| 218 | } |
| 219 | |
| 220 | class exp_UnaryCRFunction extends UnaryCRFunction { |
| 221 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 222 | return x.exp(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 223 | } |
| 224 | } |
| 225 | |
| 226 | class ln_UnaryCRFunction extends UnaryCRFunction { |
| 227 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 228 | return x.ln(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 229 | } |
| 230 | } |
| 231 | |
| 232 | class identity_UnaryCRFunction extends UnaryCRFunction { |
| 233 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 234 | return x; |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 235 | } |
| 236 | } |
| 237 | |
| 238 | class negate_UnaryCRFunction extends UnaryCRFunction { |
| 239 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 240 | return x.negate(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 241 | } |
| 242 | } |
| 243 | |
| 244 | class inverse_UnaryCRFunction extends UnaryCRFunction { |
| 245 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 246 | return x.inverse(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 247 | } |
| 248 | } |
| 249 | |
| 250 | class abs_UnaryCRFunction extends UnaryCRFunction { |
| 251 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 252 | return x.abs(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 253 | } |
| 254 | } |
| 255 | |
| 256 | class sqrt_UnaryCRFunction extends UnaryCRFunction { |
| 257 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 258 | return x.sqrt(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 259 | } |
| 260 | } |
| 261 | |
| 262 | class compose_UnaryCRFunction extends UnaryCRFunction { |
| 263 | UnaryCRFunction f1; |
| 264 | UnaryCRFunction f2; |
| 265 | compose_UnaryCRFunction(UnaryCRFunction func1, |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 266 | UnaryCRFunction func2) { |
| 267 | f1 = func1; f2 = func2; |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 268 | } |
| 269 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 270 | return f1.execute(f2.execute(x)); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 271 | } |
| 272 | } |
| 273 | |
| 274 | class inverseMonotone_UnaryCRFunction extends UnaryCRFunction { |
| 275 | // The following variables are final, so that they |
| 276 | // can be referenced from the inner class inverseIncreasingCR. |
| 277 | // I couldn't find a way to initialize these such that the |
| 278 | // compiler accepted them as final without turning them into arrays. |
| 279 | final UnaryCRFunction f[] = new UnaryCRFunction[1]; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 280 | // Monotone increasing. |
| 281 | // If it was monotone decreasing, we |
| 282 | // negate it. |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 283 | final boolean f_negated[] = new boolean[1]; |
| 284 | final CR low[] = new CR[1]; |
| 285 | final CR high[] = new CR[1]; |
| 286 | final CR f_low[] = new CR[1]; |
| 287 | final CR f_high[] = new CR[1]; |
| 288 | final int max_msd[] = new int[1]; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 289 | // Bound on msd of both f(high) and f(low) |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 290 | final int max_arg_prec[] = new int[1]; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 291 | // base**max_arg_prec is a small fraction |
| 292 | // of low - high. |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 293 | final int deriv_msd[] = new int[1]; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 294 | // Rough approx. of msd of first |
| 295 | // derivative. |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 296 | final static BigInteger BIG1023 = BigInteger.valueOf(1023); |
| 297 | static final boolean ENABLE_TRACE = false; // Change to generate trace |
| 298 | static void trace(String s) { |
| 299 | if (ENABLE_TRACE) { |
| 300 | System.out.println(s); |
| 301 | // Change to Log.v("UnaryCRFunction", s); for Android use. |
| 302 | } |
| 303 | } |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 304 | inverseMonotone_UnaryCRFunction(UnaryCRFunction func, CR l, CR h) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 305 | low[0] = l; high[0] = h; |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 306 | CR tmp_f_low = func.execute(l); |
| 307 | CR tmp_f_high = func.execute(h); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 308 | // Since func is monotone and low < high, the following test |
| 309 | // converges. |
| 310 | if (tmp_f_low.compareTo(tmp_f_high) > 0) { |
| 311 | f[0] = UnaryCRFunction.negateFunction.compose(func); |
| 312 | f_negated[0] = true; |
| 313 | f_low[0] = tmp_f_low.negate(); |
| 314 | f_high[0] = tmp_f_high.negate(); |
| 315 | } else { |
| 316 | f[0] = func; |
| 317 | f_negated[0] = false; |
| 318 | f_low[0] = tmp_f_low; |
| 319 | f_high[0] = tmp_f_high; |
| 320 | } |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 321 | max_msd[0] = low[0].abs().max(high[0].abs()).msd(); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 322 | max_arg_prec[0] = high[0].subtract(low[0]).msd() - 4; |
| 323 | deriv_msd[0] = f_high[0].subtract(f_low[0]) |
| 324 | .divide(high[0].subtract(low[0])).msd(); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 325 | } |
| 326 | class inverseIncreasingCR extends CR { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 327 | final CR arg; |
| 328 | inverseIncreasingCR(CR x) { |
| 329 | arg = f_negated[0]? x.negate() : x; |
| 330 | } |
| 331 | // Comparison with a difference of one treated as equality. |
| 332 | int sloppy_compare(BigInteger x, BigInteger y) { |
| 333 | BigInteger difference = x.subtract(y); |
| 334 | if (difference.compareTo(big1) > 0) { |
| 335 | return 1; |
| 336 | } |
| 337 | if (difference.compareTo(bigm1) < 0) { |
| 338 | return -1; |
| 339 | } |
| 340 | return 0; |
| 341 | } |
| 342 | protected BigInteger approximate(int p) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 343 | final int extra_arg_prec = 4; |
| 344 | final UnaryCRFunction fn = f[0]; |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 345 | int small_step_deficit = 0; // Number of ineffective steps not |
| 346 | // yet compensated for by a binary |
| 347 | // search step. |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 348 | int digits_needed = max_msd[0] - p; |
| 349 | if (digits_needed < 0) return big0; |
| 350 | int working_arg_prec = p - extra_arg_prec; |
| 351 | if (working_arg_prec > max_arg_prec[0]) { |
| 352 | working_arg_prec = max_arg_prec[0]; |
| 353 | } |
| 354 | int working_eval_prec = working_arg_prec + deriv_msd[0] - 20; |
| 355 | // initial guess |
| 356 | // We use a combination of binary search and something like |
| 357 | // the secant method. This always converges linearly, |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 358 | // and should converge quadratically under favorable assumptions. |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 359 | // F_l and f_h are always the approximate images of l and h. |
| 360 | // At any point, arg is between f_l and f_h, or no more than |
| 361 | // one outside [f_l, f_h]. |
| 362 | // L and h are implicitly scaled by working_arg_prec. |
| 363 | // The scaled values of l and h are strictly between low and high. |
| 364 | // If at_left is true, then l is logically at the left |
| 365 | // end of the interval. We approximate this by setting l to |
| 366 | // a point slightly inside the interval, and letting f_l |
| 367 | // approximate the function value at the endpoint. |
| 368 | // If at_right is true, r and f_r are set correspondingly. |
| 369 | // At the endpoints of the interval, f_l and f_h may correspond |
| 370 | // to the endpoints, even if l and h are slightly inside. |
| 371 | // F_l and f_u are scaled by working_eval_prec. |
| 372 | // Working_eval_prec may need to be adjusted depending |
| 373 | // on the derivative of f. |
| 374 | boolean at_left, at_right; |
| 375 | BigInteger l, f_l; |
| 376 | BigInteger h, f_h; |
| 377 | BigInteger low_appr = low[0].get_appr(working_arg_prec) |
| 378 | .add(big1); |
| 379 | BigInteger high_appr = high[0].get_appr(working_arg_prec) |
| 380 | .subtract(big1); |
| 381 | BigInteger arg_appr = arg.get_appr(working_eval_prec); |
| 382 | boolean have_good_appr = (appr_valid && min_prec < max_msd[0]); |
| 383 | if (digits_needed < 30 && !have_good_appr) { |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 384 | trace("Setting interval to entire domain"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 385 | h = high_appr; |
| 386 | f_h = f_high[0].get_appr(working_eval_prec); |
| 387 | l = low_appr; |
| 388 | f_l = f_low[0].get_appr(working_eval_prec); |
| 389 | // Check for clear out-of-bounds case. |
| 390 | // Close cases may fail in other ways. |
| 391 | if (f_h.compareTo(arg_appr.subtract(big1)) < 0 |
| 392 | || f_l.compareTo(arg_appr.add(big1)) > 0) { |
Hans Boehm | 7b2383f | 2014-11-26 18:07:25 -0800 | [diff] [blame] | 393 | throw new ArithmeticException("inverse(out-of-bounds)"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 394 | } |
| 395 | at_left = true; |
| 396 | at_right = true; |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 397 | small_step_deficit = 2; // Start with bin search steps. |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 398 | } else { |
| 399 | int rough_prec = p + digits_needed/2; |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 400 | |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 401 | if (have_good_appr && |
| 402 | (digits_needed < 30 || min_prec < p + 3*digits_needed/4)) { |
| 403 | rough_prec = min_prec; |
| 404 | } |
| 405 | BigInteger rough_appr = get_appr(rough_prec); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 406 | trace("Setting interval based on prev. appr"); |
| 407 | trace("prev. prec = " + rough_prec + " appr = " + rough_appr); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 408 | h = rough_appr.add(big1) |
| 409 | .shiftLeft(rough_prec - working_arg_prec); |
| 410 | l = rough_appr.subtract(big1) |
| 411 | .shiftLeft(rough_prec - working_arg_prec); |
| 412 | if (h.compareTo(high_appr) > 0) { |
| 413 | h = high_appr; |
| 414 | f_h = f_high[0].get_appr(working_eval_prec); |
| 415 | at_right = true; |
| 416 | } else { |
| 417 | CR h_cr = CR.valueOf(h).shiftLeft(working_arg_prec); |
| 418 | f_h = fn.execute(h_cr).get_appr(working_eval_prec); |
| 419 | at_right = false; |
| 420 | } |
| 421 | if (l.compareTo(low_appr) < 0) { |
| 422 | l = low_appr; |
| 423 | f_l = f_low[0].get_appr(working_eval_prec); |
| 424 | at_left = true; |
| 425 | } else { |
| 426 | CR l_cr = CR.valueOf(l).shiftLeft(working_arg_prec); |
| 427 | f_l = fn.execute(l_cr).get_appr(working_eval_prec); |
| 428 | at_left = false; |
| 429 | } |
| 430 | } |
| 431 | BigInteger difference = h.subtract(l); |
| 432 | for(int i = 0;; ++i) { |
| 433 | if (Thread.interrupted() || please_stop) |
Hans Boehm | 9666c57 | 2015-06-19 18:27:45 -0700 | [diff] [blame] | 434 | throw new AbortedException(); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 435 | trace("***Iteration: " + i); |
| 436 | trace("Arg prec = " + working_arg_prec |
| 437 | + " eval prec = " + working_eval_prec |
| 438 | + " arg appr. = " + arg_appr); |
| 439 | trace("l = " + l); trace("h = " + h); |
| 440 | trace("f(l) = " + f_l); trace("f(h) = " + f_h); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 441 | if (difference.compareTo(big6) < 0) { |
| 442 | // Answer is less than 1/2 ulp away from h. |
| 443 | return scale(h, -extra_arg_prec); |
| 444 | } |
| 445 | BigInteger f_difference = f_h.subtract(f_l); |
| 446 | // Narrow the interval by dividing at a cleverly |
| 447 | // chosen point (guess) in the middle. |
| 448 | { |
| 449 | BigInteger guess; |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 450 | boolean binary_step = |
| 451 | (small_step_deficit > 0 || f_difference.signum() == 0); |
| 452 | if (binary_step) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 453 | // Do a binary search step to guarantee linear |
| 454 | // convergence. |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 455 | trace("binary step"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 456 | guess = l.add(h).shiftRight(1); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 457 | --small_step_deficit; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 458 | } else { |
| 459 | // interpolate. |
| 460 | // f_difference is nonzero here. |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 461 | trace("interpolating"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 462 | BigInteger arg_difference = arg_appr.subtract(f_l); |
| 463 | BigInteger t = arg_difference.multiply(difference); |
| 464 | BigInteger adj = t.divide(f_difference); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 465 | // tentative adjustment to l to compute guess |
| 466 | // If we are within 1/1024 of either end, back off. |
| 467 | // This greatly improves the odds of bounding |
| 468 | // the answer within the smaller interval. |
| 469 | // Note that interpolation will often get us |
| 470 | // MUCH closer than this. |
| 471 | if (adj.compareTo(difference.shiftRight(10)) < 0) { |
| 472 | adj = adj.shiftLeft(8); |
| 473 | trace("adjusting left"); |
| 474 | } else if (adj.compareTo(difference.multiply(BIG1023) |
| 475 | .shiftRight(10)) > 0){ |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 476 | adj = difference.subtract(difference.subtract(adj) |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 477 | .shiftLeft(8)); |
| 478 | trace("adjusting right"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 479 | } |
| 480 | if (adj.signum() <= 0) |
| 481 | adj = big2; |
| 482 | if (adj.compareTo(difference) >= 0) |
| 483 | adj = difference.subtract(big2); |
| 484 | guess = (adj.signum() <= 0? l.add(big2) : l.add(adj)); |
| 485 | } |
| 486 | int outcome; |
| 487 | BigInteger tweak = big2; |
| 488 | BigInteger f_guess; |
| 489 | for(boolean adj_prec = false;; adj_prec = !adj_prec) { |
| 490 | CR guess_cr = CR.valueOf(guess) |
| 491 | .shiftLeft(working_arg_prec); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 492 | trace("Evaluating at " + guess_cr |
| 493 | + " with precision " + working_eval_prec); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 494 | CR f_guess_cr = fn.execute(guess_cr); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 495 | trace("fn value = " + f_guess_cr); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 496 | f_guess = f_guess_cr.get_appr(working_eval_prec); |
| 497 | outcome = sloppy_compare(f_guess, arg_appr); |
| 498 | if (outcome != 0) break; |
| 499 | // Alternately increase evaluation precision |
| 500 | // and adjust guess slightly. |
| 501 | // This should be an unlikely case. |
| 502 | if (adj_prec) { |
| 503 | // adjust working_eval_prec to get enough |
| 504 | // resolution. |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 505 | int adjustment = -f_guess.bitLength()/4; |
| 506 | if (adjustment > -20) adjustment = - 20; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 507 | CR l_cr = CR.valueOf(l) |
| 508 | .shiftLeft(working_arg_prec); |
| 509 | CR h_cr = CR.valueOf(h) |
| 510 | .shiftLeft(working_arg_prec); |
| 511 | working_eval_prec += adjustment; |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 512 | trace("New eval prec = " + working_eval_prec |
| 513 | + (at_left? "(at left)" : "") |
| 514 | + (at_right? "(at right)" : "")); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 515 | if (at_left) { |
| 516 | f_l = f_low[0].get_appr(working_eval_prec); |
| 517 | } else { |
| 518 | f_l = fn.execute(l_cr) |
| 519 | .get_appr(working_eval_prec); |
| 520 | } |
| 521 | if (at_right) { |
| 522 | f_h = f_high[0].get_appr(working_eval_prec); |
| 523 | } else { |
| 524 | f_h = fn.execute(h_cr) |
| 525 | .get_appr(working_eval_prec); |
| 526 | } |
| 527 | arg_appr = arg.get_appr(working_eval_prec); |
| 528 | } else { |
| 529 | // guess might be exactly right; tweak it |
| 530 | // slightly. |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 531 | trace("tweaking guess"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 532 | BigInteger new_guess = guess.add(tweak); |
| 533 | if (new_guess.compareTo(h) >= 0) { |
| 534 | guess = guess.subtract(tweak); |
| 535 | } else { |
| 536 | guess = new_guess; |
| 537 | } |
| 538 | // If we keep hitting the right answer, it's |
| 539 | // important to alternate which side we move it |
| 540 | // to, so that the interval shrinks rapidly. |
| 541 | tweak = tweak.negate(); |
| 542 | } |
| 543 | } |
| 544 | if (outcome > 0) { |
| 545 | h = guess; |
| 546 | f_h = f_guess; |
| 547 | at_right = false; |
| 548 | } else { |
| 549 | l = guess; |
| 550 | f_l = f_guess; |
| 551 | at_left = false; |
| 552 | } |
| 553 | BigInteger new_difference = h.subtract(l); |
Hans Boehm | a8e45a3 | 2015-02-25 14:53:13 -0800 | [diff] [blame] | 554 | if (!binary_step) { |
| 555 | if (new_difference.compareTo(difference |
| 556 | .shiftRight(1)) >= 0) { |
| 557 | ++small_step_deficit; |
| 558 | } else { |
| 559 | --small_step_deficit; |
| 560 | } |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 561 | } |
| 562 | difference = new_difference; |
| 563 | } |
| 564 | } |
| 565 | } |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 566 | } |
| 567 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 568 | return new inverseIncreasingCR(x); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 569 | } |
| 570 | } |
| 571 | |
| 572 | class monotoneDerivative_UnaryCRFunction extends UnaryCRFunction { |
| 573 | // The following variables are final, so that they |
| 574 | // can be referenced from the inner class inverseIncreasingCR. |
| 575 | final UnaryCRFunction f[] = new UnaryCRFunction[1]; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 576 | // Monotone increasing. |
| 577 | // If it was monotone decreasing, we |
| 578 | // negate it. |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 579 | final CR low[] = new CR[1]; // endpoints and mispoint of interval |
| 580 | final CR mid[] = new CR[1]; |
| 581 | final CR high[] = new CR[1]; |
| 582 | final CR f_low[] = new CR[1]; // Corresponding function values. |
| 583 | final CR f_mid[] = new CR[1]; |
| 584 | final CR f_high[] = new CR[1]; |
| 585 | final int difference_msd[] = new int[1]; // msd of interval len. |
| 586 | final int deriv2_msd[] = new int[1]; |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 587 | // Rough approx. of msd of second |
| 588 | // derivative. |
| 589 | // This is increased to be an appr. bound |
| 590 | // on the msd of |(f'(y)-f'(x))/(x-y)| |
| 591 | // for any pair of points x and y |
| 592 | // we have considered. |
| 593 | // It may be better to keep a copy per |
| 594 | // derivative value. |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 595 | |
| 596 | monotoneDerivative_UnaryCRFunction(UnaryCRFunction func, CR l, CR h) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 597 | f[0] = func; |
| 598 | low[0] = l; high[0] = h; |
| 599 | mid[0] = l.add(h).shiftRight(1); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 600 | f_low[0] = func.execute(l); |
| 601 | f_mid[0] = func.execute(mid[0]); |
| 602 | f_high[0] = func.execute(h); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 603 | CR difference = h.subtract(l); |
| 604 | // compute approximate msd of |
| 605 | // ((f_high - f_mid) - (f_mid - f_low))/(high - low) |
| 606 | // This should be a very rough appr to the second derivative. |
| 607 | // We add a little slop to err on the high side, since |
| 608 | // a low estimate will cause extra iterations. |
| 609 | CR appr_diff2 = f_high[0].subtract(f_mid[0].shiftLeft(1)).add(f_low[0]); |
| 610 | difference_msd[0] = difference.msd(); |
| 611 | deriv2_msd[0] = appr_diff2.msd() - difference_msd[0] + 4; |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 612 | } |
| 613 | class monotoneDerivativeCR extends CR { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 614 | CR arg; |
| 615 | CR f_arg; |
| 616 | int max_delta_msd; |
| 617 | monotoneDerivativeCR(CR x) { |
| 618 | arg = x; |
| 619 | f_arg = f[0].execute(x); |
| 620 | // The following must converge, since arg must be in the |
| 621 | // open interval. |
| 622 | CR left_diff = arg.subtract(low[0]); |
| 623 | int max_delta_left_msd = left_diff.msd(); |
| 624 | CR right_diff = high[0].subtract(arg); |
| 625 | int max_delta_right_msd = right_diff.msd(); |
| 626 | if (left_diff.signum() < 0 || right_diff.signum() < 0) { |
Hans Boehm | 7b2383f | 2014-11-26 18:07:25 -0800 | [diff] [blame] | 627 | throw new ArithmeticException("fn not monotone"); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 628 | } |
| 629 | max_delta_msd = (max_delta_left_msd < max_delta_right_msd? |
| 630 | max_delta_left_msd |
| 631 | : max_delta_right_msd); |
| 632 | } |
| 633 | protected BigInteger approximate(int p) { |
| 634 | final int extra_prec = 4; |
| 635 | int log_delta = p - deriv2_msd[0]; |
| 636 | // Ensure that we stay within the interval. |
| 637 | if (log_delta > max_delta_msd) log_delta = max_delta_msd; |
| 638 | log_delta -= extra_prec; |
Hans Boehm | b849a8e | 2015-05-20 17:53:27 -0700 | [diff] [blame] | 639 | CR delta = ONE.shiftLeft(log_delta); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 640 | |
| 641 | CR left = arg.subtract(delta); |
| 642 | CR right = arg.add(delta); |
| 643 | CR f_left = f[0].execute(left); |
| 644 | CR f_right = f[0].execute(right); |
| 645 | CR left_deriv = f_arg.subtract(f_left).shiftRight(log_delta); |
| 646 | CR right_deriv = f_right.subtract(f_arg).shiftRight(log_delta); |
| 647 | int eval_prec = p - extra_prec; |
| 648 | BigInteger appr_left_deriv = left_deriv.get_appr(eval_prec); |
| 649 | BigInteger appr_right_deriv = right_deriv.get_appr(eval_prec); |
| 650 | BigInteger deriv_difference = |
| 651 | appr_right_deriv.subtract(appr_left_deriv).abs(); |
| 652 | if (deriv_difference.compareTo(big8) < 0) { |
| 653 | return scale(appr_left_deriv, -extra_prec); |
| 654 | } else { |
| 655 | if (Thread.interrupted() || please_stop) |
Hans Boehm | 9666c57 | 2015-06-19 18:27:45 -0700 | [diff] [blame] | 656 | throw new AbortedException(); |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 657 | deriv2_msd[0] = |
| 658 | eval_prec + deriv_difference.bitLength() + 4/*slop*/; |
| 659 | deriv2_msd[0] -= log_delta; |
| 660 | return approximate(p); |
| 661 | } |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 662 | } |
| 663 | } |
| 664 | public CR execute(CR x) { |
Hans Boehm | 36dc537 | 2014-11-26 13:35:59 -0800 | [diff] [blame] | 665 | return new monotoneDerivativeCR(x); |
Hans Boehm | 349dbd7 | 2014-11-25 15:28:42 -0800 | [diff] [blame] | 666 | } |
| 667 | } |