use of org.matheclipse.core.eval.exception.ArgumentTypeException 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);
}
use of org.matheclipse.core.eval.exception.ArgumentTypeException in project symja_android_library by axkr.
the class Expr2LP method expr2ObjectiveFunction.
private ISignedNumber expr2ObjectiveFunction(final IExpr expr, double[] coefficients) throws ArithmeticException, ClassCastException {
if (expr instanceof IAST) {
final IAST ast = (IAST) expr;
if (ast.isPlus()) {
double constantTerm = 0.0;
for (int i = 1; i < ast.size(); i++) {
IExpr temp = ast.get(i);
ISignedNumber num = expr2ObjectiveFunction(temp, coefficients);
if (num != null) {
constantTerm += num.doubleValue();
}
}
return F.num(constantTerm);
} else if (ast.isTimes()) {
ISymbol variable = null;
double value = 1.0;
for (int i = 1; i < ast.size(); i++) {
IExpr temp = ast.get(i);
if (temp.isVariable()) {
if (variable != null) {
throw new ArgumentTypeException("conversion from expression to linear programming expression failed for " + temp.toString());
}
variable = (ISymbol) temp;
continue;
}
ISignedNumber num = temp.evalReal();
if (num != null) {
value *= num.doubleValue();
continue;
}
throw new ArgumentTypeException("conversion from expression to linear programming expression failed for " + temp.toString());
}
if (variable != null) {
for (int i = 0; i < coefficients.length; i++) {
if (variable.equals(fVariables.get(i))) {
coefficients[i] += value;
return null;
}
}
throw new ArgumentTypeException("conversion from expression to linear programming expression failed for " + ast.toString());
}
return F.num(value);
}
} else if (expr.isVariable()) {
ISymbol variable = (ISymbol) expr;
for (int i = 0; i < coefficients.length; i++) {
if (variable.equals(fVariables.get(i))) {
coefficients[i] += 1.0d;
return null;
}
}
throw new ArgumentTypeException("conversion from expression to linear programming expression failed for " + expr.toString());
}
ISignedNumber num = expr.evalReal();
if (num != null) {
return num;
}
throw new ArgumentTypeException("conversion from expression to linear programming expression failed for " + expr.toString());
}
use of org.matheclipse.core.eval.exception.ArgumentTypeException in project symja_android_library by axkr.
the class BesselJS method besselJZero.
/**
* @param n a positive value
* @param m
* @return
*/
public static double besselJZero(double n, int m) {
if (n < 0.0) {
throw new ArgumentTypeException(n + " is less than 0 (negative order)");
}
// approximations from dlmf.nist.gov/10.21#vi
double delta = Math.PI / 4.0;
double a = (m + n / 2.0 - 0.25) * Math.PI;
double e = a - (4.0 * (n * n) - 1.0) / (8.0 * a);
BisectionSolver solver = new BisectionSolver();
ISymbol x = F.Dummy("x");
UnivariateDifferentiableFunction f = new UnaryNumerical(F.BesselJ(F.num(n), x), x, EvalEngine.get(), true);
return solver.solve(100, f, e - delta, e + delta);
}
use of org.matheclipse.core.eval.exception.ArgumentTypeException in project symja_android_library by axkr.
the class BesselJS method besselYDouble.
public static double besselYDouble(double n, double x) {
// double delta = 1e-5;
if (x < 0) {
throw new ArgumentTypeException(x + " < 0.0");
}
// dlmf.nist.gov/10.2#E3
if (F.isNumIntValue(n)) {
FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(15, 0.01);
ISymbol nSymbol = F.Dummy("n");
UnivariateDifferentiableFunction f = differentiator.differentiate(new UnaryNumerical(F.BesselJ(nSymbol, F.num(x)), nSymbol, EvalEngine.get()));
DSFactory factory = new DSFactory(1, 1);
double d1 = f.value(factory.variable(0, n)).getPartialDerivative(1);
double d2 = f.value(factory.variable(0, -n)).getPartialDerivative(1);
return (d1 + Math.pow(-1.0, Math.rint(n)) * d2) / Math.PI;
// return ( diff( n => besselJ(n,x), n ) + (-1)**n * diff( n => besselJ(n,x), -n ) ) / pi;
// old without differentiation
// return (besselYDouble(n + delta, x) + besselYDouble(n - delta, x)) / 2.0;
}
return (besselJDouble(n, x) * Math.cos(n * Math.PI) - besselJDouble(-n, x)) / Math.sin(n * Math.PI);
}
use of org.matheclipse.core.eval.exception.ArgumentTypeException in project symja_android_library by axkr.
the class SparseArrayExpr method toRealVector.
/**
* {@inheritDoc}
*/
@Override
public RealVector toRealVector() {
if (fDimension.length == 1 && fDimension[0] > 0) {
try {
OpenMapRealVector result = new OpenMapRealVector(fDimension[0]);
if (!fDefaultValue.isZero()) {
double d = fDefaultValue.evalDouble();
for (int i = 0; i < fDimension[0]; i++) {
result.setEntry(i, d);
}
}
for (TrieNode<int[], IExpr> entry : fData.nodeSet()) {
int[] key = entry.getKey();
IExpr value = entry.getValue();
result.setEntry(key[0] - 1, value.evalDouble());
}
return result;
} catch (ArgumentTypeException rex) {
}
}
return null;
}
Aggregations