001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017package org.apache.commons.math4.legacy.analysis.differentiation; 018 019import java.util.Arrays; 020 021import org.apache.commons.numbers.core.Sum; 022import org.apache.commons.math4.legacy.core.Field; 023import org.apache.commons.math4.legacy.core.FieldElement; 024import org.apache.commons.math4.legacy.core.RealFieldElement; 025import org.apache.commons.math4.legacy.exception.DimensionMismatchException; 026import org.apache.commons.math4.legacy.exception.NumberIsTooLargeException; 027import org.apache.commons.math4.core.jdkmath.JdkMath; 028import org.apache.commons.math4.legacy.core.MathArrays; 029 030/** Class representing both the value and the differentials of a function. 031 * <p>This class is the workhorse of the differentiation package.</p> 032 * <p>This class is an implementation of the extension to Rall's 033 * numbers described in Dan Kalman's paper <a 034 * href="http://www1.american.edu/cas/mathstat/People/kalman/pdffiles/mmgautodiff.pdf">Doubly 035 * Recursive Multivariate Automatic Differentiation</a>, Mathematics Magazine, vol. 75, 036 * no. 3, June 2002. Rall's numbers are an extension to the real numbers used 037 * throughout mathematical expressions; they hold the derivative together with the 038 * value of a function. Dan Kalman's derivative structures hold all partial derivatives 039 * up to any specified order, with respect to any number of free parameters. Rall's 040 * numbers therefore can be seen as derivative structures for order one derivative and 041 * one free parameter, and real numbers can be seen as derivative structures with zero 042 * order derivative and no free parameters.</p> 043 * <p>{@link DerivativeStructure} instances can be used directly thanks to 044 * the arithmetic operators to the mathematical functions provided as 045 * methods by this class (+, -, *, /, %, sin, cos ...).</p> 046 * <p>Implementing complex expressions by hand using these classes is 047 * a tedious and error-prone task but has the advantage of having no limitation 048 * on the derivation order despite not requiring users to compute the derivatives by 049 * themselves. Implementing complex expression can also be done by developing computation 050 * code using standard primitive double values and to use {@link 051 * UnivariateFunctionDifferentiator differentiators} to create the {@link 052 * DerivativeStructure}-based instances. This method is simpler but may be limited in 053 * the accuracy and derivation orders and may be computationally intensive (this is 054 * typically the case for {@link FiniteDifferencesDifferentiator finite differences 055 * differentiator}.</p> 056 * <p>Instances of this class are guaranteed to be immutable.</p> 057 * @see DSCompiler 058 * @since 3.1 059 */ 060public class DerivativeStructure implements RealFieldElement<DerivativeStructure> { 061 /** Compiler for the current dimensions. */ 062 private DSCompiler compiler; 063 064 /** Combined array holding all values. */ 065 private final double[] data; 066 067 /** Build an instance with all values and derivatives set to 0. 068 * @param compiler compiler to use for computation 069 */ 070 private DerivativeStructure(final DSCompiler compiler) { 071 this.compiler = compiler; 072 this.data = new double[compiler.getSize()]; 073 } 074 075 /** Build an instance with all values and derivatives set to 0. 076 * @param parameters number of free parameters 077 * @param order derivation order 078 * @throws NumberIsTooLargeException if order is too large. 079 */ 080 public DerivativeStructure(final int parameters, final int order) { 081 this(DSCompiler.getCompiler(parameters, order)); 082 } 083 084 /** Build an instance representing a constant value. 085 * @param parameters number of free parameters 086 * @param order derivation order 087 * @param value value of the constant 088 * @throws NumberIsTooLargeException if order is too large. 089 * @see #DerivativeStructure(int, int, int, double) 090 */ 091 public DerivativeStructure(final int parameters, final int order, final double value) { 092 this(parameters, order); 093 this.data[0] = value; 094 } 095 096 /** Build an instance representing a variable. 097 * <p>Instances built using this constructor are considered 098 * to be the free variables with respect to which differentials 099 * are computed. As such, their differential with respect to 100 * themselves is +1.</p> 101 * @param parameters number of free parameters 102 * @param order derivation order 103 * @param index index of the variable (from 0 to {@code parameters - 1}) 104 * @param value value of the variable 105 * @throws NumberIsTooLargeException if {@code index ≥ parameters}. 106 * @see #DerivativeStructure(int, int, double) 107 */ 108 public DerivativeStructure(final int parameters, final int order, 109 final int index, final double value) { 110 this(parameters, order, value); 111 112 if (index >= parameters) { 113 throw new NumberIsTooLargeException(index, parameters, false); 114 } 115 116 if (order > 0) { 117 // the derivative of the variable with respect to itself is 1. 118 data[DSCompiler.getCompiler(index, order).getSize()] = 1.0; 119 } 120 } 121 122 /** Linear combination constructor. 123 * The derivative structure built will be a1 * ds1 + a2 * ds2 124 * @param a1 first scale factor 125 * @param ds1 first base (unscaled) derivative structure 126 * @param a2 second scale factor 127 * @param ds2 second base (unscaled) derivative structure 128 * @throws org.apache.commons.math4.legacy.exception.DimensionMismatchException 129 * if number of free parameters or orders are inconsistent 130 */ 131 public DerivativeStructure(final double a1, final DerivativeStructure ds1, 132 final double a2, final DerivativeStructure ds2) { 133 this(ds1.compiler); 134 compiler.checkCompatibility(ds2.compiler); 135 compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, data, 0); 136 } 137 138 /** Linear combination constructor. 139 * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 140 * @param a1 first scale factor 141 * @param ds1 first base (unscaled) derivative structure 142 * @param a2 second scale factor 143 * @param ds2 second base (unscaled) derivative structure 144 * @param a3 third scale factor 145 * @param ds3 third base (unscaled) derivative structure 146 * @throws org.apache.commons.math4.legacy.exception.DimensionMismatchException 147 * if number of free parameters or orders are inconsistent. 148 */ 149 public DerivativeStructure(final double a1, final DerivativeStructure ds1, 150 final double a2, final DerivativeStructure ds2, 151 final double a3, final DerivativeStructure ds3) { 152 this(ds1.compiler); 153 compiler.checkCompatibility(ds2.compiler); 154 compiler.checkCompatibility(ds3.compiler); 155 compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, a3, ds3.data, 0, data, 0); 156 } 157 158 /** Linear combination constructor. 159 * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 + a4 * ds4 160 * @param a1 first scale factor 161 * @param ds1 first base (unscaled) derivative structure 162 * @param a2 second scale factor 163 * @param ds2 second base (unscaled) derivative structure 164 * @param a3 third scale factor 165 * @param ds3 third base (unscaled) derivative structure 166 * @param a4 fourth scale factor 167 * @param ds4 fourth base (unscaled) derivative structure 168 * @throws DimensionMismatchException if number of free parameters or orders are inconsistent. 169 */ 170 public DerivativeStructure(final double a1, final DerivativeStructure ds1, 171 final double a2, final DerivativeStructure ds2, 172 final double a3, final DerivativeStructure ds3, 173 final double a4, final DerivativeStructure ds4) { 174 this(ds1.compiler); 175 compiler.checkCompatibility(ds2.compiler); 176 compiler.checkCompatibility(ds3.compiler); 177 compiler.checkCompatibility(ds4.compiler); 178 compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, 179 a3, ds3.data, 0, a4, ds4.data, 0, 180 data, 0); 181 } 182 183 /** Build an instance from all its derivatives. 184 * @param parameters number of free parameters 185 * @param order derivation order 186 * @param derivatives derivatives sorted according to 187 * {@link DSCompiler#getPartialDerivativeIndex(int...)} 188 * @throws DimensionMismatchException if derivatives array does not match the 189 * {@link DSCompiler#getSize() size} expected by the compiler. 190 * @throws NumberIsTooLargeException if order is too large. 191 * @see #getAllDerivatives() 192 */ 193 public DerivativeStructure(final int parameters, final int order, final double ... derivatives) { 194 this(parameters, order); 195 if (derivatives.length != data.length) { 196 throw new DimensionMismatchException(derivatives.length, data.length); 197 } 198 System.arraycopy(derivatives, 0, data, 0, data.length); 199 } 200 201 /** Copy constructor. 202 * @param ds instance to copy 203 */ 204 private DerivativeStructure(final DerivativeStructure ds) { 205 this.compiler = ds.compiler; 206 this.data = ds.data.clone(); 207 } 208 209 /** Get the number of free parameters. 210 * @return number of free parameters 211 */ 212 public int getFreeParameters() { 213 return compiler.getFreeParameters(); 214 } 215 216 /** Get the derivation order. 217 * @return derivation order 218 */ 219 public int getOrder() { 220 return compiler.getOrder(); 221 } 222 223 /** Create a constant compatible with instance order and number of parameters. 224 * <p> 225 * This method is a convenience factory method, it simply calls 226 * {@code new DerivativeStructure(getFreeParameters(), getOrder(), c)} 227 * </p> 228 * @param c value of the constant 229 * @return a constant compatible with instance order and number of parameters 230 * @see #DerivativeStructure(int, int, double) 231 * @since 3.3 232 */ 233 public DerivativeStructure createConstant(final double c) { 234 return new DerivativeStructure(getFreeParameters(), getOrder(), c); 235 } 236 237 /** {@inheritDoc} 238 * @since 3.2 239 */ 240 @Override 241 public double getReal() { 242 return data[0]; 243 } 244 245 /** Get the value part of the derivative structure. 246 * @return value part of the derivative structure 247 * @see #getPartialDerivative(int...) 248 */ 249 public double getValue() { 250 return data[0]; 251 } 252 253 /** Get a partial derivative. 254 * @param orders derivation orders with respect to each variable (if all orders are 0, 255 * the value is returned) 256 * @return partial derivative 257 * @see #getValue() 258 * @throws DimensionMismatchException if the numbers of variables does not 259 * match the instance. 260 * @throws NumberIsTooLargeException if the sum of derivation orders is larger than 261 * the instance limits. 262 */ 263 public double getPartialDerivative(final int ... orders) { 264 return data[compiler.getPartialDerivativeIndex(orders)]; 265 } 266 267 /** Get all partial derivatives. 268 * @return a fresh copy of partial derivatives, in an array sorted according to 269 * {@link DSCompiler#getPartialDerivativeIndex(int...)} 270 */ 271 public double[] getAllDerivatives() { 272 return data.clone(); 273 } 274 275 /** {@inheritDoc} 276 * @since 3.2 277 */ 278 @Override 279 public DerivativeStructure add(final double a) { 280 final DerivativeStructure ds = new DerivativeStructure(this); 281 ds.data[0] += a; 282 return ds; 283 } 284 285 /** {@inheritDoc} 286 * @throws DimensionMismatchException if number of free parameters 287 * or orders do not match. 288 */ 289 @Override 290 public DerivativeStructure add(final DerivativeStructure a) { 291 compiler.checkCompatibility(a.compiler); 292 final DerivativeStructure ds = new DerivativeStructure(this); 293 compiler.add(data, 0, a.data, 0, ds.data, 0); 294 return ds; 295 } 296 297 /** {@inheritDoc} 298 * @since 3.2 299 */ 300 @Override 301 public DerivativeStructure subtract(final double a) { 302 return add(-a); 303 } 304 305 /** {@inheritDoc} 306 * @throws DimensionMismatchException if number of free parameters 307 * or orders do not match 308 */ 309 @Override 310 public DerivativeStructure subtract(final DerivativeStructure a) { 311 compiler.checkCompatibility(a.compiler); 312 final DerivativeStructure ds = new DerivativeStructure(this); 313 compiler.subtract(data, 0, a.data, 0, ds.data, 0); 314 return ds; 315 } 316 317 /** {@inheritDoc} */ 318 @Override 319 public DerivativeStructure multiply(final int n) { 320 return multiply((double) n); 321 } 322 323 /** {@inheritDoc} 324 * @since 3.2 325 */ 326 @Override 327 public DerivativeStructure multiply(final double a) { 328 final DerivativeStructure ds = new DerivativeStructure(this); 329 for (int i = 0; i < ds.data.length; ++i) { 330 ds.data[i] *= a; 331 } 332 return ds; 333 } 334 335 /** {@inheritDoc} 336 * @throws DimensionMismatchException if number of free parameters 337 * or orders do not match 338 */ 339 @Override 340 public DerivativeStructure multiply(final DerivativeStructure a) { 341 compiler.checkCompatibility(a.compiler); 342 final DerivativeStructure result = new DerivativeStructure(compiler); 343 compiler.multiply(data, 0, a.data, 0, result.data, 0); 344 return result; 345 } 346 347 /** {@inheritDoc} 348 * @since 3.2 349 */ 350 @Override 351 public DerivativeStructure divide(final double a) { 352 final DerivativeStructure ds = new DerivativeStructure(this); 353 for (int i = 0; i < ds.data.length; ++i) { 354 ds.data[i] /= a; 355 } 356 return ds; 357 } 358 359 /** {@inheritDoc} 360 * @throws DimensionMismatchException if number of free parameters 361 * or orders do not match 362 */ 363 @Override 364 public DerivativeStructure divide(final DerivativeStructure a) { 365 compiler.checkCompatibility(a.compiler); 366 final DerivativeStructure result = new DerivativeStructure(compiler); 367 compiler.divide(data, 0, a.data, 0, result.data, 0); 368 return result; 369 } 370 371 /** {@inheritDoc} */ 372 @Override 373 public DerivativeStructure remainder(final double a) { 374 final DerivativeStructure ds = new DerivativeStructure(this); 375 ds.data[0] = JdkMath.IEEEremainder(ds.data[0], a); 376 return ds; 377 } 378 379 /** {@inheritDoc} 380 * @throws DimensionMismatchException if number of free parameters 381 * or orders do not match 382 * @since 3.2 383 */ 384 @Override 385 public DerivativeStructure remainder(final DerivativeStructure a) { 386 compiler.checkCompatibility(a.compiler); 387 final DerivativeStructure result = new DerivativeStructure(compiler); 388 compiler.remainder(data, 0, a.data, 0, result.data, 0); 389 return result; 390 } 391 392 /** {@inheritDoc} */ 393 @Override 394 public DerivativeStructure negate() { 395 final DerivativeStructure ds = new DerivativeStructure(compiler); 396 for (int i = 0; i < ds.data.length; ++i) { 397 ds.data[i] = -data[i]; 398 } 399 return ds; 400 } 401 402 /** {@inheritDoc} 403 * @since 3.2 404 */ 405 @Override 406 public DerivativeStructure abs() { 407 if (Double.doubleToLongBits(data[0]) < 0) { 408 // we use the bits representation to also handle -0.0 409 return negate(); 410 } else { 411 return this; 412 } 413 } 414 415 /** {@inheritDoc} 416 * @since 3.2 417 */ 418 @Override 419 public DerivativeStructure ceil() { 420 return new DerivativeStructure(compiler.getFreeParameters(), 421 compiler.getOrder(), 422 JdkMath.ceil(data[0])); 423 } 424 425 /** {@inheritDoc} 426 * @since 3.2 427 */ 428 @Override 429 public DerivativeStructure floor() { 430 return new DerivativeStructure(compiler.getFreeParameters(), 431 compiler.getOrder(), 432 JdkMath.floor(data[0])); 433 } 434 435 /** {@inheritDoc} 436 * @since 3.2 437 */ 438 @Override 439 public DerivativeStructure rint() { 440 return new DerivativeStructure(compiler.getFreeParameters(), 441 compiler.getOrder(), 442 JdkMath.rint(data[0])); 443 } 444 445 /** {@inheritDoc} */ 446 @Override 447 public long round() { 448 return JdkMath.round(data[0]); 449 } 450 451 /** {@inheritDoc} 452 * @since 3.2 453 */ 454 @Override 455 public DerivativeStructure signum() { 456 return new DerivativeStructure(compiler.getFreeParameters(), 457 compiler.getOrder(), 458 JdkMath.signum(data[0])); 459 } 460 461 /** {@inheritDoc} 462 * @since 3.2 463 */ 464 @Override 465 public DerivativeStructure copySign(final DerivativeStructure sign){ 466 long m = Double.doubleToLongBits(data[0]); 467 long s = Double.doubleToLongBits(sign.data[0]); 468 if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK 469 return this; 470 } 471 return negate(); // flip sign 472 } 473 474 /** {@inheritDoc} 475 * @since 3.2 476 */ 477 @Override 478 public DerivativeStructure copySign(final double sign) { 479 long m = Double.doubleToLongBits(data[0]); 480 long s = Double.doubleToLongBits(sign); 481 if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK 482 return this; 483 } 484 return negate(); // flip sign 485 } 486 487 /** 488 * Return the exponent of the instance value, removing the bias. 489 * <p> 490 * For double numbers of the form 2<sup>x</sup>, the unbiased 491 * exponent is exactly x. 492 * </p> 493 * @return exponent for instance in IEEE754 representation, without bias 494 */ 495 public int getExponent() { 496 return JdkMath.getExponent(data[0]); 497 } 498 499 /** {@inheritDoc} 500 * @since 3.2 501 */ 502 @Override 503 public DerivativeStructure scalb(final int n) { 504 final DerivativeStructure ds = new DerivativeStructure(compiler); 505 for (int i = 0; i < ds.data.length; ++i) { 506 ds.data[i] = JdkMath.scalb(data[i], n); 507 } 508 return ds; 509 } 510 511 /** {@inheritDoc} 512 * @throws DimensionMismatchException if number of free parameters 513 * or orders do not match 514 * @since 3.2 515 */ 516 @Override 517 public DerivativeStructure hypot(final DerivativeStructure y) { 518 compiler.checkCompatibility(y.compiler); 519 520 if (Double.isInfinite(data[0]) || Double.isInfinite(y.data[0])) { 521 return new DerivativeStructure(compiler.getFreeParameters(), 522 compiler.getFreeParameters(), 523 Double.POSITIVE_INFINITY); 524 } else if (Double.isNaN(data[0]) || Double.isNaN(y.data[0])) { 525 return new DerivativeStructure(compiler.getFreeParameters(), 526 compiler.getFreeParameters(), 527 Double.NaN); 528 } else { 529 530 final int expX = getExponent(); 531 final int expY = y.getExponent(); 532 if (expX > expY + 27) { 533 // y is neglectible with respect to x 534 return abs(); 535 } else if (expY > expX + 27) { 536 // x is neglectible with respect to y 537 return y.abs(); 538 } else { 539 540 // find an intermediate scale to avoid both overflow and underflow 541 final int middleExp = (expX + expY) / 2; 542 543 // scale parameters without losing precision 544 final DerivativeStructure scaledX = scalb(-middleExp); 545 final DerivativeStructure scaledY = y.scalb(-middleExp); 546 547 // compute scaled hypotenuse 548 final DerivativeStructure scaledH = 549 scaledX.multiply(scaledX).add(scaledY.multiply(scaledY)).sqrt(); 550 551 // remove scaling 552 return scaledH.scalb(middleExp); 553 } 554 } 555 } 556 557 /** 558 * Returns the hypotenuse of a triangle with sides {@code x} and {@code y} 559 * - sqrt(<i>x</i><sup>2</sup> +<i>y</i><sup>2</sup>) 560 * avoiding intermediate overflow or underflow. 561 * 562 * <ul> 563 * <li> If either argument is infinite, then the result is positive infinity.</li> 564 * <li> else, if either argument is NaN then the result is NaN.</li> 565 * </ul> 566 * 567 * @param x a value 568 * @param y a value 569 * @return sqrt(<i>x</i><sup>2</sup> +<i>y</i><sup>2</sup>) 570 * @throws DimensionMismatchException if number of free parameters 571 * or orders do not match 572 * @since 3.2 573 */ 574 public static DerivativeStructure hypot(final DerivativeStructure x, final DerivativeStructure y) { 575 return x.hypot(y); 576 } 577 578 /** Compute composition of the instance by a univariate function. 579 * @param f array of value and derivatives of the function at 580 * the current point (i.e. [f({@link #getValue()}), 581 * f'({@link #getValue()}), f''({@link #getValue()})...]). 582 * @return f(this) 583 * @throws DimensionMismatchException if the number of derivatives 584 * in the array is not equal to {@link #getOrder() order} + 1 585 */ 586 public DerivativeStructure compose(final double ... f) { 587 if (f.length != getOrder() + 1) { 588 throw new DimensionMismatchException(f.length, getOrder() + 1); 589 } 590 final DerivativeStructure result = new DerivativeStructure(compiler); 591 compiler.compose(data, 0, f, result.data, 0); 592 return result; 593 } 594 595 /** {@inheritDoc} */ 596 @Override 597 public DerivativeStructure reciprocal() { 598 final DerivativeStructure result = new DerivativeStructure(compiler); 599 compiler.pow(data, 0, -1, result.data, 0); 600 return result; 601 } 602 603 /** {@inheritDoc} 604 * @since 3.2 605 */ 606 @Override 607 public DerivativeStructure sqrt() { 608 return rootN(2); 609 } 610 611 /** {@inheritDoc} 612 * @since 3.2 613 */ 614 @Override 615 public DerivativeStructure cbrt() { 616 return rootN(3); 617 } 618 619 /** {@inheritDoc} 620 * @since 3.2 621 */ 622 @Override 623 public DerivativeStructure rootN(final int n) { 624 final DerivativeStructure result = new DerivativeStructure(compiler); 625 compiler.rootN(data, 0, n, result.data, 0); 626 return result; 627 } 628 629 /** {@inheritDoc} */ 630 @Override 631 public Field<DerivativeStructure> getField() { 632 return new Field<DerivativeStructure>() { 633 634 /** {@inheritDoc} */ 635 @Override 636 public DerivativeStructure getZero() { 637 return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 0.0); 638 } 639 640 /** {@inheritDoc} */ 641 @Override 642 public DerivativeStructure getOne() { 643 return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 1.0); 644 } 645 646 /** {@inheritDoc} */ 647 @Override 648 public Class<? extends FieldElement<DerivativeStructure>> getRuntimeClass() { 649 return DerivativeStructure.class; 650 } 651 }; 652 } 653 654 /** Compute a<sup>x</sup> where a is a double and x a {@link DerivativeStructure}. 655 * @param a number to exponentiate 656 * @param x power to apply 657 * @return a<sup>x</sup> 658 * @since 3.3 659 */ 660 public static DerivativeStructure pow(final double a, final DerivativeStructure x) { 661 final DerivativeStructure result = new DerivativeStructure(x.compiler); 662 x.compiler.pow(a, x.data, 0, result.data, 0); 663 return result; 664 } 665 666 /** {@inheritDoc} 667 * @since 3.2 668 */ 669 @Override 670 public DerivativeStructure pow(final double p) { 671 final DerivativeStructure result = new DerivativeStructure(compiler); 672 compiler.pow(data, 0, p, result.data, 0); 673 return result; 674 } 675 676 /** {@inheritDoc} 677 * @since 3.2 678 */ 679 @Override 680 public DerivativeStructure pow(final int n) { 681 final DerivativeStructure result = new DerivativeStructure(compiler); 682 compiler.pow(data, 0, n, result.data, 0); 683 return result; 684 } 685 686 /** {@inheritDoc} 687 * @throws DimensionMismatchException if number of free parameters 688 * or orders do not match 689 * @since 3.2 690 */ 691 @Override 692 public DerivativeStructure pow(final DerivativeStructure e) { 693 compiler.checkCompatibility(e.compiler); 694 final DerivativeStructure result = new DerivativeStructure(compiler); 695 compiler.pow(data, 0, e.data, 0, result.data, 0); 696 return result; 697 } 698 699 /** {@inheritDoc} 700 * @since 3.2 701 */ 702 @Override 703 public DerivativeStructure exp() { 704 final DerivativeStructure result = new DerivativeStructure(compiler); 705 compiler.exp(data, 0, result.data, 0); 706 return result; 707 } 708 709 /** {@inheritDoc} 710 * @since 3.2 711 */ 712 @Override 713 public DerivativeStructure expm1() { 714 final DerivativeStructure result = new DerivativeStructure(compiler); 715 compiler.expm1(data, 0, result.data, 0); 716 return result; 717 } 718 719 /** {@inheritDoc} 720 * @since 3.2 721 */ 722 @Override 723 public DerivativeStructure log() { 724 final DerivativeStructure result = new DerivativeStructure(compiler); 725 compiler.log(data, 0, result.data, 0); 726 return result; 727 } 728 729 /** {@inheritDoc} 730 * @since 3.2 731 */ 732 @Override 733 public DerivativeStructure log1p() { 734 final DerivativeStructure result = new DerivativeStructure(compiler); 735 compiler.log1p(data, 0, result.data, 0); 736 return result; 737 } 738 739 /** Base 10 logarithm. 740 * @return base 10 logarithm of the instance 741 */ 742 @Override 743 public DerivativeStructure log10() { 744 final DerivativeStructure result = new DerivativeStructure(compiler); 745 compiler.log10(data, 0, result.data, 0); 746 return result; 747 } 748 749 /** {@inheritDoc} 750 * @since 3.2 751 */ 752 @Override 753 public DerivativeStructure cos() { 754 final DerivativeStructure result = new DerivativeStructure(compiler); 755 compiler.cos(data, 0, result.data, 0); 756 return result; 757 } 758 759 /** {@inheritDoc} 760 * @since 3.2 761 */ 762 @Override 763 public DerivativeStructure sin() { 764 final DerivativeStructure result = new DerivativeStructure(compiler); 765 compiler.sin(data, 0, result.data, 0); 766 return result; 767 } 768 769 /** {@inheritDoc} 770 * @since 3.2 771 */ 772 @Override 773 public DerivativeStructure tan() { 774 final DerivativeStructure result = new DerivativeStructure(compiler); 775 compiler.tan(data, 0, result.data, 0); 776 return result; 777 } 778 779 /** {@inheritDoc} 780 * @since 3.2 781 */ 782 @Override 783 public DerivativeStructure acos() { 784 final DerivativeStructure result = new DerivativeStructure(compiler); 785 compiler.acos(data, 0, result.data, 0); 786 return result; 787 } 788 789 /** {@inheritDoc} 790 * @since 3.2 791 */ 792 @Override 793 public DerivativeStructure asin() { 794 final DerivativeStructure result = new DerivativeStructure(compiler); 795 compiler.asin(data, 0, result.data, 0); 796 return result; 797 } 798 799 /** {@inheritDoc} 800 * @since 3.2 801 */ 802 @Override 803 public DerivativeStructure atan() { 804 final DerivativeStructure result = new DerivativeStructure(compiler); 805 compiler.atan(data, 0, result.data, 0); 806 return result; 807 } 808 809 /** {@inheritDoc} 810 * @since 3.2 811 */ 812 @Override 813 public DerivativeStructure atan2(final DerivativeStructure x) { 814 compiler.checkCompatibility(x.compiler); 815 final DerivativeStructure result = new DerivativeStructure(compiler); 816 compiler.atan2(data, 0, x.data, 0, result.data, 0); 817 return result; 818 } 819 820 /** Two arguments arc tangent operation. 821 * @param y first argument of the arc tangent 822 * @param x second argument of the arc tangent 823 * @return atan2(y, x) 824 * @throws DimensionMismatchException if number of free parameters 825 * or orders do not match 826 * @since 3.2 827 */ 828 public static DerivativeStructure atan2(final DerivativeStructure y, final DerivativeStructure x) { 829 return y.atan2(x); 830 } 831 832 /** {@inheritDoc} 833 * @since 3.2 834 */ 835 @Override 836 public DerivativeStructure cosh() { 837 final DerivativeStructure result = new DerivativeStructure(compiler); 838 compiler.cosh(data, 0, result.data, 0); 839 return result; 840 } 841 842 /** {@inheritDoc} 843 * @since 3.2 844 */ 845 @Override 846 public DerivativeStructure sinh() { 847 final DerivativeStructure result = new DerivativeStructure(compiler); 848 compiler.sinh(data, 0, result.data, 0); 849 return result; 850 } 851 852 /** {@inheritDoc} 853 * @since 3.2 854 */ 855 @Override 856 public DerivativeStructure tanh() { 857 final DerivativeStructure result = new DerivativeStructure(compiler); 858 compiler.tanh(data, 0, result.data, 0); 859 return result; 860 } 861 862 /** {@inheritDoc} 863 * @since 3.2 864 */ 865 @Override 866 public DerivativeStructure acosh() { 867 final DerivativeStructure result = new DerivativeStructure(compiler); 868 compiler.acosh(data, 0, result.data, 0); 869 return result; 870 } 871 872 /** {@inheritDoc} 873 * @since 3.2 874 */ 875 @Override 876 public DerivativeStructure asinh() { 877 final DerivativeStructure result = new DerivativeStructure(compiler); 878 compiler.asinh(data, 0, result.data, 0); 879 return result; 880 } 881 882 /** {@inheritDoc} 883 * @since 3.2 884 */ 885 @Override 886 public DerivativeStructure atanh() { 887 final DerivativeStructure result = new DerivativeStructure(compiler); 888 compiler.atanh(data, 0, result.data, 0); 889 return result; 890 } 891 892 /** Convert radians to degrees, with error of less than 0.5 ULP. 893 * @return instance converted into degrees 894 */ 895 public DerivativeStructure toDegrees() { 896 final DerivativeStructure ds = new DerivativeStructure(compiler); 897 for (int i = 0; i < ds.data.length; ++i) { 898 ds.data[i] = JdkMath.toDegrees(data[i]); 899 } 900 return ds; 901 } 902 903 /** Convert degrees to radians, with error of less than 0.5 ULP. 904 * @return instance converted into radians 905 */ 906 public DerivativeStructure toRadians() { 907 final DerivativeStructure ds = new DerivativeStructure(compiler); 908 for (int i = 0; i < ds.data.length; ++i) { 909 ds.data[i] = JdkMath.toRadians(data[i]); 910 } 911 return ds; 912 } 913 914 /** Evaluate Taylor expansion a derivative structure. 915 * @param delta parameters offsets (Δx, Δy, ...) 916 * @return value of the Taylor expansion at x + Δx, y + Δy, ... 917 * @throws org.apache.commons.math4.legacy.exception.MathArithmeticException 918 * if factorials becomes too large 919 */ 920 public double taylor(final double ... delta) { 921 return compiler.taylor(data, 0, delta); 922 } 923 924 /** {@inheritDoc} 925 * @throws DimensionMismatchException if number of free parameters 926 * or orders do not match 927 * @since 3.2 928 */ 929 @Override 930 public DerivativeStructure linearCombination(final DerivativeStructure[] a, final DerivativeStructure[] b) { 931 // compute an accurate value, taking care of cancellations 932 final double[] aDouble = new double[a.length]; 933 for (int i = 0; i < a.length; ++i) { 934 aDouble[i] = a[i].getValue(); 935 } 936 final double[] bDouble = new double[b.length]; 937 for (int i = 0; i < b.length; ++i) { 938 bDouble[i] = b[i].getValue(); 939 } 940 final double accurateValue = Sum.ofProducts(aDouble, bDouble).getAsDouble(); 941 942 // compute a simple value, with all partial derivatives 943 DerivativeStructure simpleValue = a[0].getField().getZero(); 944 for (int i = 0; i < a.length; ++i) { 945 simpleValue = simpleValue.add(a[i].multiply(b[i])); 946 } 947 948 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 949 final double[] all = simpleValue.getAllDerivatives(); 950 all[0] = accurateValue; 951 return new DerivativeStructure(simpleValue.getFreeParameters(), simpleValue.getOrder(), all); 952 } 953 954 /** {@inheritDoc} 955 * @throws DimensionMismatchException if number of free parameters 956 * or orders do not match 957 * @since 3.2 958 */ 959 @Override 960 public DerivativeStructure linearCombination(final double[] a, final DerivativeStructure[] b) { 961 // compute an accurate value, taking care of cancellations 962 final double[] bDouble = new double[b.length]; 963 for (int i = 0; i < b.length; ++i) { 964 bDouble[i] = b[i].getValue(); 965 } 966 final double accurateValue = Sum.ofProducts(a, bDouble).getAsDouble(); 967 968 // compute a simple value, with all partial derivatives 969 DerivativeStructure simpleValue = b[0].getField().getZero(); 970 for (int i = 0; i < a.length; ++i) { 971 simpleValue = simpleValue.add(b[i].multiply(a[i])); 972 } 973 974 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 975 final double[] all = simpleValue.getAllDerivatives(); 976 all[0] = accurateValue; 977 return new DerivativeStructure(simpleValue.getFreeParameters(), simpleValue.getOrder(), all); 978 } 979 980 /** {@inheritDoc} 981 * @throws DimensionMismatchException if number of free parameters 982 * or orders do not match 983 * @since 3.2 984 */ 985 @Override 986 public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1, 987 final DerivativeStructure a2, final DerivativeStructure b2) { 988 // compute an accurate value, taking care of cancellations 989 final double accurateValue = Sum.create() 990 .addProduct(a1.getValue(), b1.getValue()) 991 .addProduct(a2.getValue(), b2.getValue()).getAsDouble(); 992 993 // compute a simple value, with all partial derivatives 994 final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)); 995 996 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 997 final double[] all = simpleValue.getAllDerivatives(); 998 all[0] = accurateValue; 999 return new DerivativeStructure(getFreeParameters(), getOrder(), all); 1000 } 1001 1002 /** {@inheritDoc} 1003 * @throws DimensionMismatchException if number of free parameters 1004 * or orders do not match 1005 * @since 3.2 1006 */ 1007 @Override 1008 public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1, 1009 final double a2, final DerivativeStructure b2) { 1010 // compute an accurate value, taking care of cancellations 1011 final double accurateValue = Sum.create() 1012 .addProduct(a1, b1.getValue()) 1013 .addProduct(a2, b2.getValue()).getAsDouble(); 1014 1015 // compute a simple value, with all partial derivatives 1016 final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)); 1017 1018 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 1019 final double[] all = simpleValue.getAllDerivatives(); 1020 all[0] = accurateValue; 1021 return new DerivativeStructure(getFreeParameters(), getOrder(), all); 1022 } 1023 1024 /** {@inheritDoc} 1025 * @throws DimensionMismatchException if number of free parameters 1026 * or orders do not match 1027 * @since 3.2 1028 */ 1029 @Override 1030 public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1, 1031 final DerivativeStructure a2, final DerivativeStructure b2, 1032 final DerivativeStructure a3, final DerivativeStructure b3) { 1033 // compute an accurate value, taking care of cancellations 1034 final double accurateValue = Sum.create() 1035 .addProduct(a1.getValue(), b1.getValue()) 1036 .addProduct(a2.getValue(), b2.getValue()) 1037 .addProduct(a3.getValue(), b3.getValue()).getAsDouble(); 1038 1039 // compute a simple value, with all partial derivatives 1040 final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3)); 1041 1042 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 1043 final double[] all = simpleValue.getAllDerivatives(); 1044 all[0] = accurateValue; 1045 return new DerivativeStructure(getFreeParameters(), getOrder(), all); 1046 } 1047 1048 /** {@inheritDoc} 1049 * @throws DimensionMismatchException if number of free parameters 1050 * or orders do not match 1051 * @since 3.2 1052 */ 1053 @Override 1054 public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1, 1055 final double a2, final DerivativeStructure b2, 1056 final double a3, final DerivativeStructure b3) { 1057 // compute an accurate value, taking care of cancellations 1058 final double accurateValue = Sum.create() 1059 .addProduct(a1, b1.getValue()) 1060 .addProduct(a2, b2.getValue()) 1061 .addProduct(a3, b3.getValue()).getAsDouble(); 1062 1063 // compute a simple value, with all partial derivatives 1064 final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3)); 1065 1066 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 1067 final double[] all = simpleValue.getAllDerivatives(); 1068 all[0] = accurateValue; 1069 return new DerivativeStructure(getFreeParameters(), getOrder(), all); 1070 } 1071 1072 /** {@inheritDoc} 1073 * @throws DimensionMismatchException if number of free parameters 1074 * or orders do not match 1075 * @since 3.2 1076 */ 1077 @Override 1078 public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1, 1079 final DerivativeStructure a2, final DerivativeStructure b2, 1080 final DerivativeStructure a3, final DerivativeStructure b3, 1081 final DerivativeStructure a4, final DerivativeStructure b4) { 1082 // compute an accurate value, taking care of cancellations 1083 final double accurateValue = Sum.create() 1084 .addProduct(a1.getValue(), b1.getValue()) 1085 .addProduct(a2.getValue(), b2.getValue()) 1086 .addProduct(a3.getValue(), b3.getValue()) 1087 .addProduct(a4.getValue(), b4.getValue()).getAsDouble(); 1088 1089 // compute a simple value, with all partial derivatives 1090 final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3)).add(a4.multiply(b4)); 1091 1092 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 1093 final double[] all = simpleValue.getAllDerivatives(); 1094 all[0] = accurateValue; 1095 return new DerivativeStructure(getFreeParameters(), getOrder(), all); 1096 } 1097 1098 /** {@inheritDoc} 1099 * @throws DimensionMismatchException if number of free parameters 1100 * or orders do not match 1101 * @since 3.2 1102 */ 1103 @Override 1104 public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1, 1105 final double a2, final DerivativeStructure b2, 1106 final double a3, final DerivativeStructure b3, 1107 final double a4, final DerivativeStructure b4) { 1108 // compute an accurate value, taking care of cancellations 1109 final double accurateValue = Sum.create() 1110 .addProduct(a1, b1.getValue()) 1111 .addProduct(a2, b2.getValue()) 1112 .addProduct(a3, b3.getValue()) 1113 .addProduct(a4, b4.getValue()).getAsDouble(); 1114 1115 // compute a simple value, with all partial derivatives 1116 final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3)).add(b4.multiply(a4)); 1117 1118 // create a result with accurate value and all derivatives (not necessarily as accurate as the value) 1119 final double[] all = simpleValue.getAllDerivatives(); 1120 all[0] = accurateValue; 1121 return new DerivativeStructure(getFreeParameters(), getOrder(), all); 1122 } 1123 1124 /** 1125 * Test for the equality of two derivative structures. 1126 * <p> 1127 * Derivative structures are considered equal if they have the same number 1128 * of free parameters, the same derivation order, and the same derivatives. 1129 * </p> 1130 * @param other Object to test for equality to this 1131 * @return true if two derivative structures are equal 1132 * @since 3.2 1133 */ 1134 @Override 1135 public boolean equals(Object other) { 1136 1137 if (this == other) { 1138 return true; 1139 } 1140 1141 if (other instanceof DerivativeStructure) { 1142 final DerivativeStructure rhs = (DerivativeStructure)other; 1143 return getFreeParameters() == rhs.getFreeParameters() && 1144 getOrder() == rhs.getOrder() && 1145 MathArrays.equals(data, rhs.data); 1146 } 1147 1148 return false; 1149 } 1150 1151 /** 1152 * Get a hashCode for the derivative structure. 1153 * @return a hash code value for this object 1154 * @since 3.2 1155 */ 1156 @Override 1157 public int hashCode() { 1158 return 227 + 229 * getFreeParameters() + 233 * getOrder() + 239 * Arrays.hashCode(data); 1159 } 1160}