diff --git a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/NumberTheory.java b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/NumberTheory.java index c2fd76c46f..44110935f1 100644 --- a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/NumberTheory.java +++ b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/NumberTheory.java @@ -133,20 +133,13 @@ private static class BellB extends AbstractFunctionEvaluator { * @param index an int number >= 0 * @return */ - private static IInteger bellNumber(int index) { + private static IExpr bellNumber(int index) { if (index < BELLB_25.length) { return AbstractIntegerSym.valueOf(BELLB_25[index]); } // Sum[StirlingS2[n, k], {k, 0, n}] - IInteger sum = F.C1; - for (int ki = 0; ki < index; ki++) { - sum = sum.add(stirlingS2(index, F.ZZ(ki), ki)); - if (sum.bitLength() > Config.MAX_BIT_LENGTH / 100) { - BigIntegerLimitExceeded.throwIt(Config.MAX_BIT_LENGTH / 100); - } - } - return sum; + return F.sum(k -> stirlingS2(index, k, k.toIntDefault()), 0, index, 1); } /** @@ -169,20 +162,21 @@ private static IExpr bellBPolynomial(int n, IExpr z) { return z; } - IASTAppendable sum = F.PlusAlloc(n + 1); - for (int k = 0; k <= n; k++) { - sum.append(F.Times(F.StirlingS2(F.ZZ(n), F.ZZ(k)), F.Power(z, k))); - } - return sum; + return F.sum(k -> F.Times(F.StirlingS2(F.ZZ(n), k), F.Power(z, k)), 0, n + 1, 1); } @Override public IExpr evaluate(final IAST ast, EvalEngine engine) { try { IExpr arg1 = ast.arg1(); + if (arg1.isNegative()) { + // Non-negative machine-sized integer expected at position `2` in `1` + return Errors.printMessage(S.BellB, "intnm", F.list(ast, F.C1), engine); + } int n = arg1.toIntDefault(); if (n < 0) { if (arg1.isNumber()) { + // Non-negative machine-sized integer expected at position `2` in `1` return Errors.printMessage(S.BellB, "intnm", F.list(ast, F.C1), engine); } } @@ -265,53 +259,48 @@ private static class BernoulliB extends AbstractFunctionEvaluator { @Override public IExpr evaluate(IAST ast, EvalEngine engine) { + if (ast.arg1().isNegative()) { + // Non-negative machine-sized integer expected at position `2` in `1`. + return Errors.printMessage(S.BernoulliB, "intnm", F.List(ast, F.C1), engine); + } if (ast.isAST1()) { - try { - int bn = ast.arg1().toIntDefault(); - if (bn >= 0) { - return bernoulliNumber(bn); - } - IExpr temp = engine.evaluate(F.Subtract(ast.arg1(), F.C3)); - if (temp.isIntegerResult() && temp.isPositiveResult() && temp.isEvenResult()) { - // http://fungrim.org/entry/a98234/ - return F.C0; - } - - } catch (RuntimeException rex) { - Errors.printMessage(S.BernoulliB, rex, engine); + int bn = ast.arg1().toIntDefault(); + if (bn >= 0) { + return bernoulliNumber(bn); + } + IExpr temp = engine.evaluate(F.Subtract(ast.arg1(), F.C3)); + if (temp.isIntegerResult() && temp.isPositiveResult() && temp.isEvenResult()) { + // http://fungrim.org/entry/a98234/ + return F.C0; } return F.NIL; } if (ast.isAST2()) { - try { - IExpr n = ast.arg1(); - IExpr x = ast.arg2(); - int xInt = x.toIntDefault(); - if (xInt != Integer.MIN_VALUE) { - if (xInt == 0) { - // http://fungrim.org/entry/a1d2d7/ - return F.BernoulliB(ast.arg1()); - } - if (xInt == 1 && n.isIntegerResult()) { - // http://fungrim.org/entry/829185/ - return F.Times(F.Power(F.CN1, n), F.BernoulliB(n)); - } + IExpr n = ast.arg1(); + IExpr x = ast.arg2(); + int xInt = x.toIntDefault(); + if (xInt != Integer.MIN_VALUE) { + if (xInt == 0) { + // http://fungrim.org/entry/a1d2d7/ + return F.BernoulliB(ast.arg1()); } - if (n.isInteger() && n.isNonNegativeResult()) { - if (x.isNumEqualRational(F.C1D2)) { - // http://fungrim.org/entry/03ee0b/ - return F.Times(F.Subtract(F.Power(F.C2, F.Subtract(F.C1, n)), F.C1), F.BernoulliB(n)); - } - int bn = n.toIntDefault(); - if (bn >= 0) { - // http://fungrim.org/entry/555e10/ - return F.sum( - k -> F.Times(F.Binomial(n, k), F.BernoulliB(F.Subtract(n, k)), F.Power(x, k)), 0, - bn); - } + if (xInt == 1 && n.isIntegerResult()) { + // http://fungrim.org/entry/829185/ + return F.Times(F.Power(F.CN1, n), F.BernoulliB(n)); + } + } + if (n.isInteger() && n.isNonNegativeResult()) { + if (x.isNumEqualRational(F.C1D2)) { + // http://fungrim.org/entry/03ee0b/ + return F.Times(F.Subtract(F.Power(F.C2, F.Subtract(F.C1, n)), F.C1), F.BernoulliB(n)); + } + int bn = n.toIntDefault(); + if (bn >= 0) { + // http://fungrim.org/entry/555e10/ + return F.sum( + k -> F.Times(F.Binomial(n, k), F.BernoulliB(F.Subtract(n, k)), F.Power(x, k)), 0, + bn); } - } catch (RuntimeException rex) { - Errors.printMessage(S.BernoulliB, rex, engine); } } return F.NIL; @@ -1709,10 +1698,6 @@ private static IExpr divisorSigma(IExpr arg1, IInteger n) { // general formula IASTAppendable sum = F.PlusAlloc(size); return sum.appendArgs(size, i -> F.Power(list.get(i), arg1)); - // for (int i = 1; i < size; i++) { - // sum.append(F.Power(list.get(i), arg1)); - // } - // return sum; } return F.NIL; } @@ -5089,9 +5074,9 @@ public IExpr evaluate(final IAST ast, EvalEngine engine) { } // try to convert into a fractional number return rationalize(arg1, epsilon, useConvergenceMethod).orElse(arg1); - } catch (Exception ex) { + } catch (RuntimeException rex) { // ex.printStackTrace(); - Errors.printMessage(S.Rationalize, ex, engine); + Errors.printMessage(S.Rationalize, rex, engine); } return F.NIL; diff --git a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/PolynomialFunctions.java b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/PolynomialFunctions.java index 13c3145b60..dd8c9b94a0 100644 --- a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/PolynomialFunctions.java +++ b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/builtin/PolynomialFunctions.java @@ -1516,9 +1516,14 @@ public IExpr evaluate(final IAST ast, EvalEngine engine) { if (arg1.isInteger() && ast.arg2().isInteger()) { int n = arg1.toIntDefault(); int k = ast.arg2().toIntDefault(); - if (n < 0 || k < 0 || !ast.arg3().isList() || ast.arg3().isMatrix() != null) { + if (n < 0 || k < 0) { return F.NIL; } + IExpr list = ast.arg3().normal(false); + if (!list.isList()) { + return F.NIL; + } + if (n == 0 && k == 0) { return F.C1; } @@ -1533,7 +1538,7 @@ public IExpr evaluate(final IAST ast, EvalEngine engine) { } int max = n - k + 2; if (max >= 0) { - return bellY(n, k, (IAST) ast.arg3(), ast, engine); + return bellY(n, k, (IAST) list, ast, engine); } } @@ -2259,13 +2264,24 @@ private static IAST monomialListModulus(IExpr polynomial, List variablesL } } - public static IExpr bellY(int n, int k, IAST symbols, IAST ast, EvalEngine engine) { + /** + * Calculate the partial Bell polynomial recursively. + * + * @param n + * @param k + * @param symbols + * @param listOfVariables + * @param engine + * @return + * + */ + public static IExpr bellY(int n, int k, IAST symbols, IAST listOfVariables, EvalEngine engine) { final int recursionLimit = engine.getRecursionLimit(); try { if (recursionLimit > 0) { int counter = engine.incRecursionCounter(); if (counter > recursionLimit) { - RecursionLimitExceeded.throwIt(counter, ast); + RecursionLimitExceeded.throwIt(counter, listOfVariables); } } if (n == 0 && k == 0) { @@ -2280,11 +2296,11 @@ public static IExpr bellY(int n, int k, IAST symbols, IAST ast, EvalEngine engin int iterationLimit = engine.getIterationLimit(); if (iterationLimit >= 0 && iterationLimit <= max) { - IterationLimitExceeded.throwIt(max, ast); + IterationLimitExceeded.throwIt(max, listOfVariables); } for (int m = 1; m < max; m++) { if ((m < symbols.size()) && !symbols.get(m).isZero()) { - IExpr bellY = bellY(n - m, k - 1, symbols, ast, engine); + IExpr bellY = bellY(n - m, k - 1, symbols, listOfVariables, engine); if (bellY.isPlus()) { bellY = bellY.mapThread(F.Times(a, F.Slot1, symbols.get(m)), 2); } else { diff --git a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/expression/F.java b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/expression/F.java index ed010e9368..5fbb8ccdbf 100644 --- a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/expression/F.java +++ b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/expression/F.java @@ -109,6 +109,7 @@ import org.matheclipse.core.eval.EvalAttributes; import org.matheclipse.core.eval.EvalEngine; import org.matheclipse.core.eval.exception.ASTElementLimitExceeded; +import org.matheclipse.core.eval.exception.BigIntegerLimitExceeded; import org.matheclipse.core.eval.interfaces.AbstractCoreFunctionEvaluator; import org.matheclipse.core.eval.interfaces.ICoreFunctionEvaluator; import org.matheclipse.core.eval.util.BiIntFunction; @@ -8155,7 +8156,8 @@ public static IAST Product(final IExpr expr, final IExpr iterationSpecification) * @param to * @return */ - public static IExpr product(final Function function, final int from, final int to) { + public static IExpr product(final Function function, final int from, + final int to) { return intProduct(function, from, to, 1); } @@ -9565,7 +9567,8 @@ public static IRational sumRational(final IntFunction function, final * @param iMax * @return */ - public static IExpr sum(final Function function, final int iMin, final int iMax) { + public static IExpr sum(final Function function, final int iMin, + final int iMax) { return intSum(function, iMin, iMax, 1); } @@ -9580,8 +9583,11 @@ public static IExpr sum(final Function function, final int iMin, f * @param step * @return */ - public static IExpr intProduct(final Function function, final int from, + public static IExpr intProduct(final Function function, final int from, final int to, final int step) { + if (from > to && step > 0) { + return F.C1; + } IASTAppendable result = ast(S.Times, 15); long numberOfLeaves = 0; INumber number = F.C1; @@ -9592,6 +9598,10 @@ public static IExpr intProduct(final Function function, final int IExpr temp = engine.evaluate(function.apply(ZZ(i))); if (temp.isNumber()) { number = number.times((INumber) temp); + if (number instanceof IInteger // + && ((IInteger) number).bitLength() > Config.MAX_BIT_LENGTH / 100) { + BigIntegerLimitExceeded.throwIt(Config.MAX_BIT_LENGTH / 100); + } } else { numberOfLeaves += temp.leafCount() + 1; if (numberOfLeaves >= Config.MAX_AST_SIZE / 2) { @@ -9614,6 +9624,9 @@ public static IExpr intProduct(final Function function, final int * @return */ public static IExpr intSum(final IntFunction function, final int iMin, final int iMax) { + if (iMin > iMax) { + return F.C0; + } IASTAppendable result = ast(S.Plus, 15); int numberOfLeaves = 0; EvalEngine engine = EvalEngine.get(); @@ -9624,6 +9637,10 @@ public static IExpr intSum(final IntFunction function, final int iMin, fi IExpr temp = engine.evaluate(function.apply(i)); if (temp.isNumber()) { number = number.plus((INumber) temp); + if (number instanceof IInteger // + && ((IInteger) number).bitLength() > Config.MAX_BIT_LENGTH / 100) { + BigIntegerLimitExceeded.throwIt(Config.MAX_BIT_LENGTH / 100); + } } else { numberOfLeaves += temp.leafCount() + 1; if (numberOfLeaves >= Config.MAX_AST_SIZE / 2) { @@ -9647,7 +9664,7 @@ public static IExpr intSum(final IntFunction function, final int iMin, fi * @param step * @return */ - public static IExpr intSum(final Function function, final int from, final int to, + public static IExpr intSum(final Function function, final int from, final int to, final int step) { return intSum(function, from, to, step, false); } @@ -9663,8 +9680,11 @@ public static IExpr intSum(final Function function, final int from * @param expand expand {@link S#Plus}, {@link S#Times}, {@link S#Power} subexpressions * @return */ - public static IExpr intSum(final Function function, final int from, final int to, + public static IExpr intSum(final Function function, final int from, final int to, final int step, boolean expand) { + if (from > to && step > 0) { + return F.C0; + } IASTAppendable result = ast(S.Plus, 15); long numberOfLeaves = 0; INumber number = F.C0; @@ -9675,6 +9695,10 @@ public static IExpr intSum(final Function function, final int from IExpr temp = engine.evaluate(function.apply(ZZ(i))); if (temp.isNumber()) { number = number.plus((INumber) temp); + if (number instanceof IInteger // + && ((IInteger) number).bitLength() > Config.MAX_BIT_LENGTH / 100) { + BigIntegerLimitExceeded.throwIt(Config.MAX_BIT_LENGTH / 100); + } } else { numberOfLeaves += temp.leafCount() + 1; if (numberOfLeaves >= Config.MAX_AST_SIZE / 2) { @@ -9701,7 +9725,7 @@ public static IExpr intSum(final Function function, final int from * @param iStep * @return */ - public static IExpr sum(final Function function, final int iMin, final int iMax, + public static IExpr sum(final Function function, final int iMin, final int iMax, final int iStep) { return intSum(function, iMin, iMax, iStep); } @@ -9716,7 +9740,7 @@ public static IExpr sum(final Function function, final int iMin, f * @param expand expand {@link S#Plus}, {@link S#Times}, {@link S#Power} subexpressions * @return */ - public static IExpr sum(final Function function, final int iMin, final int iMax, + public static IExpr sum(final Function function, final int iMin, final int iMax, final int iStep, boolean expand) { return intSum(function, iMin, iMax, iStep, expand); } diff --git a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/interfaces/IAST.java b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/interfaces/IAST.java index 952529ae8e..580b175b2f 100644 --- a/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/interfaces/IAST.java +++ b/symja_android_library/matheclipse-core/src/main/java/org/matheclipse/core/interfaces/IAST.java @@ -1511,10 +1511,11 @@ default IAST most() { } /** - * Get the argument at index 1, if the size() == 2 or the complete ast if the - * size() > 2 (useful for ASTs with attribute OneIdentity for example for - * Plus[] you can call getOneIdentity(F.C0) or for Times[]) - * you can call getOneIdentity(F.C1). + * Return the argument at index 1, if the argSize() == 1. Or return the complete ast + * if the argSize() > 1 If the argSize() == 0 return + * defaultValue (useful for ASTs with attribute OneIdentity for example + * for Plus[] you can call getOneIdentity(F.C0) or for + * Times[]) you can call getOneIdentity(F.C1). * * @param defaultValue default value, if size() < 2. * @return @@ -1522,8 +1523,8 @@ default IAST most() { public IExpr oneIdentity(IExpr defaultValue); /** - * Get the argument at index 1, if the size() == 2 or the complete ast if the - * size() > 2 If the size() == 1 return 0. + * Return the argument at index 1, if the argSize() == 1. Or return the complete ast + * if the argSize() > 1 If the argSize() == 0 return 0. * * @return */ @@ -1532,8 +1533,8 @@ default IExpr oneIdentity0() { } /** - * Get the argument at index 1, if the size() == 2 or the complete ast if the - * size() > 2 If the size() == 1 return 1. + * Return the argument at index 1, if the argSize() == 1. Or return the complete ast + * if the argSize() > 1 If the argSize() == 0 return 1. * * @return */ diff --git a/symja_android_library/matheclipse-core/src/test/java/org/matheclipse/core/system/LowercaseTestCase.java b/symja_android_library/matheclipse-core/src/test/java/org/matheclipse/core/system/LowercaseTestCase.java index 7a47963097..79c60bad3e 100644 --- a/symja_android_library/matheclipse-core/src/test/java/org/matheclipse/core/system/LowercaseTestCase.java +++ b/symja_android_library/matheclipse-core/src/test/java/org/matheclipse/core/system/LowercaseTestCase.java @@ -1566,6 +1566,8 @@ public void testBellB() { "BellB(n)"); check("BellB(5,x)", // "x+15*x^2+25*x^3+10*x^4+x^5"); + check("BellB(5,x+y^2)", // + "x+y^2+15*(x+y^2)^2+25*(x+y^2)^3+10*(x+y^2)^4+(x+y^2)^5"); check("Table(BellB(k), {k, 0, 14})", // "{1,1,2,5,15,52,203,877,4140,21147,115975,678570,4213597,27644437,190899322}"); check("BellB(10)", // @@ -1593,6 +1595,8 @@ public void testBellY() { check("BellY(2,1,{1/2,0})", // "0"); + check("BellY(5, 2, {})", // + "0"); check("BellY(5, 2, {1})", // "0"); check("BellY(5, 2, {1,2})", // @@ -1616,11 +1620,21 @@ public void testBellY() { public void testBernoulliB() { // slow + check("BernoulliB(2009,-1+Sqrt(2))", // - "BernoulliB(2009,-1+Sqrt(2))"); + "Hold(BernoulliB(2009,-1+Sqrt(2)))"); + // message: Non-negative machine-sized integer expected at position 1 in + // BernoulliB(-2147483648,1/2). check("BernoulliB(-2147483648,1/2)", // "BernoulliB(-2147483648,1/2)"); + // message: Non-negative machine-sized integer expected at position 1 in + // BernoulliB(-3). + check("BernoulliB(-3)", // + "BernoulliB(-3)"); + check("BernoulliB(3,-2)", // + "-15"); + check("BernoulliB(4, 9)", // "155519/30"); check("BernoulliB(4, -9)", //