Java Code Examples for org.apache.commons.math3.complex.Complex#divide()

private static ComplexPair transform(final Complex in, final double a, final double b) {
    if (in.isInfinite()) {
        return new ComplexPair(new Complex(-1), new Complex(1));
    }

    final Complex c = new Complex(1).add(in).divide(new Complex(1).subtract(in)); // bilinear

    final double a2 = a * a;
    final double b2 = b * b;
    final double ab = a * b;
    final double ab2 = 2 * ab;
    Complex v = new Complex(0).add(c.multiply(4 * (b2 * (a2 - 1) + 1)));
    v = v.add(8 * (b2 * (a2 - 1) - 1));
    v = v.multiply(c);
    v = v.add(4 * (b2 * (a2 - 1) + 1));
    v = v.sqrt();

    final Complex u = v.multiply(-1).add(c.multiply(ab2)).add(ab2);

    v = v.add(c.multiply(ab2)).add(ab2);

    final Complex d = new Complex(0).add(c.multiply(2 * (b - 1))).add(2 * (1 + b));

    return new ComplexPair(u.divide(d), v.divide(d));
}
 

public Complex response(final double normalizedFrequency) {
    final double w = 2 * Math.PI * normalizedFrequency;
    final Complex czn1 = ComplexUtils.polar2Complex(1., -w);
    final Complex czn2 = ComplexUtils.polar2Complex(1., -2 * w);
    Complex ch = new Complex(1);
    Complex cbot = new Complex(1);

    for (int i = 0; i < mNumBiquads; i++) {
        final Biquad stage = mBiquads[i];

        Complex ct = new Complex(stage.getB0() / stage.getA0()); // NOPMD
        ct = ct.add(czn1.multiply(stage.getB1() / stage.getA0()));
        ct = ct.add(czn2.multiply(stage.getB2() / stage.getA0()));

        Complex cb = new Complex(1); // NOPMD
        cb = cb.add(czn1.multiply(stage.getA1() / stage.getA0()));
        cb = cb.add(czn2.multiply(stage.getA2() / stage.getA0()));

        ch = ch.multiply(ct);
        cbot = cbot.multiply(cb);
    }

    return ch.divide(cbot);
}
 

public SV otherSide(double r, double x, double g, double b, double ratio) {
    Complex z = new Complex(r, x); // z=r+jx
    Complex y = new Complex(g, b); // y=g+jb
    Complex s1 = new Complex(p, q); // s1=p1+jq1
    Complex u1 = ComplexUtils.polar2Complex(u, Math.toRadians(a));
    Complex v1 = u1.divide(Math.sqrt(3f)); // v1=u1/sqrt(3)

    Complex v1p = v1.multiply(ratio); // v1p=v1*rho
    Complex i1 = s1.divide(v1.multiply(3)).conjugate(); // i1=conj(s1/(3*v1))
    Complex i1p = i1.divide(ratio); // i1p=i1/rho
    Complex i2 = i1p.subtract(y.multiply(v1p)); // i2=i1p-y*v1p
    Complex v2 = v1p.subtract(z.multiply(i2)); // v2=v1p-z*i2
    Complex s2 = v2.multiply(3).multiply(i2.conjugate()); // s2=3*v2*conj(i2)

    Complex u2 = v2.multiply(Math.sqrt(3f));
    return new SV(-s2.getReal(), -s2.getImaginary(), u2.abs(), Math.toDegrees(u2.getArgument()));
}
 

public SV otherSide(double r, double x, double g1, double b1, double g2, double b2, double ratio) {
    Complex z = new Complex(r, x); // z=r+jx
    Complex y1 = new Complex(g1, b1); // y1=g1+jb1
    Complex y2 = new Complex(g2, b2); // y2=g2+jb2
    Complex s1 = new Complex(p, q); // s1=p1+jq1
    Complex u1 = ComplexUtils.polar2Complex(u, Math.toRadians(a));
    Complex v1 = u1.divide(Math.sqrt(3f)); // v1=u1/sqrt(3)

    Complex v1p = v1.multiply(ratio); // v1p=v1*rho
    Complex i1 = s1.divide(v1.multiply(3)).conjugate(); // i1=conj(s1/(3*v1))
    Complex i1p = i1.divide(ratio); // i1p=i1/rho
    Complex i2p = i1p.subtract(y1.multiply(v1p)); // i2p=i1p-y1*v1p
    Complex v2 = v1p.subtract(z.multiply(i2p)); // v2p=v1p-z*i2
    Complex i2 = i2p.subtract(y2.multiply(v2)); // i2=i2p-y2*v2
    Complex s2 = v2.multiply(3).multiply(i2.conjugate()); // s2=3*v2*conj(i2)

    Complex u2 = v2.multiply(Math.sqrt(3f));
    return new SV(-s2.getReal(), -s2.getImaginary(), u2.abs(), Math.toDegrees(u2.getArgument()));
}
 

