Examples with Complex org.hipparchus.complex.Complex used on opensource projects

Example 26 with Complex

use of org.hipparchus.complex.Complex in project symja_android_library by axkr.

the class GammaJS method gamma.

/**
 * Incomplete gamma function.
 *
 * @param x
 * @param y
 * @return
 */
public static Complex gamma(Complex x, Complex y) {
    // patch lower end or evaluate exponential integral independently
    if (F.isZero(x)) {
        if (F.isZero(y)) {
            throw new ArgumentTypeException("Gamma function pole");
        }
        // if (Complex.equals(x, Complex.ZERO, Config.SPECIAL_FUNCTIONS_TOLERANCE)) {
        // // taylorSeries => (-EulerGamma - Log(y)) + x - 1/4 * x^2 + 1/18 * x^3 - 1/96 * x^4 + 1/600
        // * x^5
        // Complex result = y.log().add(ConstantDefinitions.EULER_GAMMA).negate();
        // double[] coeff = new double[] { 1.0, -0.25, 1.0 / 18.0, -1.0 / 96.0, 1.0 / 600.0 };
        // Complex yPow = y;
        // for (int i = 0; i < coeff.length; i++) {
        // result = result.add(yPow.multiply(coeff[i]));
        // yPow = yPow.multiply(y);
        // }
        // return result;
        Complex yLogNegate = y.log().negate();
        Complex result = // 
        expIntegralEi(y.negate()).negate().add(// 
        y.negate().log().subtract(y.reciprocal().negate().log()).multiply(0.5)).add(yLogNegate);
        if (F.isZero(y.getImaginary()) && y.getReal() > 0.0) {
            return new Complex(result.getReal());
        }
        return result;
    }
    double delta = 1e-5;
    if (x.norm() < delta) {
    // TODO implement case for abs value near 0
    // return taylorSeries( t => gamma(t,y), mul( x, delta/x.abs( ) ), 2.0)(x);
    }
    // dlmf.nist.gov/8.4.15
    double xRe = x.getReal();
    if (xRe < 0.0 && F.isNumIntValue(xRe) && F.isZero(x.getImaginary())) {
        // x is a negative integer
        final double n = -xRe;
        int iterationLimit = EvalEngine.get().getIterationLimit();
        final Complex t = // 
        y.negate().exp().multiply(ZetaJS.complexSummation(// 
        k -> new Complex(Math.pow(-1.0, k) * factorialInt(k)).divide(y.pow(k + 1)), 0.0, n - 1.0, iterationLimit));
        // dlmf.nist.gov/8.4.4
        final double plusMinusOne = Math.pow(-1.0, n);
        return gamma(Complex.ZERO, y).subtract(t).multiply(plusMinusOne / factorialInt(n));
    }
    return Arithmetic.lanczosApproxGamma(x).subtract(gamma(x, 0.0, y));
}

Example 27 with Complex

use of org.hipparchus.complex.Complex in project symja_android_library by axkr.

the class ZetaJS method hurwitzZeta.

/**
 * @param x
 * @param a
 * @return
 */
