/*
    * Licensed to the Hipparchus project under one or more
    * contributor license agreements.  See the NOTICE file distributed with
    * this work for additional information regarding copyright ownership.
    * The Hipparchus project licenses this file to You under the Apache License, Version 2.0
    * (the "License"); you may not use this file except in compliance with
    * the License.  You may obtain a copy of the License at
    *
    *      https://www.apache.org/licenses/LICENSE-2.0
   *
   * Unless required by applicable law or agreed to in writing, software
   * distributed under the License is distributed on an "AS IS" BASIS,
   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
   * See the License for the specific language governing permissions and
   * limitations under the License.
   */
  
  package org.hipparchus.ode.nonstiff;
  
  
  import org.hipparchus.analysis.UnivariateFunction;
  import org.hipparchus.analysis.polynomials.PolynomialFunction;
  import org.hipparchus.linear.Array2DRowRealMatrix;
  import org.junit.jupiter.api.Test;
  
  import static org.junit.jupiter.api.Assertions.assertEquals;
  
  class AdamsNordsieckTransformerTest {
  
      @Test
      void testPolynomialExtraDerivative() {
          checkNordsieckStart(new PolynomialFunction(new double[] { 6, 5, 4, 3, 2, 1 }),
                              5, 0.0, 0.125, 3.2e-16);
      }
  
      @Test
      void testPolynomialRegular() {
          checkNordsieckStart(new PolynomialFunction(new double[] { 6, 5, 4, 3, 2, 1 }),
                              4, 0.0, 0.125, 3.1e-16);
      }
  
      @Test
      void testPolynomialMissingLastDerivative() {
          // this test intentionally uses not enough start points,
          // the Nordsieck vector is therefore not expected to match the exact scaled derivatives
          checkNordsieckStart(new PolynomialFunction(new double[] { 6, 5, 4, 3, 2, 1 }),
                              3, 0.0, 0.125, 1.6e-4);
      }
  
      @Test
      void testTransformExact() {
          // a 5 steps transformer handles a degree 5 polynomial exactly
          // the Nordsieck vector holds the full information about the function
          // transforming the vector from t0 to t0+h or recomputing it from scratch
          // at t0+h yields the same result
          checkTransform(new PolynomialFunction(new double[] { 6, 5, 4, 3, 2, 1 }), 5, 2.567e-15);
      }
  
      @Test
      void testTransformInexact() {
          // a 4 steps transformer cannot handle a degree 5 polynomial exactly
          // the Nordsieck vector lacks some high degree information about the function
          // transforming the vector from t0 to t0+h or recomputing it from scratch
          // at t0+h yields different results
          checkTransform(new PolynomialFunction(new double[] { 6, 5, 4, 3, 2, 1 }), 4, 5.658e-4);
      }
  
      private void checkNordsieckStart(final PolynomialFunction polynomial, final int nbSteps, final double t0,
                                       final double h, final double epsilon) {
  
          final AdamsNordsieckTransformer transformer = AdamsNordsieckTransformer.getInstance(nbSteps);
          PolynomialFunction derivative = polynomial.polynomialDerivative();
          final Array2DRowRealMatrix nordsieck = start(transformer, nbSteps, t0, h, polynomial, derivative);
  
          assertEquals(nbSteps - 1, nordsieck.getRowDimension());
          double coeff = h;
          for (int i = 0; i < nordsieck.getRowDimension(); ++i) {
              coeff *= h / (i + 2);
              derivative = derivative.polynomialDerivative();
              assertEquals(derivative.value(t0) * coeff, nordsieck.getEntry(i, 0), epsilon);
          }
  
      }
  
      private void checkTransform(final PolynomialFunction polynomial, final int nbSteps, final double expectedError) {
  
          final AdamsNordsieckTransformer transformer = AdamsNordsieckTransformer.getInstance(nbSteps);
          final PolynomialFunction derivative = polynomial.polynomialDerivative();
  
          final double t0 = 0.0;
          final double h  = 0.125;
          final Array2DRowRealMatrix n0 = start(transformer, nbSteps, t0, h, polynomial, derivative);
          final Array2DRowRealMatrix n1 = transformer.updateHighOrderDerivativesPhase1(n0);
          transformer.updateHighOrderDerivativesPhase2(new double[] { h * derivative.value(t0)     },
                                                       new double[] { h * derivative.value(t0 + h) },
                                                       n1);
          final Array2DRowRealMatrix n2 = start(transformer, nbSteps, t0 + h, h, polynomial, derivative);
  
          assertEquals(expectedError, n2.subtract(n1).getNorm1(), expectedError * 0.001);
 
     }
 
     private Array2DRowRealMatrix start(final AdamsNordsieckTransformer transformer, final int nbSteps,
                                        final double t0, final double h,
                                        final UnivariateFunction f0, final UnivariateFunction f1) {
 
         final int        nbStartPoints = (nbSteps + 3) / 2;
         final double[]   t             = new double[nbStartPoints];
         final double[][] y             = new double[nbStartPoints][1];
         final double[][] yDot          = new double[nbStartPoints][1];
         for (int i = 0; i < nbStartPoints; ++i) {
             t[i]       = t0 + i * h;
             y[i][0]    = f0.value(t[i]);
             yDot[i][0] = f1.value(t[i]);
         }
 
         return transformer.initializeHighOrderDerivatives(h, t, y, yDot);
 
     }
 
 }