private ComplexPair transform(Complex c) {
	if (c.isInfinite()) {
		return new ComplexPair(new Complex(-1), new Complex(1));
	}

	c = ((new Complex(1)).add(c)).divide((new Complex(1)).subtract(c)); // bilinear

	Complex v = new Complex(0);
	v = MathSupplement.addmul(v, 4 * (b2 * (a2 - 1) + 1), c);
	v = v.add(8 * (b2 * (a2 - 1) - 1));
	v = v.multiply(c);
	v = v.add(4 * (b2 * (a2 - 1) + 1));
	v = v.sqrt();

	Complex u = v.multiply(-1);
	u = MathSupplement.addmul(u, ab_2, c);
	u = u.add(ab_2);

	v = MathSupplement.addmul(v, ab_2, c);
	v = v.add(ab_2);

	Complex d = new Complex(0);
	d = MathSupplement.addmul(d, 2 * (b - 1), c).add(2 * (1 + b));

	return new ComplexPair(u.divide(d), v.divide(d));
}
 

public Complex response(double normalizedFrequency) {
	double w = 2 * Math.PI * normalizedFrequency;
	Complex czn1 = ComplexUtils.polar2Complex(1., -w);
	Complex czn2 = ComplexUtils.polar2Complex(1., -2 * w);
	Complex ch = new Complex(1);
	Complex cbot = new Complex(1);

	for (int i = 0; i < m_numBiquads; i++) {
		Biquad stage = m_biquads[i];
		Complex cb = new Complex(1);
		Complex ct = new Complex(stage.getB0() / stage.getA0());
		ct = MathSupplement.addmul(ct, stage.getB1() / stage.getA0(), czn1);
		ct = MathSupplement.addmul(ct, stage.getB2() / stage.getA0(), czn2);
		cb = MathSupplement.addmul(cb, stage.getA1() / stage.getA0(), czn1);
		cb = MathSupplement.addmul(cb, stage.getA2() / stage.getA0(), czn2);
		ch = ch.multiply(ct);
		cbot = cbot.multiply(cb);
	}

	return ch.divide(cbot);
}
 

public Complex response(double normalizedFrequency) {
    double a0 = getA0();
    double a1 = getA1();
    double a2 = getA2();
    double b0 = getB0();
    double b1 = getB1();
    double b2 = getB2();

    double w = 2 * Math.PI * normalizedFrequency;
    Complex czn1 = ComplexUtils.polar2Complex(1., -w);
    Complex czn2 = ComplexUtils.polar2Complex(1., -2 * w);
    Complex ch = new Complex(1);
    Complex cbot = new Complex(1);

    Complex ct = new Complex(b0 / a0);
    Complex cb = new Complex(1);
    ct = MathSupplement.addmul(ct, b1 / a0, czn1);
    ct = MathSupplement.addmul(ct, b2 / a0, czn2);
    cb = MathSupplement.addmul(cb, a1 / a0, czn1);
    cb = MathSupplement.addmul(cb, a2 / a0, czn2);
    ch = ch.multiply(ct);
    cbot = cbot.multiply(cb);

    return ch.divide(cbot);
}
 

public Complex response(final double normalizedFrequency) {
    final double a0 = getA0();
    final double a1 = getA1();
    final double a2 = getA2();
    final double b0 = getB0();
    final double b1 = getB1();
    final double b2 = getB2();

    final double w = 2 * Math.PI * normalizedFrequency;
    final Complex czn1 = ComplexUtils.polar2Complex(1., -w);
    final Complex czn2 = ComplexUtils.polar2Complex(1., -2 * w);
    Complex ch = new Complex(1);
    Complex cbot = new Complex(1);

    Complex ct = new Complex(b0 / a0);

    ct = ct.add(czn1.multiply(b1 / a0));
    ct = ct.add(czn2.multiply(b2 / a0));

    Complex cb = new Complex(1);
    cb = cb.add(czn1.multiply(a1 / a0));
    cb = cb.add(czn2.multiply(a2 / a0));

    ch = ch.multiply(ct);
    cbot = cbot.multiply(cb);

    return ch.divide(cbot);
}
 

