/*
 * Licensed to the Apache Software Foundation (ASF) under one or more
 * contributor license agreements.  See the NOTICE file distributed with
 * this work for additional information regarding copyright ownership.
 * The ASF 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
 *
 *      http://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.apache.commons.math3.ode;

import java.util.ArrayList;
import java.util.List;

import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.MaxCountExceededException;


This class represents a combined set of first order differential equations,
with at least a primary set of equations expandable by some sets of secondary
equations.

One typical use case is the computation of the Jacobian matrix for some ODE.
In this case, the primary set of equations corresponds to the raw ODE, and we
add to this set another bunch of secondary equations which represent the Jacobian
matrix of the primary set.


We want the integrator to use only the primary set to estimate the
errors and hence the step sizes. It should not use the secondary equations in this computation. The integrator will be able to know where the primary set ends and so where the secondary sets begin. 
See Also:
FirstOrderDifferentialEquations
JacobianMatrices
Since:
3.0/**
 * This class represents a combined set of first order differential equations,
 * with at least a primary set of equations expandable by some sets of secondary
 * equations.
 * <p>
 * One typical use case is the computation of the Jacobian matrix for some ODE.
 * In this case, the primary set of equations corresponds to the raw ODE, and we
 * add to this set another bunch of secondary equations which represent the Jacobian
 * matrix of the primary set.
 * </p>
 * <p>
 * We want the integrator to use <em>only</em> the primary set to estimate the
 * errors and hence the step sizes. It should <em>not</em> use the secondary
 * equations in this computation. The {@link AbstractIntegrator integrator} will
 * be able to know where the primary set ends and so where the secondary sets begin.
 * </p>
 *
 * @see FirstOrderDifferentialEquations
 * @see JacobianMatrices
 *
 * @since 3.0
 */

