Package org.apache.commons.math.complex

Examples of org.apache.commons.math.complex.Complex


        int iterationCount = 0;
        if (n < 1) {
            throw MathRuntimeException.createIllegalArgumentException(
                  NON_POSITIVE_DEGREE_MESSAGE, n);
        }
        Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
        for (int i = 0; i <= n; i++) {
            c[i] = coefficients[i];
        }

        // solve individual root successively
        Complex root[] = new Complex[n];
        for (int i = 0; i < n; i++) {
            Complex subarray[] = new Complex[n-i+1];
            System.arraycopy(c, 0, subarray, 0, subarray.length);
            root[i] = solve(subarray, initial);
            // polynomial deflation using synthetic division
            Complex newc = c[n-i];
            Complex oldc = null;
            for (int j = n-i-1; j >= 0; j--) {
                oldc = c[j];
                c[j] = newc;
                newc = oldc.add(newc.multiply(root[i]));
            }
            iterationCount += this.iterationCount;
        }

        resultComputed = true;
View Full Code Here


        int n = coefficients.length - 1;
        if (n < 1) {
            throw MathRuntimeException.createIllegalArgumentException(
                  NON_POSITIVE_DEGREE_MESSAGE, n);
        }
        Complex N  = new Complex(n,     0.0);
        Complex N1 = new Complex(n - 1, 0.0);

        int i = 1;
        Complex pv = null;
        Complex dv = null;
        Complex d2v = null;
        Complex G = null;
        Complex G2 = null;
        Complex H = null;
        Complex delta = null;
        Complex denominator = null;
        Complex z = initial;
        Complex oldz = new Complex(Double.POSITIVE_INFINITY, Double.POSITIVE_INFINITY);
        while (i <= maximalIterationCount) {
            // Compute pv (polynomial value), dv (derivative value), and
            // d2v (second derivative value) simultaneously.
            pv = coefficients[n];
            dv = Complex.ZERO;
            d2v = Complex.ZERO;
            for (int j = n-1; j >= 0; j--) {
                d2v = dv.add(z.multiply(d2v));
                dv = pv.add(z.multiply(dv));
                pv = coefficients[j].add(z.multiply(pv));
            }
            d2v = d2v.multiply(new Complex(2.0, 0.0));

            // check for convergence
            double tolerance = Math.max(relativeAccuracy * z.abs(),
                                        absoluteAccuracy);
            if ((z.subtract(oldz)).abs() <= tolerance) {
                resultComputed = true;
                iterationCount = i;
                return z;
            }
            if (pv.abs() <= functionValueAccuracy) {
                resultComputed = true;
                iterationCount = i;
                return z;
            }

            // now pv != 0, calculate the new approximation
            G = dv.divide(pv);
            G2 = G.multiply(G);
            H = G2.subtract(d2v.divide(pv));
            delta = N1.multiply((N.multiply(H)).subtract(G2));
            // choose a denominator larger in magnitude
            Complex deltaSqrt = delta.sqrt();
            Complex dplus = G.add(deltaSqrt);
            Complex dminus = G.subtract(deltaSqrt);
            denominator = dplus.abs() > dminus.abs() ? dplus : dminus;
            // Perturb z if denominator is zero, for instance,
            // p(x) = x^3 + 1, z = 0.
            if (denominator.equals(new Complex(0.0, 0.0))) {
                z = z.add(new Complex(absoluteAccuracy, absoluteAccuracy));
                oldz = new Complex(Double.POSITIVE_INFINITY,
                                   Double.POSITIVE_INFINITY);
            } else {
                oldz = z;
                z = z.subtract(N.divide(denominator));
            }
View Full Code Here

            final double b = 0.5 * (f[i] - f[n-i]);
            x[i]     = a + b;
            x[n - i] = a - b;
        }
        FastFourierTransformer transformer = new FastFourierTransformer();
        Complex y[] = transformer.transform(x);

        // reconstruct the FST result for the original array
        transformed[0] = 0.0;
        transformed[1] = 0.5 * y[0].getReal();
        for (int i = 1; i < (n >> 1); i++) {
View Full Code Here

            x[i] = a - b;
            x[n-i] = a + b;
            t1 += c;
        }
        FastFourierTransformer transformer = new FastFourierTransformer();
        Complex y[] = transformer.transform(x);

        // reconstruct the FCT result for the original array
        transformed[0] = y[0].getReal();
        transformed[1] = t1;
        for (int i = 1; i < (n >> 1); i++) {
View Full Code Here

     */
    protected Complex[] fft(double f[], boolean isInverse)
        throws IllegalArgumentException {

        verifyDataSet(f);
        Complex F[] = new Complex[f.length];
        if (f.length == 1) {
            F[0] = new Complex(f[0], 0.0);
            return F;
        }

        // Rather than the naive real to complex conversion, pack 2N
        // real numbers into N complex numbers for better performance.
        int N = f.length >> 1;
        Complex c[] = new Complex[N];
        for (int i = 0; i < N; i++) {
            c[i] = new Complex(f[2*i], f[2*i+1]);
        }
        roots.computeOmega(isInverse ? -N : N);
        Complex z[] = fft(c);

        // reconstruct the FFT result for the original array
        roots.computeOmega(isInverse ? -2*N : 2*N);
        F[0] = new Complex(2 * (z[0].getReal() + z[0].getImaginary()), 0.0);
        F[N] = new Complex(2 * (z[0].getReal() - z[0].getImaginary()), 0.0);
        for (int i = 1; i < N; i++) {
            Complex A = z[N-i].conjugate();
            Complex B = z[i].add(A);
            Complex C = z[i].subtract(A);
            //Complex D = roots.getOmega(i).multiply(Complex.I);
            Complex D = new Complex(-roots.getOmegaImaginary(i),
                                    roots.getOmegaReal(i));
            F[i] = B.subtract(C.multiply(D));
            F[2*N-i] = F[i].conjugate();
        }

View Full Code Here

     */
    protected Complex[] fft(Complex data[])
        throws IllegalArgumentException {

        final int n = data.length;
        final Complex f[] = new Complex[n];

        // initial simple cases
        verifyDataSet(data);
        if (n == 1) {
            f[0] = data[0];
            return f;
        }
        if (n == 2) {
            f[0] = data[0].add(data[1]);
            f[1] = data[0].subtract(data[1]);
            return f;
        }

        // permute original data array in bit-reversal order
        int ii = 0;
        for (int i = 0; i < n; i++) {
            f[i] = data[ii];
            int k = n >> 1;
            while (ii >= k && k > 0) {
                ii -= k; k >>= 1;
            }
            ii += k;
        }

        // the bottom base-4 round
        for (int i = 0; i < n; i += 4) {
            final Complex a = f[i].add(f[i+1]);
            final Complex b = f[i+2].add(f[i+3]);
            final Complex c = f[i].subtract(f[i+1]);
            final Complex d = f[i+2].subtract(f[i+3]);
            final Complex e1 = c.add(d.multiply(Complex.I));
            final Complex e2 = c.subtract(d.multiply(Complex.I));
            f[i] = a.add(b);
            f[i+2] = a.subtract(b);
            // omegaCount indicates forward or inverse transform
            f[i+1] = roots.isForward() ? e2 : e1;
            f[i+3] = roots.isForward() ? e1 : e2;
        }

        // iterations from bottom to top take O(N*logN) time
        for (int i = 4; i < n; i <<= 1) {
            final int m = n / (i<<1);
            for (int j = 0; j < n; j += i<<1) {
                for (int k = 0; k < i; k++) {
                    //z = f[i+j+k].multiply(roots.getOmega(k*m));
                    final int k_times_m = k*m;
                    final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                    final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                    //z = f[i+j+k].multiply(omega[k*m]);
                    final Complex z = new Complex(
                        f[i+j+k].getReal() * omega_k_times_m_real -
                        f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                        f[i+j+k].getReal() * omega_k_times_m_imaginary +
                        f[i+j+k].getImaginary() * omega_k_times_m_real);

 
View Full Code Here

     * @param d the real scaling coefficient
     * @return a reference to the scaled array
     */
    public static Complex[] scaleArray(Complex f[], double d) {
        for (int i = 0; i < f.length; i++) {
            f[i] = new Complex(d * f[i].getReal(), d * f[i].getImaginary());
        }
        return f;
    }
View Full Code Here

            Object[] lastDimension = (Object[]) multiDimensionalComplexArray;
            for (int i = 0; i < dimensionSize.length - 1; i++) {
                lastDimension = (Object[]) lastDimension[vector[i]];
            }

            Complex lastValue = (Complex) lastDimension[vector[dimensionSize.length - 1]];
            lastDimension[vector[dimensionSize.length - 1]] = magnitude;

            return lastValue;
        }
View Full Code Here

    /**
     * Test of solver for the quintic function using solveAll().
     */
    public void testQuinticFunction2() throws MathException {
        double initial = 0.0, tolerance;
        Complex expected, result[];

        // p(x) = x^5 + 4x^3 + x^2 + 4 = (x+1)(x^2-x+1)(x^2+4)
        double coefficients[] = { 4.0, 0.0, 1.0, 4.0, 0.0, 1.0 };
        PolynomialFunction f = new PolynomialFunction(coefficients);
        LaguerreSolver solver = new LaguerreSolver(f);
        result = solver.solveAll(coefficients, initial);

        expected = new Complex(0.0, -2.0);
        tolerance = Math.max(solver.getAbsoluteAccuracy(),
                    Math.abs(expected.abs() * solver.getRelativeAccuracy()));
        TestUtils.assertContains(result, expected, tolerance);

        expected = new Complex(0.0, 2.0);
        tolerance = Math.max(solver.getAbsoluteAccuracy(),
                    Math.abs(expected.abs() * solver.getRelativeAccuracy()));
        TestUtils.assertContains(result, expected, tolerance);

        expected = new Complex(0.5, 0.5 * Math.sqrt(3.0));
        tolerance = Math.max(solver.getAbsoluteAccuracy(),
                    Math.abs(expected.abs() * solver.getRelativeAccuracy()));
        TestUtils.assertContains(result, expected, tolerance);

        expected = new Complex(-1.0, 0.0);
        tolerance = Math.max(solver.getAbsoluteAccuracy(),
                    Math.abs(expected.abs() * solver.getRelativeAccuracy()));
        TestUtils.assertContains(result, expected, tolerance);
       
        expected = new Complex(0.5, -0.5 * Math.sqrt(3.0));
        tolerance = Math.max(solver.getAbsoluteAccuracy(),
                    Math.abs(expected.abs() * solver.getRelativeAccuracy()));
        TestUtils.assertContains(result, expected, tolerance);
    }
View Full Code Here

        if (p.value(min) == 0.0) { return min; }
        if (p.value(max) == 0.0) { return max; }
        verifyBracketing(min, max, p);

        double coefficients[] = p.getCoefficients();
        Complex c[] = new Complex[coefficients.length];
        for (int i = 0; i < coefficients.length; i++) {
            c[i] = new Complex(coefficients[i], 0.0);
        }
        Complex initial = new Complex(0.5 * (min + max), 0.0);
        Complex z = solve(c, initial);
        if (isRootOK(min, max, z)) {
            setResult(z.getReal(), iterationCount);
            return result;
        }

        // solve all roots and select the one we're seeking
        Complex[] root = solveAll(c, initial);
View Full Code Here

TOP

Related Classes of org.apache.commons.math.complex.Complex

Copyright © 2018 www.massapicom. All rights reserved.
All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.