public StateVariable toSv2(StateVariable sv1) {
    Complex s1 = new Complex(-sv1.p, -sv1.q); // s1=p1+jq1
    Complex u1 = ComplexUtils.polar2Complex(sv1.u, Math.toRadians(sv1.theta));
    Complex v1 = u1.divide(SQUARE_3); // v1=u1/sqrt(3)
    Complex v1p = v1.multiply(ratio); // v1p=v1*rho
    Complex i1 = s1.divide(v1.multiply(3)).conjugate(); // i1=conj(s1/(3*v1))
    Complex i1p = i1.divide(ratio); // i1p=i1/rho
    Complex i2 = i1p.subtract(y.multiply(v1p)).negate(); // i2=-(i1p-y*v1p)
    Complex v2 = v1p.subtract(z.multiply(i2)); // v2=v1p-z*i2
    Complex s2 = v2.multiply(3).multiply(i2.conjugate()); // s2=3*v2*conj(i2)
    Complex u2 = v2.multiply(SQUARE_3);
    return new StateVariable(-s2.getReal(), -s2.getImaginary(), u2.abs(), Math.toDegrees(u2.getArgument()));
}
 

public StateVariable toSv1(StateVariable sv2) {
    Complex s2 = new Complex(-sv2.p, -sv2.q); // s2=p2+jq2
    Complex u2 = ComplexUtils.polar2Complex(sv2.u, Math.toRadians(sv2.theta));
    Complex v2 = u2.divide(SQUARE_3); // v2=u2/sqrt(3)
    Complex i2 = s2.divide(v2.multiply(3)).conjugate(); // i2=conj(s2/(3*v2))
    Complex v1p = v2.add(z.multiply(i2)); // v1'=v2+z*i2
    Complex i1p = i2.negate().add(y.multiply(v1p)); // i1'=-i2+v1'*y
    Complex i1 = i1p.multiply(ratio); // i1=i1p*ration
    Complex v1 = v1p.divide(ratio); // v1=v1p/ration
    Complex s1 = v1.multiply(3).multiply(i1.conjugate()); // s1=3*v1*conj(i1)
    Complex u1 = v1.multiply(SQUARE_3);
    return new StateVariable(-s1.getReal(), -s1.getImaginary(), u1.abs(), Math.toDegrees(u1.getArgument()));
}
 

@Override
public Complex apply(final Complex argument) {
	final Complex iargument = argument.multiply(Complex.I);
	final Complex exponent = iargument.add(1.0);
	final Complex numerator = (new Complex(strike)).pow(exponent.subtract(1.0)).multiply(exponent);
	final Complex denominator = (argument.multiply(argument)).subtract(iargument);

	return numerator.divide(denominator);
}
 

@Override
public Complex apply(final Complex argument) {
	final Complex iargument = argument.multiply(Complex.I);
	final Complex exponent = iargument.add(1.0);
	final Complex numerator = (new Complex(strike)).pow(exponent.subtract(1.0)).multiply(exponent);
	final Complex denominator = (argument.multiply(argument)).subtract(iargument);

	return numerator.divide(denominator);
}
 

@Override
public Map<String, Function<Double, Double>> getValue(final double evaluationTime, final CharacteristicFunctionModel model) throws CalculationException {

	final CharacteristicFunction modelCF = model.apply(getMaturity());

	final double lineOfIntegration = 0.5 * (getIntegrationDomainImagUpperBound()+getIntegrationDomainImagLowerBound());

	final double lambda = 2*Math.PI/(numberOfPoints*gridSpacing); //Equation 23 Carr and Madan
	final double upperBound = (numberOfPoints * lambda)/2.0; //Equation 20 Carr and Madan

	final Complex[] integrandEvaluations = new Complex[numberOfPoints];

	for(int i = 0; i<numberOfPoints; i++) {

		final double u = gridSpacing * i;

		//Integration over a line parallel to the real axis
		final Complex z = new Complex(u,-lineOfIntegration);

		//The characteristic function is already discounted
		final Complex numerator = modelCF.apply(z.subtract(Complex.I));

		final Complex denominator = apply(z);
		Complex ratio = numerator.divide(denominator);
		ratio = (ratio.multiply(((Complex.I).multiply(upperBound*u)).exp())).multiply(gridSpacing);

		double delta;
		if (i==0){
			delta=1.0;
		}else{
			delta = 0.0;
		}
		final double simpsonWeight = (3+Math.pow(-1,i+1)-delta)/3;

		integrandEvaluations[i] = ratio.multiply(simpsonWeight);
	}

	//Compute the FFT
	Complex[] transformedVector = new Complex[numberOfPoints];
	final FastFourierTransformer fft=new FastFourierTransformer(DftNormalization.STANDARD);
	transformedVector=	fft.transform(integrandEvaluations,TransformType.FORWARD);

	//Find relevant prices via interpolation
	final double[] logStrikeVector = new double[numberOfPoints];
	final double[] strikeVector = new double[numberOfPoints];
	final double[] optionPriceVector = new double[numberOfPoints];

	for(int j = 0; j<numberOfPoints; j++) {
		logStrikeVector[j] = -upperBound+lambda*j;
		strikeVector[j] = Math.exp(logStrikeVector[j]);
		optionPriceVector[j] = (transformedVector[j].multiply(Math.exp(-lineOfIntegration * logStrikeVector[j]))).getReal()/Math.PI;
	}

	final RationalFunctionInterpolation interpolation = new RationalFunctionInterpolation(strikeVector, optionPriceVector,intMethod, extMethod);

	final Complex minusI = new Complex(0,-1);
	final double residueTerm = (modelCF.apply(minusI)).getReal();

	final Function<Double, Double> strikeToPrice = new Function<Double, Double>(){

		@Override
		public Double apply(final Double t) {
			return residueTerm + interpolation.getValue(t);
		}

	};

	final HashMap<String, Function<Double, Double>> results = new HashMap<>();
	results.put("valuePerStrike", strikeToPrice);
	return results;
}
 

