001package models;
002
003import static java.lang.Math.*;
004
005import org.apache.commons.math3.analysis.UnivariateFunction;
006import org.apache.commons.math3.analysis.integration.SimpsonIntegrator;
007import org.apache.commons.math3.analysis.integration.UnivariateIntegrator;
008import org.apache.commons.math3.complex.Complex;
009
010/**
011 * This class provides helper functionality for computing the minimal solutions
012 * of various models, given by the exponential of a certain integral expression.
013 * 
014 */
015public class SinhTools
016{
017    /**
018     * Computes an integral expression that occurs in the minimal solution of sinh-Gordon-type models.
019     * Explicitly, we compute
020     * 
021     *   \[ \exp 8\int_{0}^{\infty} \frac{dx}{x}\frac{\sinh\Big(\frac{x B}{4}\Big)\sinh\Big(\frac{x}{2}(1-\frac{B}{2}) \Big)
022     *   \sinh\frac{x}{2}}{\sinh^2 x} \sin^2\Big( \frac{x\theta}{2\pi}\Big) \]
023     *  
024     *  where \(\theta \in \mathbb{R} \) and \( B \in (0,2) \) are given as a parameter.  
025     *  
026     *  The integral is numerically approximated using the composite Simpson's rule.
027     *  
028     * @param theta the parameter \( \theta \) in the integrand
029     * @param b the parameter \( B \) in the integrand
030     * @return the value of the integral expression as described above
031     */
032    public static double sinhIntegral(final double theta, final double b) {
033
034        UnivariateIntegrator integrator = new SimpsonIntegrator(1e-12, 1e-15, 
035                        SimpsonIntegrator.DEFAULT_MIN_ITERATIONS_COUNT, SimpsonIntegrator.SIMPSON_MAX_ITERATIONS_COUNT);
036        
037        UnivariateFunction integrand = new UnivariateFunction() {
038                public double value(double x) {
039                    return 1.0/x*sinh(x*b/4)*sinh((2-b)/4*x)*sinh(x/2)/pow(sinh(x), 2) * pow(sin(x*theta/2/PI), 2);
040                }
041            };
042        double intval = integrator.integrate(1000000, integrand, 1e-10, 30);
043        return exp(8*intval);
044    }
045
046    /**
047     * Computes the limit of the integral expression {@link #sinhIntegral} as \(\theta \to \infty\).
048     * By the Riemann-Lebesgue lemma, this limit is given by
049     * 
050     *   \[ \exp 4\int_{0}^{\infty} \frac{dx}{x}\frac{\sinh\Big(\frac{x B}{4}\Big)\sinh\Big(\frac{x}{2}(1-\frac{B}{2}) \Big)
051     *   \sinh\frac{x}{2}}{\sinh^2 x}  \]
052     *  
053     *  The integral is numerically approximated using the composite Simpson's rule.
054     *  
055     * @param b the parameter \( B \) in the integrand
056     * @return the value of the integral expression as described above
057     */
058    public static double sinhIntegralLimit(final double b) {
059
060        UnivariateIntegrator integrator = new SimpsonIntegrator(1e-12, 1e-15,
061                        SimpsonIntegrator.DEFAULT_MIN_ITERATIONS_COUNT, SimpsonIntegrator.SIMPSON_MAX_ITERATIONS_COUNT);
062        
063        UnivariateFunction integrand = new UnivariateFunction() {
064                public double value(double x) {
065                    return 1.0/x*sinh(x*b/4)*sinh((2-b)/4*x)*sinh(x/2)/pow(sinh(x), 2);
066                }
067            };
068        double intval = integrator.integrate(1000000, integrand, 1e-10, 30);
069        return exp(4*intval);
070    }
071
072    /**
073     * Computes the product of two integral expression as in {@link #sinhIntegral} but for a complex value \(B\) 
074     * and for its complex conjugate.
075     * 
076     * The expression is explicitly given by 
077     *   \[ \exp 8\int_{0}^{\infty} \frac{dx}{x}\frac{ 2 \operatorname{Re} \Big(\sinh\Big(\frac{x B}{4}\Big)\sinh\Big(\frac{x}{2}(1-\frac{B}{2})  \Big) \Big)
078     *   \sinh\frac{x}{2}}{\sinh^2 x} \; \sin^2\Big( \frac{x\theta}{2\pi}\Big) \]
079     *  
080     *  where \(\theta \in \mathbb{R} \) and \( B \in (0,2) + i \mathbb{R} \) are given as parameters.  
081     * 
082     * The integral is numerically approximated using the composite Simpson's rule.
083     * @param theta the parameter \( \theta \) in the integrand
084     * @param b the complex parameter \( B \) in the integrand
085     * @return the value of the integral expression as described above
086     */ 
087    public static double sinhIntegralCCPair(final double theta, final Complex b) {
088
089        UnivariateIntegrator integrator = new SimpsonIntegrator(1e-12, 1e-15, 
090                        SimpsonIntegrator.DEFAULT_MIN_ITERATIONS_COUNT, SimpsonIntegrator.SIMPSON_MAX_ITERATIONS_COUNT);
091
092        UnivariateFunction integrand = new UnivariateFunction() {
093                public double value(double x) {
094                        Complex sinhFact1 = b.multiply(x/4).sinh();
095                        Complex sinhFact2 = new Complex(2, 0).subtract(b).multiply(x/4).sinh();
096                        double sinhFact = sinhFact1.multiply(sinhFact2).getReal()*2;
097                    return 1.0/x*sinhFact*sinh(x/2)/pow(sinh(x), 2) * pow(sin(x*theta/2/PI), 2);
098                }
099            };
100        double intval = integrator.integrate(1000000, integrand, 1e-10, 30);
101        return exp(8*intval);
102    }
103
104    /**
105     * Computes the limit of the integral expression {@link #sinhIntegralCCPair} as \(\theta \to \infty\).
106     * By the Riemann-Lebesgue lemma, this limit is given by
107     *  \[ \exp 4\int_{0}^{\infty} \frac{dx}{x}\frac{ 2 \operatorname{Re} \Big(\sinh\Big(\frac{x B}{4}\Big)\sinh\Big(\frac{x}{2}(1-\frac{B}{2})  \Big) \Big)
108     *     \sinh\frac{x}{2}}{\sinh^2 x}  \]
109     *  
110     *  
111     *  The integral is numerically approximated using the composite Simpson's rule.
112     *  
113     * @param b the parameter \( B \) in the integrand
114     * @return the value of the integral expression as described above
115     */
116    public static double sinhIntegralCCPairLimit(final Complex b) {
117
118        UnivariateIntegrator integrator = new SimpsonIntegrator(1e-12, 1e-15, 
119                        SimpsonIntegrator.DEFAULT_MIN_ITERATIONS_COUNT, SimpsonIntegrator.SIMPSON_MAX_ITERATIONS_COUNT);
120
121        UnivariateFunction integrand = new UnivariateFunction() {
122                public double value(double x) {
123                        Complex sinhFact1 = b.multiply(x/4).sinh();
124                        Complex sinhFact2 = new Complex(2, 0).subtract(b).multiply(x/4).sinh();
125                        double sinhFact = sinhFact1.multiply(sinhFact2).getReal()*2;
126                    return 1.0/x*sinhFact*sinh(x/2)/pow(sinh(x), 2) ;
127                }
128            };
129        double intval = integrator.integrate(1000000, integrand, 1e-10, 30);
130        return exp(4*intval);
131    }
132
133}