public class ExpandableStatefulODE {

    
Primary differential equation. /** Primary differential equation. */
    private final FirstOrderDifferentialEquations primary;

    
Mapper for primary equation. /** Mapper for primary equation. */
    private final EquationsMapper primaryMapper;

    
Time. /** Time. */
    private double time;

    
State. /** State. */
    private final double[] primaryState;

    
State derivative. /** State derivative. */
    private final double[] primaryStateDot;

    
Components of the expandable ODE. /** Components of the expandable ODE. */
    private List<SecondaryComponent> components;

    
Build an expandable set from its primary ODE set.
Params:
primary – the primary set of differential equations to be integrated./** Build an expandable set from its primary ODE set.
     * @param primary the primary set of differential equations to be integrated.
     */
    public ExpandableStatefulODE(final FirstOrderDifferentialEquations primary) {
        final int n          = primary.getDimension();
        this.primary         = primary;
        this.primaryMapper   = new EquationsMapper(0, n);
        this.time            = Double.NaN;
        this.primaryState    = new double[n];
        this.primaryStateDot = new double[n];
        this.components      = new ArrayList<ExpandableStatefulODE.SecondaryComponent>();
    }

    
Get the primary set of differential equations.
Returns:
primary set of differential equations/** Get the primary set of differential equations.
     * @return primary set of differential equations
     */
    public FirstOrderDifferentialEquations getPrimary() {
        return primary;
    }

    
Return the dimension of the complete set of equations.

The complete set of equations correspond to the primary set plus all secondary sets.

Returns:
dimension of the complete set of equations/** Return the dimension of the complete set of equations.
     * <p>
     * The complete set of equations correspond to the primary set plus all secondary sets.
     * </p>
     * @return dimension of the complete set of equations
     */
    public int getTotalDimension() {
        if (components.isEmpty()) {
            // there are no secondary equations, the complete set is limited to the primary set
            return primaryMapper.getDimension();
        } else {
            // there are secondary equations, the complete set ends after the last set
            final EquationsMapper lastMapper = components.get(components.size() - 1).mapper;
            return lastMapper.getFirstIndex() + lastMapper.getDimension();
        }
    }

    
Get the current time derivative of the complete state vector.
Params:
t – current value of the independent time variable
y – array containing the current value of the complete state vector
yDot – placeholder array where to put the time derivative of the complete state vector
Throws:
MaxCountExceededException – if the number of functions evaluations is exceeded
DimensionMismatchException – if arrays dimensions do not match equations settings/** Get the current time derivative of the complete state vector.
     * @param t current value of the independent <I>time</I> variable
     * @param y array containing the current value of the complete state vector
     * @param yDot placeholder array where to put the time derivative of the complete state vector
     * @exception MaxCountExceededException if the number of functions evaluations is exceeded
     * @exception DimensionMismatchException if arrays dimensions do not match equations settings
     */
    public void computeDerivatives(final double t, final double[] y, final double[] yDot)
        throws MaxCountExceededException, DimensionMismatchException {

        // compute derivatives of the primary equations
        primaryMapper.extractEquationData(y, primaryState);
        primary.computeDerivatives(t, primaryState, primaryStateDot);

        // Add contribution for secondary equations
        for (final SecondaryComponent component : components) {
            component.mapper.extractEquationData(y, component.state);
            component.equation.computeDerivatives(t, primaryState, primaryStateDot,
                                                  component.state, component.stateDot);
            component.mapper.insertEquationData(component.stateDot, yDot);
        }

        primaryMapper.insertEquationData(primaryStateDot, yDot);

    }

    
Add a set of secondary equations to be integrated along with the primary set.
Params:
secondary – secondary equations set
Returns:
index of the secondary equation in the expanded state/** Add a set of secondary equations to be integrated along with the primary set.
     * @param secondary secondary equations set
     * @return index of the secondary equation in the expanded state
     */
    public int addSecondaryEquations(final SecondaryEquations secondary) {

        final int firstIndex;
        if (components.isEmpty()) {
            // lazy creation of the components list
            components = new ArrayList<ExpandableStatefulODE.SecondaryComponent>();
            firstIndex = primary.getDimension();
        } else {
            final SecondaryComponent last = components.get(components.size() - 1);
            firstIndex = last.mapper.getFirstIndex() + last.mapper.getDimension();
        }

        components.add(new SecondaryComponent(secondary, firstIndex));

        return components.size() - 1;

    }

    
Get an equations mapper for the primary equations set.
See Also:
getSecondaryMappers()
Returns:
mapper for the primary set/** Get an equations mapper for the primary equations set.
     * @return mapper for the primary set
     * @see #getSecondaryMappers()
     */
    public EquationsMapper getPrimaryMapper() {
        return primaryMapper;
    }

    
Get the equations mappers for the secondary equations sets.
See Also:
getPrimaryMapper()
Returns:
equations mappers for the secondary equations sets/** Get the equations mappers for the secondary equations sets.
     * @return equations mappers for the secondary equations sets
     * @see #getPrimaryMapper()
     */
    public EquationsMapper[] getSecondaryMappers() {
        final EquationsMapper[] mappers = new EquationsMapper[components.size()];
        for (int i = 0; i < mappers.length; ++i) {
            mappers[i] = components.get(i).mapper;
        }
        return mappers;
    }

    
Set current time.
Params:
time – current time/** Set current time.
     * @param time current time
     */
    public void setTime(final double time) {
        this.time = time;
    }

    
Get current time.
Returns:
current time/** Get current time.
     * @return current time
     */
    public double getTime() {
        return time;
    }

    
Set primary part of the current state.
Params:
primaryState – primary part of the current state
Throws:
DimensionMismatchException – if the dimension of the array does not
match the primary set/** Set primary part of the current state.
     * @param primaryState primary part of the current state
     * @throws DimensionMismatchException if the dimension of the array does not
     * match the primary set
     */
    public void setPrimaryState(final double[] primaryState) throws DimensionMismatchException {

        // safety checks
        if (primaryState.length != this.primaryState.length) {
            throw new DimensionMismatchException(primaryState.length, this.primaryState.length);
        }

        // set the data
        System.arraycopy(primaryState, 0, this.primaryState, 0, primaryState.length);

    }

    
Get primary part of the current state.
Returns:
primary part of the current state/** Get primary part of the current state.
     * @return primary part of the current state
     */
    public double[] getPrimaryState() {
        return primaryState.clone();
    }

    
Get primary part of the current state derivative.
Returns:
primary part of the current state derivative/** Get primary part of the current state derivative.
     * @return primary part of the current state derivative
     */
    public double[] getPrimaryStateDot() {
        return primaryStateDot.clone();
    }

    
Set secondary part of the current state.
Params:
index – index of the part to set as returned by addSecondaryEquations(SecondaryEquations)
secondaryState – secondary part of the current state
Throws:
DimensionMismatchException – if the dimension of the partial state does not
match the selected equations set dimension/** Set secondary part of the current state.
     * @param index index of the part to set as returned by {@link
     * #addSecondaryEquations(SecondaryEquations)}
     * @param secondaryState secondary part of the current state
     * @throws DimensionMismatchException if the dimension of the partial state does not
     * match the selected equations set dimension
     */
    public void setSecondaryState(final int index, final double[] secondaryState)
        throws DimensionMismatchException {

        // get either the secondary state
        double[] localArray = components.get(index).state;

        // safety checks
        if (secondaryState.length != localArray.length) {
            throw new DimensionMismatchException(secondaryState.length, localArray.length);
        }

        // set the data
        System.arraycopy(secondaryState, 0, localArray, 0, secondaryState.length);

    }

    
Get secondary part of the current state.
Params:
index – index of the part to set as returned by addSecondaryEquations(SecondaryEquations)
Returns:
secondary part of the current state/** Get secondary part of the current state.
     * @param index index of the part to set as returned by {@link
     * #addSecondaryEquations(SecondaryEquations)}
     * @return secondary part of the current state
     */
    public double[] getSecondaryState(final int index) {
        return components.get(index).state.clone();
    }

    
Get secondary part of the current state derivative.
Params:
index – index of the part to set as returned by addSecondaryEquations(SecondaryEquations)
Returns:
secondary part of the current state derivative/** Get secondary part of the current state derivative.
     * @param index index of the part to set as returned by {@link
     * #addSecondaryEquations(SecondaryEquations)}
     * @return secondary part of the current state derivative
     */
    public double[] getSecondaryStateDot(final int index) {
        return components.get(index).stateDot.clone();
    }

    
Set the complete current state.
Params:
completeState – complete current state to copy data from
Throws:
DimensionMismatchException – if the dimension of the complete state does not
match the complete equations sets dimension/** Set the complete current state.
     * @param completeState complete current state to copy data from
     * @throws DimensionMismatchException if the dimension of the complete state does not
     * match the complete equations sets dimension
     */
    public void setCompleteState(final double[] completeState)
        throws DimensionMismatchException {

        // safety checks
        if (completeState.length != getTotalDimension()) {
            throw new DimensionMismatchException(completeState.length, getTotalDimension());
        }

        // set the data
        primaryMapper.extractEquationData(completeState, primaryState);
        for (final SecondaryComponent component : components) {
            component.mapper.extractEquationData(completeState, component.state);
        }

    }

    
Get the complete current state.
Throws:
DimensionMismatchException – if the dimension of the complete state does not
match the complete equations sets dimension
Returns:
complete current state/** Get the complete current state.
     * @return complete current state
     * @throws DimensionMismatchException if the dimension of the complete state does not
     * match the complete equations sets dimension
     */
    public double[] getCompleteState() throws DimensionMismatchException {

        // allocate complete array
        double[] completeState = new double[getTotalDimension()];

        // set the data
        primaryMapper.insertEquationData(primaryState, completeState);
        for (final SecondaryComponent component : components) {
            component.mapper.insertEquationData(component.state, completeState);
        }

        return completeState;

    }

    
Components of the compound stateful ODE. /** Components of the compound stateful ODE. */
    private static class SecondaryComponent {

        
Secondary differential equation. /** Secondary differential equation. */
        private final SecondaryEquations equation;

        
Mapper between local and complete arrays. /** Mapper between local and complete arrays. */
        private final EquationsMapper mapper;

        
State. /** State. */
        private final double[] state;

        
State derivative. /** State derivative. */
        private final double[] stateDot;

        
Simple constructor.
Params:
equation – secondary differential equation
firstIndex – index to use for the first element in the complete arrays/** Simple constructor.
         * @param equation secondary differential equation
         * @param firstIndex index to use for the first element in the complete arrays
         */
        SecondaryComponent(final SecondaryEquations equation, final int firstIndex) {
            final int n   = equation.getDimension();
            this.equation = equation;
            mapper        = new EquationsMapper(firstIndex, n);
            state         = new double[n];
            stateDot      = new double[n];
        }

    }

}