io.yields.math.concepts.operator.CurveLength.java Source code

Introduction

Here is the source code for io.yields.math.concepts.operator.CurveLength.java

Source

/*
 * Copyright 2014 by Yields.
 *
 * Licensed 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 io.yields.math.concepts.operator;

import org.apache.commons.math3.analysis.UnivariateFunction;
import org.apache.commons.math3.analysis.integration.RombergIntegrator;
import org.apache.commons.math3.analysis.integration.UnivariateIntegrator;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;

import java.util.function.Function;
import java.util.stream.IntStream;

/**
 * Curve Length is an Operator that measures the length of a one-dimensional RealFunction.
 */
public class CurveLength implements Function<Function<Double, Double>, Double> {

    private static final Logger logger = LoggerFactory.getLogger(CurveLength.class);

    private static final double DELTA = 1.e-4;

    private double min, max;
    private UnivariateIntegrator integrator;
    protected static final int MAX_EVAL = 1000;

    protected CurveLength() {
        this.integrator = new RombergIntegrator();
    }

    /**
     * Helper function for use in the testing framework. This method sets the start and end point of the length measurement.
     */
    public CurveLength between(double min, double max) {
        this.min = min;
        this.max = max;
        return this;
    }

    @Override
    public Double apply(Function<Double, Double> function) {
        try {
            return integrator.integrate(MAX_EVAL, getLengthFunction(function), min, max);
        } catch (Exception e) {
            logger.warn("Failed to integrate using Romberg algorithm - switching to Simpson algorithm.", e);
            return simpsonIntegration(function);
        }
    }

    protected Double simpsonIntegration(Function<Double, Double> originalFunction) {
        UnivariateFunction function = getLengthFunction(originalFunction);
        double delta = (max - min) / MAX_EVAL;
        double integral = IntStream.rangeClosed(0, MAX_EVAL).mapToDouble(index -> function.value(index * delta))
                .sum() * delta;
        return integral - .5 * delta * (function.value(min) + function.value(max));
    }

    private UnivariateFunction getLengthFunction(Function<Double, Double> function) {
        double epsilon = Math.min((max - min) / MAX_EVAL, DELTA);
        UnivariateFunction dFunction = x -> (function.apply(x + epsilon) - function.apply(x - epsilon))
                / (2 * epsilon);
        return x -> Math.sqrt(1 + Math.pow(dFunction.value(x), 2));
    }

}