[poincare] Arithmetic::PrimeFactorization return the number of prime

factors found
This commit is contained in:
Émilie Feral
2018-09-20 13:53:04 +02:00
parent 9ae6a3a9d2
commit 0e1fd1bad7
5 changed files with 22 additions and 35 deletions

View File

@@ -19,7 +19,7 @@ public:
* small integer optimization. Thereby, the tables will not allocate any
* Integer on the static Integer table that limit the number of Integer
* simultaneously alive. */
static void PrimeFactorization(const Integer & i, Integer outputFactors[], Integer outputCoefficients[], int outputLength);
static int PrimeFactorization(const Integer & i, Integer outputFactors[], Integer outputCoefficients[], int outputLength);
constexpr static int k_numberOfPrimeFactors = 1000;
constexpr static int k_maxNumberOfPrimeFactors = 32;
private:

View File

@@ -40,7 +40,7 @@ int primeFactors[Arithmetic::k_numberOfPrimeFactors] = {2, 3, 5, 7, 11, 13, 17,
3659, 3671, 3673, 3677, 3691, 3697, 3701, 3709, 3719, 3727, 3733, 3739, 3761, 3767, 3769, 3779, 3793, 3797, 3803, 3821, 3823, 3833, 3847, 3851, 3853, 3863, 3877, 3881, 3889, 3907, 3911, 3917, 3919, 3923, 3929, 3931, 3943, 3947, 3967, 3989, 4001, 4003, 4007, 4013, 4019, 4021, 4027, 4049, 4051, 4057, 4073, 4079, 4091, 4093, 4099, 4111, 4127, 4129, 4133, 4139, 4153, 4157, 4159, 4177, 4201, 4211, 4217, 4219, 4229, 4231, 4241, 4243, 4253, 4259, 4261, 4271, 4273, 4283, 4289, 4297, 4327, 4337, 4339, 4349, 4357, 4363, 4373, 4391, 4397, 4409, 4421, 4423, 4441, 4447, 4451, 4457, 4463, 4481, 4483, 4493, 4507, 4513, 4517, 4519, 4523, 4547, 4549, 4561, 4567, 4583, 4591, 4597, 4603, 4621, 4637, 4639, 4643, 4649, 4651, 4657, 4663, 4673, 4679, 4691, 4703, 4721, 4723, 4729, 4733, 4751, 4759, 4783, 4787, 4789, 4793, 4799, 4801, 4813, 4817, 4831, 4861, 4871, 4877, 4889, 4903, 4909, 4919, 4931, 4933, 4937, 4943, 4951, 4957, 4967, 4969, 4973, 4987, 4993, 4999, 5003, 5009, 5011, 5021, 5023, 5039, 5051, 5059, 5077, 5081, 5087, 5099, 5101, 5107, 5113, 5119, 5147, 5153, 5167, 5171, 5179, 5189, 5197, 5209, 5227, 5231, 5233, 5237, 5261, 5273, 5279, 5281, 5297, 5303, 5309, 5323, 5333, 5347, 5351, 5381, 5387, 5393, 5399, 5407, 5413, 5417, 5419, 5431, 5437, 5441, 5443, 5449, 5471, 5477, 5479, 5483, 5501, 5503, 5507, 5519, 5521, 5527, 5531, 5557, 5563, 5569, 5573, 5581, 5591, 5623, 5639, 5641, 5647, 5651, 5653, 5657, 5659, 5669, 5683, 5689, 5693, 5701, 5711, 5717, 5737, 5741, 5743, 5749, 5779, 5783, 5791, 5801, 5807, 5813, 5821, 5827, 5839, 5843, 5849, 5851, 5857, 5861, 5867, 5869, 5879, 5881, 5897, 5903, 5923, 5927, 5939, 5953, 5981, 5987, 6007, 6011, 6029, 6037, 6043, 6047, 6053, 6067, 6073, 6079, 6089, 6091, 6101, 6113, 6121, 6131, 6133, 6143, 6151, 6163, 6173, 6197, 6199, 6203, 6211, 6217, 6221, 6229, 6247, 6257, 6263, 6269, 6271, 6277, 6287, 6299, 6301, 6311, 6317, 6323, 6329, 6337, 6343, 6353, 6359, 6361, 6367, 6373, 6379, 6389, 6397, 6421, 6427, 6449, 6451, 6469, 6473, 6481, 6491, 6521, 6529, 6547, 6551, 6553, 6563, 6569, 6571, 6577, 6581, 6599, 6607, 6619, 6637, 6653, 6659, 6661, 6673, 6679, 6689, 6691, 6701, 6703, 6709, 6719, 6733, 6737, 6761, 6763, 6779, 6781, 6791, 6793, 6803, 6823, 6827, 6829, 6833, 6841, 6857, 6863, 6869, 6871, 6883, 6899, 6907, 6911, 6917, 6947, 6949, 6959, 6961, 6967, 6971, 6977, 6983, 6991, 6997, 7001, 7013, 7019, 7027, 7039, 7043, 7057, 7069, 7079, 7103, 7109, 7121, 7127, 7129, 7151, 7159, 7177, 7187, 7193, 7207, 7211, 7213, 7219, 7229, 7237, 7243, 7247, 7253, 7283, 7297, 7307, 7309, 7321, 7331, 7333, 7349, 7351, 7369, 7393, 7411, 7417, 7433, 7451, 7457, 7459, 7477, 7481, 7487, 7489, 7499, 7507, 7517, 7523, 7529, 7537, 7541, 7547, 7549, 7559, 7561, 7573, 7577, 7583, 7589, 7591, 7603, 7607, 7621, 7639, 7643, 7649, 7669, 7673, 7681, 7687, 7691, 7699, 7703, 7717, 7723, 7727, 7741, 7753, 7757, 7759, 7789, 7793, 7817, 7823, 7829, 7841, 7853, 7867, 7873, 7877, 7879, 7883, 7901, 7907, 7919};
// we can go to 7907*7907 = 62 520 649
void Arithmetic::PrimeFactorization(const Integer & n, Integer outputFactors[], Integer outputCoefficients[], int outputLength) {
int Arithmetic::PrimeFactorization(const Integer & n, Integer outputFactors[], Integer outputCoefficients[], int outputLength) {
assert(!n.isInfinity());
// Compute the absolute value of n
@@ -51,8 +51,7 @@ void Arithmetic::PrimeFactorization(const Integer & n, Integer outputFactors[],
* the prime factorization for low numbers). When k_numberOfPrimeFactors is
* overflow, try every number as divisor. */
if (Integer::NaturalOrder(m, Integer(1)) == 0) {
outputCoefficients[0] = -1;
return;
return 0;
}
const char * primorial = "525896479052627740771371797072411912900610967452630";
const Integer primorial32(primorial, strlen(primorial), false);
@@ -60,11 +59,7 @@ void Arithmetic::PrimeFactorization(const Integer & n, Integer outputFactors[],
/* Special case 1: We do not want to break i in prime factor because it
* might take too many factors... More than k_maxNumberOfPrimeFactors.
* outputCoefficients[0] is set to -1 to indicate a special case. */
outputCoefficients[0] = -1;
return;
}
for (int index = 0; index < outputLength; index++) {
outputCoefficients[index] = Integer(0);
return -1;
}
int t = 0; // n prime factor index
@@ -80,7 +75,7 @@ void Arithmetic::PrimeFactorization(const Integer & n, Integer outputFactors[],
outputCoefficients[t] = Integer::Addition(outputCoefficients[t], Integer(1));
m = d.quotient;
if (m.isOne()) {
return;
return t+1;
}
continue;
}
@@ -96,11 +91,11 @@ void Arithmetic::PrimeFactorization(const Integer & n, Integer outputFactors[],
* take too much time: the prime factor that should be tested is above
* k_biggestPrimeFactor.
* outputCoefficients[0] is set to -1 to indicate a special case. */
outputCoefficients[0] = -1;
return;
return -2;
}
outputFactors[t] = m;
outputCoefficients[t] = Integer::Addition(outputCoefficients[t], Integer(1));
return t+1;
}
}

