All Downloads are FREE. Search and download functionalities are using the official Maven repository.

cc.redberry.core.number.Complex Maven / Gradle / Ivy

Go to download

Redberry is an open source computer algebra system designed for tensor manipulation. It implements basic computer algebra system routines as well as complex tools for real computations in physics. This is Redberry core, which contains the implementation of the basic computer algebra routines and general-purpose transformations.

There is a newer version: 1.1.10
Show newest version
/*
 * Redberry: symbolic tensor computations.
 *
 * Copyright (c) 2010-2014:
 *   Stanislav Poslavsky   
 *   Bolotin Dmitriy       
 *
 * This file is part of Redberry.
 *
 * Redberry is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * Redberry is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with Redberry. If not, see .
 */
package cc.redberry.core.number;

import cc.redberry.core.context.OutputFormat;
import cc.redberry.core.indices.Indices;
import cc.redberry.core.indices.IndicesFactory;
import cc.redberry.core.tensor.*;
import org.apache.commons.math3.Field;
import org.apache.commons.math3.exception.NullArgumentException;
import org.apache.commons.math3.fraction.BigFraction;
import org.apache.commons.math3.util.FastMath;

import java.io.Serializable;
import java.math.BigInteger;

/**
 * Implementation of complex numbers.
 *
 * @author Stanislav Poslavsky
 */