public static Complex hurwitzZeta(final Complex x, final Complex a) {
    // TODO INVALID at the moment
    if (x.getReal() == 1.0 && x.getImaginary() == 0.0) {
        throw new ArgumentTypeException("Hurwitz zeta pole");
    }
    // direct summation more accurate than dlmf.nist.gov/25.11.4 for positive a
    int iterationLimit = EvalEngine.get().getIterationLimit();
    if (a.getReal() < 0.0) {
        double m = -Math.floor(a.getReal());
        return hurwitzZeta(x, a.add(m)).add(summation(i -> a.add(i).pow(x.negate()), 0, m - 1.0, iterationLimit));
    }
    // Johansson arxiv.org/abs/1309.2877
    // recommendation of Vepstas, Efficient Algorithm, p.12
    int n = 15;
    // Euler-Maclaurin has differences of large values in left-hand plane
    boolean useArbitrary = x.getReal() < 0;
    // if ( useArbitrary ) {
    // 
    // }
    double switchForms = -5.0;
    if (x.getReal() < switchForms) {
        throw new ArgumentTypeException("Currently unsuppported complex Hurwitz zeta");
    }
    Complex S = summation(i -> a.add(i).pow(x.negate()), 0, n - 1, iterationLimit);
    Complex I = a.add(n).pow(Complex.ONE.subtract(x)).divide(x.subtract(1.0));
    Complex p = x.multiply(0.5).multiply(a.add(n).reciprocal());
    Complex t = p.multiply(bernoulliInt(2));
    int i = 1;
    // converges rather quickly
    while (Math.abs(p.getReal()) > Config.SPECIAL_FUNCTIONS_TOLERANCE || Math.abs(p.getImaginary()) > Config.SPECIAL_FUNCTIONS_TOLERANCE) {
        if (i++ > iterationLimit && iterationLimit > 0) {
            IterationLimitExceeded.throwIt(i, org.matheclipse.core.expression.S.HurwitzZeta);
        }
        if (i > MAX_VALUE_HALF) {
            throw new ArgumentTypeException("Hurwitz zeta: i > MAX_VALUE_HALF");
        }
        int iPlusi = i + i;
        p = p.multiply(x.add(iPlusi - 2.0).multiply(x.add(iPlusi - 3.0)).multiply(a.add(n).pow(2.0).multiply(iPlusi * (iPlusi - 1)).reciprocal()));
        t = t.add(p.multiply(bernoulliInt(iPlusi)));
    }
    Complex T = t.add(0.5).divide(a.add(n).pow(x));
    return S.add(I).add(T);
}

Example 28 with Complex

use of org.hipparchus.complex.Complex in project symja_android_library by axkr.

the class ZetaJS method complexAverage.

public static Complex complexAverage(java.util.function.Function<Complex, Complex> f, Complex x) {
    double offset = 1e-5;
    Complex arg1 = f.apply(x.add(offset));
    Complex arg2 = f.apply(x.subtract(offset));
    return arg1.add(arg2).divide(2.0);
}

Example 29 with Complex

use of org.hipparchus.complex.Complex in project symja_android_library by axkr.

the class Convert method complexMatrix2List.

public static IASTAppendable complexMatrix2List(final FieldMatrix<Complex> matrix, boolean matrixFormat) {
    if (matrix == null) {
        return F.NIL;
    }
    final int rowSize = matrix.getRowDimension();
    final int colSize = matrix.getColumnDimension();
    final IASTAppendable out = F.ListAlloc(rowSize);
    IASTAppendable currOutRow;
    for (int i = 0; i < rowSize; i++) {
        currOutRow = F.ListAlloc(colSize);
        out.append(currOutRow);
        for (int j = 0; j < colSize; j++) {
            final Complex expr = matrix.getEntry(i, j);
            currOutRow.append(F.complexNum(expr));
        }
    }
    if (matrixFormat) {
        // because the rows can contain sub lists the IAST.IS_MATRIX flag cannot be set directly.
        // isMatrix() must be
        // used!
        out.isMatrix(true);
    }
    return out;
}

Example 30 with Complex

use of org.hipparchus.complex.Complex in project symja_android_library by axkr.

the class BesselJS method besselI.

public static Complex besselI(Complex n, Complex x) {
    if (F.isNumIntValue(n.getReal()) && n.getReal() < 0 && n.getImaginary() == 0.0) {
        return besselI(n.negate(), x);
    }
    Complex product = x.divide(2).pow(n).divide(Arithmetic.lanczosApproxGamma(n.add(1)));
    Complex sqrX = x.multiply(x);
    return product.multiply(HypergeometricJS.hypergeometric0F1(n.add(1), sqrX.multiply(.25)));
}