View File

@@ -62,25 +62,23 @@ Multiplication Factor::createMultiplicationOfIntegerPrimeDecomposition(Integer i
assert(!i.isZero());
assert(!i.isNegative());
Multiplication m;
if (i.isOne()) {
Integer factors[Arithmetic::k_maxNumberOfPrimeFactors];
Integer coefficients[Arithmetic::k_maxNumberOfPrimeFactors];
int numberOfPrimeFactors = Arithmetic::PrimeFactorization(i, factors, coefficients, Arithmetic::k_maxNumberOfPrimeFactors);
if (numberOfPrimeFactors == 0) {
m.addChildAtIndexInPlace(Rational(i), 0, 0);
return m;
}
Integer factors[Arithmetic::k_maxNumberOfPrimeFactors];
Integer coefficients[Arithmetic::k_maxNumberOfPrimeFactors];
Arithmetic::PrimeFactorization(i, factors, coefficients, Arithmetic::k_maxNumberOfPrimeFactors);
int index = 0;
if (coefficients[0].isMinusOne()) {
if (numberOfPrimeFactors < 0) {
// Exception: the decomposition failed
return m;
}
while (!coefficients[index].isZero() && index < Arithmetic::k_maxNumberOfPrimeFactors) {
for (int index = 0; index < numberOfPrimeFactors; index++) {
Expression factor = Rational(factors[index]);
if (!coefficients[index].isOne()) {
factor = Power(factor, Rational(coefficients[index]));
}
m.addChildAtIndexInPlace(factor, m.numberOfChildren(), m.numberOfChildren());
index++;
}
return m;
}

View File

@@ -194,14 +194,13 @@ bool Logarithm::parentIsAPowerOfSameBase() const {
Expression Logarithm::splitInteger(Integer i, bool isDenominator, Context & context, Preferences::AngleUnit angleUnit) {
assert(!i.isZero());
assert(!i.isNegative());
if (i.isOne()) {
return Rational(0);
}
assert(!i.isOne());
Integer factors[Arithmetic::k_maxNumberOfPrimeFactors];
Integer coefficients[Arithmetic::k_maxNumberOfPrimeFactors];
Arithmetic::PrimeFactorization(i, factors, coefficients, Arithmetic::k_maxNumberOfPrimeFactors);
if (coefficients[0].isMinusOne()) {
int numberOfPrimeFactors = Arithmetic::PrimeFactorization(i, factors, coefficients, Arithmetic::k_maxNumberOfPrimeFactors);
if (numberOfPrimeFactors == 0) {
return Rational(0);
}
if (numberOfPrimeFactors < 0) {
/* We could not break i in prime factor (either it might take too many
* factors or too much time). */
Expression e = clone();
@@ -213,8 +212,7 @@ Expression Logarithm::splitInteger(Integer i, bool isDenominator, Context & cont
return m;
}
Addition a;
int index = 0;
while (!coefficients[index].isZero() && index < Arithmetic::k_maxNumberOfPrimeFactors) {
for (int index = 0; index < numberOfPrimeFactors; index++) {
if (isDenominator) {
coefficients[index].setNegative(true);
}
@@ -224,7 +222,6 @@ Expression Logarithm::splitInteger(Integer i, bool isDenominator, Context & cont
e.simpleShallowReduce(context, angleUnit);
a.addChildAtIndexInPlace(m, a.numberOfChildren(), a.numberOfChildren());
m.shallowReduce(context, angleUnit);
index++;
}
return a;
}

View File

@@ -717,9 +717,8 @@ Expression Power::CreateSimplifiedIntegerRationalPower(Integer i, Rational r, bo
}
Integer factors[Arithmetic::k_maxNumberOfPrimeFactors];
Integer coefficients[Arithmetic::k_maxNumberOfPrimeFactors];
Arithmetic::PrimeFactorization(i, factors, coefficients, Arithmetic::k_maxNumberOfPrimeFactors);
if (coefficients[0].isMinusOne()) {
int numberOfPrimeFactors = Arithmetic::PrimeFactorization(i, factors, coefficients, Arithmetic::k_maxNumberOfPrimeFactors);
if (numberOfPrimeFactors <= 0) {
/* We could not break i in prime factors (it might take either too many
* factors or too much time). */
Expression rClone = r.clone().setSign(isDenominator ? ExpressionNode::Sign::Negative : ExpressionNode::Sign::Positive, context, angleUnit);
@@ -728,13 +727,11 @@ Expression Power::CreateSimplifiedIntegerRationalPower(Integer i, Rational r, bo
Integer r1(1);
Integer r2(1);
int index = 0;
while (!coefficients[index].isZero() && index < Arithmetic::k_maxNumberOfPrimeFactors) {
for (int index = 0; index < numberOfPrimeFactors; index++) {
Integer n = Integer::Multiplication(coefficients[index], r.signedIntegerNumerator());
IntegerDivision div = Integer::Division(n, r.integerDenominator());
r1 = Integer::Multiplication(r1, Integer::Power(factors[index], div.quotient));
r2 = Integer::Multiplication(r2, Integer::Power(factors[index], div.remainder));
index++;
}
if (r2.isInfinity() || r1.isInfinity()) {
// we overflow Integer at one point: we abort