public final class Complex extends Tensor
        implements Number,
        Serializable {

    public static final Complex ComplexNaN =
            new Complex(Numeric.NaN, Numeric.NaN);

    public static final Complex REAL_POSITIVE_INFINITY =
            new Complex(Numeric.POSITIVE_INFINITY, Numeric.ZERO);

    public static final Complex REAL_NEGATIVE_INFINITY =
            new Complex(Numeric.NEGATIVE_INFINITY, Numeric.ZERO);

    public static final Complex IMAGINARY_POSITIVE_INFINITY =
            new Complex(Numeric.ZERO, Numeric.POSITIVE_INFINITY);

    public static final Complex IMAGINARY_NEGATIVE_INFINITY =
            new Complex(Numeric.ZERO, Numeric.NEGATIVE_INFINITY);

    public static final Complex COMPLEX_NEGATIVE_INFINITY =
            new Complex(Numeric.NEGATIVE_INFINITY, Numeric.NEGATIVE_INFINITY);

    public static final Complex COMPLEX_POSITIVE_INFINITY =
            new Complex(Numeric.POSITIVE_INFINITY, Numeric.POSITIVE_INFINITY);

    public static final Complex COMPLEX_INFINITY = COMPLEX_POSITIVE_INFINITY;

    public static final Complex ZERO =
            new Complex(Rational.ZERO, Rational.ZERO);

    public static final Complex ONE =
            new Complex(Rational.ONE, Rational.ZERO);

    public static final Complex ONE_HALF =
            new Complex(Rational.ONE_HALF, Rational.ZERO);

    public static final Complex TWO =
            new Complex(Rational.TWO, Rational.ZERO);

    public static final Complex FOUR =
            new Complex(Rational.FOUR, Rational.ZERO);

    public static final Complex MINUS_ONE =
            new Complex(Rational.MINUS_ONE, Rational.ZERO);

    public static final Complex MINUS_ONE_HALF =
            new Complex(Rational.MINUSE_ONE_HALF, Rational.ZERO);

    public static final Complex MINUS_TWO =
            new Complex(Rational.MINUS_TWO, Rational.ZERO);

    public static final Complex IMAGINARY_UNIT =
            new Complex(Rational.ZERO, Rational.ONE);

    public static final Complex NEGATIVE_IMAGINARY_UNIT =
            new Complex(Rational.ZERO, Rational.MINUS_ONE);

    private final Real real;
    private final Real imaginary;

    public Complex(Real real, Real imaginary) {
        if (real instanceof Numeric || imaginary instanceof Numeric) {
            this.real = real.getNumericValue();
            this.imaginary = imaginary.getNumericValue();
        } else {
            this.real = real;
            this.imaginary = imaginary;
        }
    }

    public Complex(Real real) {
        if (real instanceof Numeric) {
            this.real = real.getNumericValue();
            this.imaginary = Numeric.ZERO;
        } else {
            this.real = real;
            this.imaginary = Rational.ZERO;
        }
    }

    public Complex(org.apache.commons.math3.complex.Complex complex) {
        this(complex.getReal(), complex.getImaginary());
    }

    public Complex(int real, int imaginary) {
        this(new Rational(real), new Rational(imaginary));
    }

    public Complex(int real) {
        this(new Rational(real), Rational.ZERO);
    }

    public Complex(double real, double imaginary) {
        this(new Numeric(real), new Numeric(imaginary));
    }

    public Complex(int real, double imaginary) {
        this(new Numeric(real), new Numeric(imaginary));
    }

    public Complex(double real, int imaginary) {
        this(new Numeric(real), new Numeric(imaginary));
    }

    public Complex(double real) {
        this(new Numeric(real), Numeric.ZERO);
    }

    public Complex(long real) {
        this(new Rational(real), Rational.ZERO);
    }

    public Complex(long real, long imaginary) {
        this(new Rational(real), new Rational(imaginary));
    }

    public Complex(BigInteger real, BigInteger imaginary) {
        this(new Rational(real), new Rational(imaginary));
    }

    public Complex(BigInteger real) {
        this(new Rational(real), Rational.ZERO);
    }

    @Override
    public Tensor get(int i) {
        throw new IndexOutOfBoundsException();
    }

    @Override
    public Indices getIndices() {
        return IndicesFactory.EMPTY_INDICES;
    }

    @Override
    protected int hash() {
        return 47 * (329 + real.hashCode()) + this.imaginary.hashCode();
    }

    @Override
    public int size() {
        return 0;
    }

    @Override
    public String toString(OutputFormat mode) {
        if (real.isZero())
            if (imaginary.isZero())
                return "0";
            else if (imaginary.isOne())
                return "I";
            else if (imaginary.isMinusOne())
                return "-I";
            else
                return imaginary.toString() + "*I";
        int is = imaginary.signum();
        if (is == 0)
            return real.toString();
        Real abs = imaginary.abs();
        if (is < 0)
            if (abs.isOne())
                return real.toString() + "-I";
            else
                return real.toString() + "-I*" + imaginary.abs();
        if (abs.isOne())
            return real.toString() + "+I";
        else
            return real.toString() + "+I*" + imaginary.abs();

    }

    @Override
    protected String toString(OutputFormat mode, Class clazz) {
        if (clazz == Product.class || clazz == Power.class) {
            if (!imaginary.isZero() || real.signum() < 0 || !real.isInteger())
                return "(" + toString(mode) + ")";
        }
        return toString(mode);
    }

    @Override
    public TensorBuilder getBuilder() {
        return new ComplexBuilder(this);
    }

    @Override
    public TensorFactory getFactory() {
        return null;
    }

    private static class ComplexBuilder implements TensorBuilder {

        private final Complex complex;

        public ComplexBuilder(Complex complex) {
            this.complex = complex;
        }

        @Override
        public Tensor build() {
            return complex;
        }

        @Override
        public void put(Tensor tensor) {
            throw new IllegalStateException("Can not put to Complex tensor builder!");
        }

        @Override
        public TensorBuilder clone() {
            return this;
        }
    }

    public Real getImaginary() {
        return imaginary;
    }

    public Real getReal() {
        return real;
    }

    public Complex getImaginaryPart() {
        return new Complex(imaginary.getField().getZero(), imaginary);
    }

    public Complex getRealPart() {
        return new Complex(real);
    }

    /**
     * Returns {@code true} if this number is real.
     *
     * @return {@code true} if this number is real
     */
    public boolean isReal() {
        return imaginary.isZero();
    }

    /**
     * Returns {@code true} if this number has only imaginary part.
     *
     * @return {@code true} if this number has only imaginary part
     */
    public boolean isImaginary() {
        return real.isZero();
    }

    public Complex getImaginaryAsComplex() {
        if (imaginary.isOne())
            return ONE;
        if (imaginary.isZero())
            return ZERO;
        return new Complex(Rational.ZERO, imaginary);
    }

    public Complex getRealAsComplex() {
        if (real.isOne())
            return ONE;
        if (real.isZero())
            return ZERO;
        return new Complex(real);
    }

    @Override
    public boolean isNaN() {
        return real.isNaN() || imaginary.isNaN();
    }

    @Override
    public boolean isNumeric() {
        return real.isNumeric();
    }

    @Override
    public boolean isInfinite() {
        return real.isInfinite() || imaginary.isInfinite();
    }

    @Override
    public boolean isZero() {
        return real.isZero() && imaginary.isZero();
    }

    @Override
    public boolean isOne() {
        return real.isOne() && imaginary.isZero();
    }

    @Override
    public boolean isMinusOne() {
        return real.isMinusOne() && imaginary.isZero();
    }

    public boolean isOneOrMinusOne() {
        return imaginary.isZero() && (real.isOne() || real.isMinusOne());
    }

    /**
     * Returns double value of the real part.
     *
     * @return double value of the real part
     */
    @Override
    public double doubleValue() {
        return real.doubleValue();
    }

    /**
     * Returns float value of the real part.
     *
     * @return float value of the real part
     */
    @Override
    public float floatValue() {
        return real.floatValue();
    }

    /**
     * Returns int value of the real part.
     *
     * @return int value of the real part
     */
    @Override
    public int intValue() {
        return real.intValue();
    }

    @Override
    public BigInteger bigIntValue() {
        return real.bigIntValue();
    }

    /**
     * Returns long value of the real part.
     *
     * @return long value of the real part
     */
    @Override
    public long longValue() {
        return real.longValue();
    }

    /**
     * Return the conjugate of this complex number. The conjugate of {@code a + bi} is {@code a - bi}. 
{@link * #ComplexNaN} is returned if either the real or imaginary part of this Complex number equals {@code Double.NaN}. *
If the imaginary part is infinite, and the real part is not {@code NaN}, the returned value has infinite * imaginary part of the opposite sign, e.g. the conjugate of {@code 1 + POSITIVE_INFINITY i} is {@code 1 - * NEGATIVE_INFINITY i}. * * @return the conjugate of this Complex object. */ public Complex conjugate() { return new Complex(real, imaginary.negate()); } /** * Returns a {@code Complex} whose value is {@code (this + addend)}. Uses the definitional formula *
     * 
     *   (a + bi) + (c + di) = (a+c) + (b+d)i
     * 
     * 

If either {@code this} or {@code addend} has a {@code NaN} * value in either part, {@link #ComplexNaN} is returned; otherwise {@code Infinite} and {@code NaN} values are * returned in the parts of the result according to the rules for {@link java.lang.Double} arithmetic. * * @param a value to be added to this {@code Complex}. * @return {@code this + addend}. * @throws NullArgumentException if {@code addend} is {@code null}. */ @Override public Complex add(Complex a) { NumberUtils.checkNotNull(a); return a.isZero() ? a.isNumeric() ? this.getNumericValue() : this : new Complex(real.add(a.real), imaginary.add(a.imaginary)); } /** * Returns a {@code Complex} whose value is {@code (this / divisor)}. Implements the definitional formula *
     * 
     *    a + bi          ac + bd + (bc - ad)i
     *    ----------- = -------------------------
     *    c + di         c2 + d2
     * 
     * 
but uses * prescaling of operands to limit the effects of overflows and underflows in the computation.
{@code * Infinite} and {@code NaN} values are handled according to the following rules, applied in the order presented: *
  • If either {@code this} or {@code divisor} has a {@code NaN} value in either part, {@link #ComplexNaN} is * returned.
  • If {@code divisor} equals {@link #ZERO}, {@link #ComplexNaN} is returned.
  • If {@code * this} and {@code divisor} are both infinite, {@link #ComplexNaN} is returned.
  • If {@code this} is finite * (i.e., has no {@code Infinite} or {@code NaN} parts) and {@code divisor} is infinite (one or both parts * infinite), {@link #ZERO} is returned.
  • If {@code this} is infinite and {@code divisor} is finite, {@code * NaN} values are returned in the parts of the result if the {@link java.lang.Double} rules applied to the * definitional formula force {@code NaN} results.
* * @param divisor Value by which this {@code Complex} is to be divided. * @return {@code this / divisor}. * @throws NullArgumentException if {@code divisor} is {@code null}. */ @Override public Complex divide(Complex divisor) { NumberUtils.checkNotNull(divisor); if (divisor.isOne()) return divisor.isNumeric() ? this.getNumericValue() : this; if (divisor.isNaN() || isNaN()) return ComplexNaN; final Real c = divisor.real; final Real d = divisor.imaginary; // Real denominator = (c.multiply(c)).add(d.multiply(d)); // return new Complex(((real.multiply(c)).add(imaginary.multiply(d))).divide(denominator), // ((imaginary.multiply(c)).subtract(real.multiply(d))).divide(denominator)); if (c.abs().compareTo(d.abs()) < 0) { Real q = c.divide(d); Real denominator = c.multiply(q).add(d); return new Complex((real.multiply(q).add(imaginary)).divide(denominator), (imaginary.multiply(q).subtract(real)).divide(denominator)); } else { Real q = d.divide(c); Real denominator = d.multiply(q).add(c); return new Complex(((imaginary.multiply(q)).add(real)).divide(denominator), (imaginary.subtract(real.multiply(q))).divide(denominator)); } } @Override public Field getField() { return ComplexField.getInstance(); } @Override public Complex multiply(int n) { return n == 1 ? this : new Complex(real.multiply(n), imaginary.multiply(n)); } /** * Returns a {@code Complex} whose value is {@code this * factor}. Implements preliminary checks for {@code NaN} and * infinity followed by the definitional formula: *
     * 
     *   (a + bi)(c + di) = (ac - bd) + (ad + bc)i
     * 
     * 
Returns {@link #ComplexNaN} if either {@code this} or {@code factor} has * one or more {@code NaN} parts.
Returns {@link #COMPLEX_INFINITY} if neither {@code this} nor {@code factor} * has one or more {@code NaN} parts and if either {@code this} or {@code factor} has one or more infinite parts * (same result is returned regardless of the sign of the components).
Returns finite values in components of * the result per the definitional formula in all remaining cases. * * @param factor value to be multiplied by this {@code Complex}. * @return {@code this * factor}. * @throws NullArgumentException if {@code factor} is {@code null}. */ @Override public Complex multiply(Complex factor) { NumberUtils.checkNotNull(factor); if (factor.isNaN()) return ComplexNaN; return new Complex(real.multiply(factor.real).subtract(imaginary.multiply(factor.imaginary)), real.multiply(factor.imaginary).add(imaginary.multiply(factor.real))); } @Override public Complex negate() { return new Complex(real.negate(), imaginary.negate()); } @Override public Complex reciprocal() { if (isNaN()) return ComplexNaN; // if (FastMath.abs(real) < FastMath.abs(imaginary)) { // double q = real / imaginary; // double scale = 1. / (real * q + imaginary); // return createComplex(scale * q, -scale); // } else { // double q = imaginary / real; // double scale = 1. / (imaginary * q + real); // return createComplex(scale, -scale * q); // } if (real.abs().compareTo(imaginary.abs()) < 0) { Real q = real.divide(imaginary); Real scale = (real.multiply(q).add(imaginary)).reciprocal(); return new Complex(scale.multiply(q), scale.negate()); } else { Real q = imaginary.divide(real); Real scale = (imaginary.multiply(q).add(real)).reciprocal(); return new Complex(scale, scale.multiply(q).negate()); } } @Override public Complex subtract(Complex a) { NumberUtils.checkNotNull(a); return a.isZero() ? a.isNumeric() ? this.getNumericValue() : this : new Complex(real.subtract(a.real), imaginary.subtract(a.imaginary)); } @Override public Complex getNumericValue() { return isNumeric() ? this : new Complex(real.getNumericValue(), imaginary.getNumericValue()); } @Override public Complex abs() { if (isZero() || isOne() || isInfinite() || isNaN()) return this; Real abs2 = real.multiply(real).add(imaginary.multiply(imaginary)); if (isNumeric()) return new Complex(abs2.pow(0.5)); Rational abs2r = (Rational) abs2; BigInteger num = abs2r.getNumerator(); BigInteger den = abs2r.getDenominator(); BigInteger nR = NumberUtils.sqrt(num); if (!NumberUtils.isSqrt(num, nR)) throw new IllegalStateException(); BigInteger dR = NumberUtils.sqrt(den); if (!NumberUtils.isSqrt(den, dR)) throw new IllegalStateException(); return new Complex(new Rational(nR, dR)); } public double absNumeric() { if (isNaN()) return Double.NaN; if (isInfinite()) return Double.POSITIVE_INFINITY; final double real = this.real.doubleValue(); final double imaginary = this.imaginary.doubleValue(); return absNumeric(real, imaginary); } public static double absNumeric(final double real, final double imaginary) { if (FastMath.abs(real) < FastMath.abs(imaginary)) { if (imaginary == 0.0) { return FastMath.abs(real); } double q = real / imaginary; return FastMath.abs(imaginary) * FastMath.sqrt(1 + q * q); } else { if (real == 0.0) { return FastMath.abs(imaginary); } double q = imaginary / real; return FastMath.abs(real) * FastMath.sqrt(1 + q * q); } } @Override public Complex add(BigFraction fraction) { return fraction.compareTo(BigFraction.ZERO) == 0 ? this : new Complex(real.add(fraction), imaginary); } @Override public Complex add(double d) { return d == 0.0 ? this.getNumericValue() : new Complex(real.add(d), imaginary); } @Override public Complex add(long d) { return d == 0 ? this : new Complex(real.add(d), imaginary); } @Override public Complex add(int d) { return d == 0 ? this : new Complex(real.add(d), imaginary); } @Override public Complex add(BigInteger bg) { NumberUtils.checkNotNull(bg); return bg.equals(BigInteger.ZERO) ? this : new Complex(real.add(bg), imaginary); } @Override public Complex subtract(BigFraction fraction) { return fraction.compareTo(BigFraction.ZERO) == 0 ? this : new Complex(real.subtract(fraction), imaginary); } @Override public Complex subtract(double d) { return d == 0.0 ? this.getNumericValue() : new Complex(real.subtract(d), imaginary); } @Override public Complex subtract(long l) { return l == 0 ? this : new Complex(real.subtract(l), imaginary); } @Override public Complex subtract(int i) { return i == 0 ? this : new Complex(real.subtract(i), imaginary); } @Override public Complex subtract(BigInteger bg) { NumberUtils.checkNotNull(bg); return bg.equals(BigInteger.ZERO) ? this : new Complex(real.subtract(bg), imaginary); } @Override public Complex multiply(BigFraction fraction) { return fraction.compareTo(BigFraction.ONE) == 0 ? this : new Complex(real.multiply(fraction), imaginary.multiply(fraction)); } @Override public Complex multiply(double d) { return d == 1.0 ? this.getNumericValue() : Double.isNaN(d) ? ComplexNaN : new Complex(real.multiply(d), imaginary.multiply(d)); } @Override public Complex multiply(BigInteger bg) { return bg.compareTo(BigInteger.ONE) == 0 ? this : new Complex(real.multiply(bg), imaginary.multiply(bg)); } @Override public Complex multiply(long d) { return d == 1 ? this : new Complex(real.multiply(d), imaginary.multiply(d)); } @Override public Complex divide(BigFraction fraction) { return fraction.compareTo(BigFraction.ONE) == 0 ? this : new Complex(real.divide(fraction), imaginary.divide(fraction)); } @Override public Complex divide(long l) { return l == 1 ? this : new Complex(real.divide(l), imaginary.divide(l)); } @Override public Complex divide(int i) { return i == 1 ? this : new Complex(real.divide(i), imaginary.divide(i)); } @Override public Complex divide(BigInteger bg) { return bg.compareTo(BigInteger.ONE) == 0 ? this : new Complex(real.divide(bg), imaginary.divide(bg)); } @Override public Complex divide(double d) { return d == 1.0 ? this : new Complex(real.divide(d), imaginary.divide(d)); } public Complex add(Real d) { NumberUtils.checkNotNull(d); return d.isZero() ? (d.isNumeric() ? this.getNumericValue() : this) : new Complex(real.add(d), imaginary.add(d)); } public Complex subtract(Real d) { NumberUtils.checkNotNull(d); return d.isZero() ? (d.isNumeric() ? this.getNumericValue() : this) : new Complex(real.subtract(d), imaginary.add(d)); } public Complex multiply(Real d) { NumberUtils.checkNotNull(d); return d.isOne() ? (d.isNumeric() ? this.getNumericValue() : this) : new Complex(real.multiply(d), imaginary.multiply(d)); } public Complex divide(Real d) { NumberUtils.checkNotNull(d); return d.isOne() ? (d.isNumeric() ? this.getNumericValue() : this) : new Complex(real.divide(d), imaginary.divide(d)); } @Override public boolean equals(Object obj) { if (obj == null) return false; if (getClass() != obj.getClass()) return false; final Complex other = (Complex) obj; return real.equals(other.real) && imaginary.equals(other.imaginary); } /** * Returns natural logarithm of this * * @return natural logarithm of this */ public Complex logNumeric() { if (isNaN()) { return ComplexNaN; } final double real = this.real.doubleValue(); final double imaginary = this.imaginary.doubleValue(); return new Complex(FastMath.log(absNumeric(real, imaginary)), FastMath.atan2(imaginary, real)); } @Override public Complex pow(double exponent) { return this.logNumeric().multiply(exponent).expNumeric(); } @Override public Complex pow(BigInteger exponent) { if (exponent.compareTo(BigInteger.ZERO) < 0) return reciprocal().pow(exponent.negate()); Complex result = Complex.ONE; Complex k2p = this; while (!BigInteger.ZERO.equals(exponent)) { if (exponent.testBit(0)) result = result.multiply(k2p); k2p = k2p.multiply(k2p); exponent = exponent.shiftRight(1); } return result; } @Override public Complex pow(long exponent) { if (exponent < 0) return reciprocal().pow(-exponent); Complex result = Complex.ONE; Complex k2p = this; while (exponent != 0) { if (exponent != 0) result = result.multiply(k2p); k2p = k2p.multiply(k2p); exponent = exponent >> 1; } return result; } @Override public Complex pow(int exponent) { if (exponent < 0) return reciprocal().pow(-exponent); Complex result = Complex.ONE; Complex k2p = this; while (exponent != 0) { if ((exponent & 0x1) != 0) result = result.multiply(k2p); k2p = k2p.multiply(k2p); exponent = exponent >> 1; } return result; } public Complex expNumeric() { if (isNaN()) { return ComplexNaN; } final double real = this.real.doubleValue(); final double imaginary = this.imaginary.doubleValue(); double expReal = FastMath.exp(real); return new Complex(expReal * FastMath.cos(imaginary), expReal * FastMath.sin(imaginary)); } public Complex powNumeric(Complex exponent) { return this.logNumeric().multiply(exponent).expNumeric(); } @Override public boolean isInteger() { return imaginary.isZero() && real.isInteger(); } @Override public boolean isNatural() { return imaginary.isZero() && real.isNatural(); } public boolean isNegativeNatural() { return imaginary.isZero() && real.isInteger() && real.signum() < 0; } public boolean isPositiveNatural() { return imaginary.isZero() && real.isInteger() && real.signum() >= 0; } }




© 2015 - 2025 Weber Informatics LLC | Privacy Policy