@Override
public Map<String, Function<Double, Double>> getValue(final double evaluationTime, final CharacteristicFunctionModel model) throws CalculationException {

	final CharacteristicFunction modelCF = model.apply(getMaturity());

	final double lineOfIntegration = 0.5 * (getIntegrationDomainImagUpperBound()+getIntegrationDomainImagLowerBound());

	final double lambda = 2*Math.PI/(numberOfPoints*gridSpacing); //Equation 23 Carr and Madan
	final double upperBound = (numberOfPoints * lambda)/2.0; //Equation 20 Carr and Madan

	final Complex[] integrandEvaluations = new Complex[numberOfPoints];

	for(int i = 0; i<numberOfPoints; i++) {

		final double u = gridSpacing * i;

		//Integration over a line parallel to the real axis
		final Complex z = new Complex(u,-lineOfIntegration);

		//The characteristic function is already discounted
		final Complex numerator = modelCF.apply(z.subtract(Complex.I));

		final Complex denominator = apply(z);
		Complex ratio = numerator.divide(denominator);
		ratio = (ratio.multiply(((Complex.I).multiply(upperBound*u)).exp())).multiply(gridSpacing);

		double delta;
		if (i==0){
			delta=1.0;
		}else{
			delta = 0.0;
		}
		final double simpsonWeight = (3+Math.pow(-1,i+1)-delta)/3;

		integrandEvaluations[i] = ratio.multiply(simpsonWeight);
	}

	//Compute the FFT
	Complex[] transformedVector = new Complex[numberOfPoints];
	final FastFourierTransformer fft=new FastFourierTransformer(DftNormalization.STANDARD);
	transformedVector=	fft.transform(integrandEvaluations,TransformType.FORWARD);

	//Find relevant prices via interpolation
	final double[] logStrikeVector = new double[numberOfPoints];
	final double[] strikeVector = new double[numberOfPoints];
	final double[] optionPriceVector = new double[numberOfPoints];

	for(int j = 0; j<numberOfPoints; j++) {
		logStrikeVector[j] = -upperBound+lambda*j;
		strikeVector[j] = Math.exp(logStrikeVector[j]);
		optionPriceVector[j] = (transformedVector[j].multiply(Math.exp(-lineOfIntegration * logStrikeVector[j]))).getReal()/Math.PI;
	}

	final RationalFunctionInterpolation interpolation = new RationalFunctionInterpolation(strikeVector, optionPriceVector,intMethod, extMethod);

	final Complex minusI = new Complex(0,-1);
	final double residueTerm = (modelCF.apply(minusI)).getReal();

	final Function<Double, Double> strikeToPrice = new Function<Double, Double>(){

		@Override
		public Double apply(final Double t) {
			return residueTerm + interpolation.getValue(t);
		}

	};

	final HashMap<String, Function<Double, Double>> results = new HashMap<>();
	results.put("valuePerStrike", strikeToPrice);
	return results;
}
 

public Complex characteristicFunctionByMonteCarlo(final Complex zeta, final RandomVariable processAtTime) {

		final int states = processAtTime.getRealizations().length;

		Complex runningSum = new Complex(0.0,0.0);

		for(int i = 0; i< states; i++) {
			runningSum = runningSum.add((Complex.I.multiply(zeta.multiply(processAtTime.get(i)))).exp());
		}

		return runningSum.divide(states);
	}
 

public Complex characteristicFunctionByMonteCarlo(final Complex zeta, final RandomVariable processAtTime) {

		final int states = processAtTime.getRealizations().length;

		Complex runningSum = new Complex(0.0,0.0);

		for(int i = 0; i< states; i++) {
			runningSum = runningSum.add((Complex.I.multiply(zeta.multiply(processAtTime.get(i)))).exp());
		}

		return runningSum.divide(states);
	}