(***********************************************************************
Mathematica-Compatible Notebook
This notebook can be used on any computer system with Mathematica 4.0,
MathReader 4.0, or any compatible application. The data for the notebook
starts with the line containing stars above.
To get the notebook into a Mathematica-compatible application, do one of
the following:
* Save the data starting with the line of stars above into a file
with a name ending in .nb, then open the file inside the application;
* Copy the data starting with the line of stars above to the
clipboard, then use the Paste menu command inside the application.
Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode. Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).
NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing the
word CacheID, otherwise Mathematica-compatible applications may try to
use invalid cache data.
For more information on notebooks and Mathematica-compatible
applications, contact Wolfram Research:
web: http://www.wolfram.com
email: info@wolfram.com
phone: +1-217-398-0700 (U.S.)
Notebook reader applications are available free of charge from
Wolfram Research.
***********************************************************************)
(*CacheID: 232*)
(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[ 33722, 958]*)
(*NotebookOutlinePosition[ 34396, 982]*)
(* CellTagsIndexPosition[ 34352, 978]*)
(*WindowFrame->Normal*)
Notebook[{
Cell[TextData[{
"Here we will implement an algorithm for computing the 'inverse' of Euler \
totient function ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](n)\)]],
". ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](n)\)]],
" is defined as the number of positive integers less than ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" and relatively prime to ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
". The inverse function for a given ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
" will return the list of all the numbers ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" such that ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](n) = m\)]],
"."
}], "Text",
CellMargins->{{12, 10}, {Inherited, Inherited}}],
Cell[TextData[{
"How can ",
StyleBox["all",
FontSlant->"Italic"],
" the solutions to the equation ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](n) = m\)]],
" be found? Naturally, let's start by looking at the decomposition of ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
" and ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" into primes. Let\n\n \
",
Cell[BoxData[
\(TraditionalForm\`n = \[Product]\+\(k = 1\)\%r p\_k\%\(a\_k\), \
m = \[Product]\+\(k = 1\)\%s q\_k\%\(b\_k\), \[Phi](n) = m\)]],
"\n\nEuler function can be expressed explicitly via the prime factors of ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
":\n\n \
",
Cell[BoxData[
\(TraditionalForm\`\[Phi](
n) = \(\[Product]\+\(k = 1\)\%r\( p\_k\%\(a\_k - 1\)\)(
p\_k - 1) = \(m = \[Product]\+\(k = 1\)\%s \
q\_k\%\(b\_k\)\)\)\)]],
" ",
StyleBox["( 1 )",
FontWeight->"Bold"]
}], "Text"],
Cell[TextData[{
"It immediately follows from (1) that ",
Cell[BoxData[
\(TraditionalForm\`p\_k\%\(a\_k - 1\)\)]],
" and ",
Cell[BoxData[
\(TraditionalForm\`\((p\_k - 1)\)\)]],
" cannot contain any prime factors that are not in the set ",
Cell[BoxData[
\(TraditionalForm\`Q = {q\_k}\_\(k = 1\)\%s\)]],
". Hence we can rewrite their decompositions thus:\n\n \
",
Cell[BoxData[
\(TraditionalForm\`\[ForAll] k = \(\(1, r\)\&_\ \ a\_k = \(1 \[Or] \
\[Exists] q\_i = p\_k\)\)\)]],
"\nand\n \
",
Cell[BoxData[
\(TraditionalForm\`\[ForAll] k = \(\(1, r\)\&_\ \ p\_k = \(2 \[Or]
p\_k = \[Product]\+\(i = 1\)\%\(t\_k\)q\_\(k, i\)\%\(b\_\(k, \
i\)\) + 1\)\)\)]],
"\n\nwhere each ",
Cell[BoxData[
\(TraditionalForm\`q\_\(k, i\)\)]],
" belongs to ",
Cell[BoxData[
\(TraditionalForm\`Q\)]],
". Allowing the powers ",
Cell[BoxData[
\(TraditionalForm\`b\_\(k, i\)\)]],
" to be zero, we can rewrite the last formula including all the ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" in the product:\n\n \
",
Cell[BoxData[
\(TraditionalForm\`\[ForAll] k = \(\(1, r\)\&_\ \ p\_k = \
\[Product]\+\(i = 1\)\%s q\_i\%\(b\_\(k, i\)\) + 1\)\)]]
}], "Text"],
Cell[TextData[{
"Another fact that follows from (1) is that ",
Cell[BoxData[
\(TraditionalForm\`\[Product]\+1\%r\((p\_k - 1)\)\)]],
" divides ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
". Therefore, for each ",
Cell[BoxData[
\(TraditionalForm\`i, k\)]],
" we have ",
Cell[BoxData[
\(TraditionalForm\`b\_\(k, i\) \[LessEqual] b\_i\)]],
". This means that we can find all the ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" that can possibly occur in ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" by simply trying all possible products of ",
Cell[BoxData[
\(TraditionalForm\`q\_i\)]],
" raised to any power not greater than ",
Cell[BoxData[
\(TraditionalForm\`b\_i\)]],
". And by transforming (1) we get\n\n \
",
Cell[BoxData[
\(TraditionalForm\`\[Phi](
n) = \(\[Product]\+1\%r\(\( p\_k\%\(a\_k\)\)(p\_k - 1)\)\/p\_k = \
\(n \(\[Product]\+1\%r\((1 -
1\/p\_k)\)\) = \(\(m\)\(\ \)\(\[Implies]\)\)\)\)\)]],
"\n \n \
\
",
Cell[BoxData[
\(TraditionalForm\`\(\(\[Implies]\)\(\ \)\(n\)\) =
m \(\[Product]\+1\%r\((1 - 1\/p\_k)\)\^\(-1\)\)\)]],
" \
",
StyleBox["( 2 )",
FontWeight->"Bold"],
"\n\nso we can get all the solutions ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" by substituting all possible combinations of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" into (2). Not all ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" that we get this way will satisfy the equation ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](n) = m\)]],
", but right now it is important to us that no solutions will be lost."
}], "Text"],
Cell[TextData[{
"Another simple way to find ",
Cell[BoxData[
\(TraditionalForm\`\(\[Phi]\^\(-1\)\)(m)\)]],
" would be to estimate the range containing all the solutions. First, ",
Cell[BoxData[
\(TraditionalForm\`\((1 - 1\/p\_k)\)\^\(-1\) > 1\)]],
", and second, if ",
Cell[BoxData[
\(TraditionalForm\`p\_i < p\_j\)]],
", then ",
Cell[BoxData[
\(TraditionalForm\`\((1 - 1\/p\_i)\)\^\(-1\) > \((1 - \
1\/p\_j)\)\^\(-1\)\)]],
". If the total number of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" is ",
Cell[BoxData[
\(TraditionalForm\`R\)]],
" and they are sorted from smallest to largest: ",
Cell[BoxData[
\(TraditionalForm\`\(p\_1 < p\_2 < ... \) < p\_R\)]],
", then \n\n \
",
Cell[BoxData[
\(TraditionalForm\`\(m(1 - 1\/p\_R)\)\^\(-1\) \[LessEqual]
n \[LessEqual]
m \(\[Product]\+\(k = 1\)\%R\((1 - 1\/p\_k)\)\^\(-1\)\)\)]],
"\n\nNow we need to estimate the values of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
". The largest is not greater than ",
Cell[BoxData[
\(TraditionalForm\`\((m + 1)\)\)]],
" because each ",
Cell[BoxData[
\(TraditionalForm\`\((p\_k - 1)\)\)]],
" divides ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
" (see above). And if we denote the ",
Cell[BoxData[
\(TraditionalForm\`k\)]],
"th prime in the row of natural numbers as ",
Cell[BoxData[
\(TraditionalForm\`\[Pi]\_k\)]],
", then ",
Cell[BoxData[
\(TraditionalForm\`p\_k \[GreaterEqual] \[Pi]\_k\)]],
". Additionally, for all ",
Cell[BoxData[
\(TraditionalForm\`k\)]],
" this inequality holds: ",
Cell[BoxData[
\(TraditionalForm\`\[Pi]\_k > k\ \(ln(k)\)\)]],
". For ",
Cell[BoxData[
\(TraditionalForm\`k\)]],
" equal to one the right side is zero, so we use this inequality for ",
Cell[BoxData[
\(TraditionalForm\`k \[GreaterEqual] 2\)]],
" (and explicitly substitute ",
Cell[BoxData[
\(TraditionalForm\`\[Pi]\_1 = 2\)]],
"), and get\n\n ",
Cell[BoxData[
\(TraditionalForm\`m \(\[Product]\+\(k = 1\)\%R\((1 - \
1\/p\_k)\)\^\(-1\)\) \[LessEqual]
m \(\[Product]\+\(k = 1\)\%R\((1 - 1\/\[Pi]\_k)\)\^\(-1\)\) \
\[LessEqual]
2 m \(\[Product]\+\(k = 2\)\%R\((1 - 1\/\(k\ \
\(ln(k)\)\))\)\^\(-1\)\)\)]],
" \n\nThe constant ",
Cell[BoxData[
\(TraditionalForm\`R\)]],
" is not greater than ",
Cell[BoxData[
\(TraditionalForm\`\[Chi](m)\)]],
", the number of divisors of ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
" (again, because ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
" is divisible by each ",
Cell[BoxData[
\(TraditionalForm\`\((p\_k - 1)\)\)]],
"), and finally we get\n\n \
",
Cell[BoxData[
\(TraditionalForm\`m + 1 \[LessEqual] n \[LessEqual]
2 m \(\[Product]\+\(k = 2\)\%\(\[Chi](m)\)\((1 - 1\/\(k\ \
\(ln(k)\)\))\)\^\(-1\)\)\)]]
}], "Text"],
Cell[BoxData[
\(invphidemo1[m_] :=
Module[\[IndentingNewLine]{Lp, Lq, Lb, Lpow, Lfact,
n}, \[IndentingNewLine]{Lq, Lb} =
Transpose[FactorInteger[m]]; \[IndentingNewLine]Lpow =
MapThread[\[IndentingNewLine]Table[#1^i, {i,
0, #2}] &, \[IndentingNewLine]{Lq,
Lb}]; \[IndentingNewLine]Lp =
Outer[Times, Sequence @@ Lpow]; \[IndentingNewLine]Lp =
Select[Flatten[Lp] + 1, PrimeQ]; \[IndentingNewLine]Lfact =
Table[\[IndentingNewLine]\((Lp[\([i]\)]/\((Lp[\([i]\)] - 1)\))\)^
j, \[IndentingNewLine]{i, Length[Lp]}, {j, 0,
1}]; \[IndentingNewLine]n =
m\ Outer[Times, Sequence @@ Lfact]; \[IndentingNewLine]Select[
Flatten[n], \[IndentingNewLine]EulerPhi[#] \[Equal]
m &]\[IndentingNewLine]]\)], "Input"],
Cell[BoxData[
\(invphidemo2[m_] :=
Select[\[IndentingNewLine]Range[m + 1,
2 m \(\[Product]\+\(k = 2\)\%\(DivisorSigma[0, m]\)\((1 - 1\/\(k\ \
Log[k]\))\)\^\(-1\)\)], \[IndentingNewLine]EulerPhi[#] \[Equal]
m &]\)], "Input"],
Cell[CellGroupData[{
Cell[BoxData[{
\(invphidemo1[12]\), "\[IndentingNewLine]",
\(invphidemo2[12]\)}], "Input"],
Cell[BoxData[
\({13, 21, 26, 28, 36, 42}\)], "Output"],
Cell[BoxData[
\({13, 21, 26, 28, 36, 42}\)], "Output"]
}, Open ]],
Cell[TextData[{
"Of course, these are bad solutions. They are very slow, and the first one \
requires prodigious amounts of memory. The reason is that if the powers ",
Cell[BoxData[
\(TraditionalForm\`b\_k\)]],
" are large, the set of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" can be many times larger than the set ",
Cell[BoxData[
\(TraditionalForm\`Q\)]],
", and creating the table of all the possible products of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" is nearly impossible. The function ",
Cell[BoxData[
\(TraditionalForm\`invphi\)]],
" still will involve the search for all combinations of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" that yield solutions but our goal is to make this search more efficient."
}], "Text"],
Cell[TextData[{
"Let's denote the set of all the ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" that can possibly occur in a solution as ",
Cell[BoxData[
\(TraditionalForm\`P = {p\_k}\_\(k = 1\)\%R\)]],
". We want to create list ",
Cell[BoxData[
\(TraditionalForm\`L\_r = {p\_\(i\_1\)\ , p\_\(i\_2\)\ , ... ,
p\_\(i\_r\)}\)]],
", with elements from ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
", corresponding to one of the solutions ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
". This means that we'll get a solution by applying formula (2) to ",
Cell[BoxData[
\(TraditionalForm\`L\_r\)]],
". Note that ",
Cell[BoxData[
\(TraditionalForm\`L\_r\)]],
" denotes the list with ",
Cell[BoxData[
\(TraditionalForm\`r\)]],
" elements and not the ",
Cell[BoxData[
\(TraditionalForm\`r\)]],
"th element. This notation will be handy for us because ",
Cell[BoxData[
\(TraditionalForm\`r\)]],
" is different for each ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" and is not known in advance. Let on some step in the algorithm we already \
have ",
Cell[BoxData[
\(TraditionalForm\`t\)]],
" elements ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
". Then the problem is reduced to deciding adding which ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" to ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" can yield a valid solution, and which ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" needn't be tried, being 'incompatible' with those ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" that are in the list already."
}], "Text"],
Cell[TextData[{
"First, as was already shown, the product ",
Cell[BoxData[
\(TraditionalForm\`\[Product]\+\(k = 1\)\%r\((p\_\(i\_k\) - 1)\)\)]],
" divides ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
", so their quotient\n\n \
",
Cell[BoxData[
\(TraditionalForm\`d\_r =
m \(\[Product]\+\(k = 1\)\%r\((p\_\(i\_k\) - 1)\)\^\(-1\)\)\)]],
"\n\nis integer. Therefore, before adding ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" to ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" we should check that the resulting quotient ",
Cell[BoxData[
\(TraditionalForm\`d\_\(t + 1\)\)]],
" is integer, and if it's not, ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" is not added."
}], "Text"],
Cell[TextData[{
"Second, ",
Cell[BoxData[
\(TraditionalForm\`d\_r\)]],
" contains only prime factors from ",
Cell[BoxData[
\(TraditionalForm\`Q\)]],
" (because ",
Cell[BoxData[
\(TraditionalForm\`d\_r\)]],
" is ",
Cell[BoxData[
\(TraditionalForm\`m\)]],
" divided by something). But by formula (1), \n\n \
",
Cell[BoxData[
\(TraditionalForm\`d\_r = \[Product]\+\(k = 1\)\%r \
p\_\(i\_k\)\%\(a\_\(i\_k\) - 1\)\)]],
"\n\nIt follows that only those primes can occur in ",
Cell[BoxData[
\(TraditionalForm\`d\_r\)]],
" which belong to ",
Cell[BoxData[
\(TraditionalForm\`P \[Intersection] Q\)]],
" and are also members of ",
Cell[BoxData[
\(TraditionalForm\`L\_r\)]],
". So on the step where we have ",
Cell[BoxData[
\(TraditionalForm\`t\)]],
" elements we should check if the 'current' quotient ",
Cell[BoxData[
\(TraditionalForm\`d\_t\)]],
" contains a factor ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" which is not a member of ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
". If such ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" exists, then it cannot appear in the 'final' quotient ",
Cell[BoxData[
\(TraditionalForm\`d\_r\)]],
" unless it belongs to ",
Cell[BoxData[
\(TraditionalForm\`L\_r\)]],
". Therefore, we can add ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" to ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" (possible only if ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" belongs to ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
") or eliminate ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" from the quotient on the next steps, if we want to arrive to a solution. \
The elimination can be done by selecting ",
Cell[BoxData[
\(TraditionalForm\`p\_\(i\_\(t + 1\)\)\)]],
" such that ",
Cell[BoxData[
\(TraditionalForm\`\((p\_\(i\_\(t + 1\)\) - 1)\)\)]],
" is divisible by ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
", so in the next quotient ",
Cell[BoxData[
\(TraditionalForm\`d\_\(t + 1\) =
d\_t/\((p\_\(i\_\(t + 1\)\) - 1)\)\)]],
" the power to which ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" appears in it will be lowered (though ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" may still be present in ",
Cell[BoxData[
\(TraditionalForm\`d\_\(t + 1\)\)]],
"). How this all benefits the algorithm is that instead of trying to add \
all the possible ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" that are not in the list yet, we can find out that there is a much \
smaller subset of ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
" such that at least one of its elements should necessarily be appended to \
",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
", in order to eliminate a ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
". Then we may try building new lists from ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" by adding only ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" from that subset first (adding them separately, thus creating several new \
'candidates' from ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
"). This way we can significantly reduce the amount of search without \
losing any of the solutions."
}], "Text"],
Cell[TextData[{
"Another important consequence is that if want to know whether ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" corresponds to a valid solution or not, it can be checked without \
computing ",
Cell[BoxData[
\(TraditionalForm\`EulerPhi\)]],
". The criterion of its validity is exactly the condition discussed above, \
namely that ",
Cell[BoxData[
\(TraditionalForm\`d\_t\)]],
" can contain only those ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" that appear in ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
"."
}], "Text"],
Cell[TextData[{
"The program ",
Cell[BoxData[
\(TraditionalForm\`invphi\)]],
" first finds ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
". ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
" and ",
Cell[BoxData[
\(TraditionalForm\`Q\)]],
" are rearranged so that all their common elements, the number of which is \
denoted ",
Cell[BoxData[
\(TraditionalForm\`r\_0\)]],
", come in the beginning and have the same positions in both lists. The \
main structure ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
" is a list of solution candidates. Each element of ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
" has in turn two elements. The first one defines ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
". It contains an indicator for each ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" from ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
" showing if that ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" is already in ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
", can, or can not be added to it (1, 0, -1 respectively). The second \
element is ",
Cell[BoxData[
\(TraditionalForm\`d\_t\)]],
". On each step we check if there are any factors ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" in the quotient ",
Cell[BoxData[
\(TraditionalForm\`d\_t\)]],
" that should be eliminated. If there are, then we see how many ",
Cell[BoxData[
\(TraditionalForm\`p\_i\)]],
" can reduce the power to which ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" occurs in ",
Cell[BoxData[
\(TraditionalForm\`d\_t\)]],
", ie how many ",
Cell[BoxData[
\(TraditionalForm\`p\_i\)]],
" have the property that ",
Cell[BoxData[
\(TraditionalForm\`\((p\_i - 1)\)\)]],
" is divisible by ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
". (As discussed above, another possibility is ",
Cell[BoxData[
\(TraditionalForm\`p\_i = q\_k\)]],
"). We always select the ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" for which the number of those ",
Cell[BoxData[
\(TraditionalForm\`p\_i\)]],
" is the smallest, thus trying to minimize the number of new candidates \
that will be added. If no ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" to be eliminated are found, then all members of ",
Cell[BoxData[
\(TraditionalForm\`P\)]],
" that are allowed for the current ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" are added to it. For each new element of ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
" thus created we find the new quotient ",
Cell[BoxData[
\(TraditionalForm\`d\_\(t + 1\)\)]],
", and mark as 'incompatible' those ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" that would make ",
Cell[BoxData[
\(TraditionalForm\`d\_\(t + 2\)\)]],
" fractional. Also, if we are adding several ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" to ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
", say ",
Cell[BoxData[
\(TraditionalForm\`p\_1\)]],
" and ",
Cell[BoxData[
\(TraditionalForm\`p\_2\)]],
", creating two new lists, ",
Cell[BoxData[
\(TraditionalForm\`L\_\(t + 1\)\%\((1)\) = Append[L\_t, p\_1]\)]],
" and ",
Cell[BoxData[
\(TraditionalForm\`L\_\(t + 1\)\%\((2)\) = Append[L\_t, p\_2]\)]],
", then in ",
Cell[BoxData[
\(TraditionalForm\`L\_\(t + 1\)\%\((2)\)\)]],
" ",
Cell[BoxData[
\(TraditionalForm\`p\_1\)]],
" is marked as 'incompatible' too in order to avoid possible repetitions of \
the same sets of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
" later. When nothing further can be added to ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
", the final list of solutions is generated by selecting those elements \
from ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
" that yield correct values of ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
". This is done in a similar way; we consider the corresponding quotient, \
and if no ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" need to be eliminated, then we have a solution."
}], "Text"],
Cell[TextData[{
"There are some other ways to improve the performance. It can be easily \
noticed that if we add ",
Cell[BoxData[
\(TraditionalForm\`p\_k = 2\)]],
" to the list, then the quotient doesn't change on the next step. We can \
speed up the search by never considering ",
Cell[BoxData[
\(TraditionalForm\`p\_k = 2\)]],
", and only taking into consideration the possible presence of this factor \
later, when generating solutions from the lists of ",
Cell[BoxData[
\(TraditionalForm\`p\_k\)]],
". Also, for selecting the ",
Cell[BoxData[
\(TraditionalForm\`p\_i\)]],
" that can eliminate a given ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
", the structure ",
Cell[BoxData[
\(TraditionalForm\`Mdiv\)]],
" is created in advance, in which ",
Cell[BoxData[
\(TraditionalForm\`Mdiv\[LeftDoubleBracket]k\[RightDoubleBracket]\)]],
" is the list of indices ",
Cell[BoxData[
\(TraditionalForm\`i\)]],
" such that ",
Cell[BoxData[
\(TraditionalForm\`q\_k\)]],
" divides ",
Cell[BoxData[
\(TraditionalForm\`\((p\_i - 1)\)\)]],
". Finally, after we have created new lists from ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
", we'll never return to ",
Cell[BoxData[
\(TraditionalForm\`L\_t\)]],
" again. When the structure ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
" grows large, we can convert those elements to the final result ",
Cell[BoxData[
\(TraditionalForm\`n\)]],
" and remove them from ",
Cell[BoxData[
\(TraditionalForm\`wrk\)]],
". This not only reduces the memory usage, but makes large computations \
significantly faster, because operations will be performed with smaller \
lists."
}], "Text"],
Cell[TextData[{
"In ",
StyleBox["Mathematica",
FontSlant->"Italic"],
", ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](0) = 0, \[Phi](1) = 1\)]],
". ",
Cell[BoxData[
\(TraditionalForm\`invphi\)]],
" complies with this conventions. Also note that only the list of \
nonnegative solutions is returned, while in ",
StyleBox["Mathematica",
FontSlant->"Italic"],
" ",
Cell[BoxData[
\(TraditionalForm\`\[Phi](\(-n\)) = \[Phi](n)\)]],
"."
}], "Text"],
Cell[BoxData[
\(invphi[m_Integer] :=
Module[\[IndentingNewLine]{main, init, genp, bestcand, gencand,
addcand, genans, \[IndentingNewLine]Lp, Lq, r, s, r0,
Mdiv}, \[IndentingNewLine]\[IndentingNewLine]main[] :=
Module[\[IndentingNewLine]{ans = {}, wrk,
threshold = 100, \[IndentingNewLine]Lstate, Ladd, quo, indx,
i}, \[IndentingNewLine]wrk = {{Table[0, {r}],
m}}; \[IndentingNewLine]wrk[\([1, 1,
1]\)] = \(-1\); \[IndentingNewLine]For[i = 1,
i \[LessEqual] Length[wrk], \(i++\), \[IndentingNewLine]If[
i \[Equal]
threshold + 1, \[IndentingNewLine]ans = {ans,
genans[Take[wrk,
threshold]\ ]\ }; \[IndentingNewLine]wrk =
Drop[wrk, threshold]; \[IndentingNewLine]i =
1\[IndentingNewLine]]; \[IndentingNewLine]{Lstate, quo} =
wrk[\([i]\)]; \[IndentingNewLine]indx =
bestcand[Lstate, quo]; \[IndentingNewLine]Ladd =
gencand[Lstate, indx]; \[IndentingNewLine]wrk =
Join[wrk,
addcand[Lstate, quo,
Ladd]\ ]\[IndentingNewLine]]; \
\[IndentingNewLine]Flatten[{ans,
genans[wrk]}]\[IndentingNewLine]]; \[IndentingNewLine]\
\[IndentingNewLine]init[] :=
Module[\[IndentingNewLine]{Lb, Lpq}, \[IndentingNewLine]{Lq, Lb} =
Transpose[FactorInteger[m]]; \[IndentingNewLine]genp[
Lb]; \[IndentingNewLine]Lpq =
Intersection[Lp, Lq]; \[IndentingNewLine]{Lp,
Lq} = \(Join[
Lpq, \[IndentingNewLine]Complement[#,
Lpq]\ ] &\) /@ \[IndentingNewLine]{Lp,
Lq}; \[IndentingNewLine]{r, s, r0} =
Length /@ {Lp, Lq,
Lpq}; \[IndentingNewLine]Mdiv = \
\(Cases[\[IndentingNewLine]Range[r], \[IndentingNewLine]x_ /;
Mod[Lp[\([x]\)] - 1, #] \[Equal]
0\[IndentingNewLine]] &\) /@
Lq;\[IndentingNewLine]]; \[IndentingNewLine]\
\[IndentingNewLine]genp[Lb_] :=
Module[\[IndentingNewLine]{Lpow, tmp}, \[IndentingNewLine]Lpow =
MapThread[\[IndentingNewLine]Table[#1^i, {i,
0, #2}] &, \[IndentingNewLine]{Lq,
Lb}]; \[IndentingNewLine]Lp = {}; \[IndentingNewLine]Outer[\
\[IndentingNewLine]If[
PrimeQ[tmp = Times[##] + 1], \[IndentingNewLine]Lp = {Lp,
tmp}\[IndentingNewLine]] &, \[IndentingNewLine]Sequence \
@@ Lpow]; \[IndentingNewLine]Lp =
Flatten[Lp];\[IndentingNewLine]]; \[IndentingNewLine]\
\[IndentingNewLine]bestcand[Lstate_, quo_] :=
Module[\[IndentingNewLine]{len = Infinity, indx = 0, cur,
i}, \[IndentingNewLine]For[i = 1,
i \[LessEqual]
s, \(i++\), \[IndentingNewLine]If[\((\((i \[LessEqual]
r0 && \[IndentingNewLine]Lstate[\([i]\)] \
\[NotEqual] 1)\) || \[IndentingNewLine]i > r0)\) && \[IndentingNewLine]Mod[
quo, Lq[\([i]\)]\ ] \[Equal]
0, \[IndentingNewLine]cur =
Length[Mdiv[\([i]\)]\ ]; \[IndentingNewLine]If[
cur < len, \[IndentingNewLine]len =
cur; \[IndentingNewLine]indx =
i\[IndentingNewLine]]\[IndentingNewLine]]\
\[IndentingNewLine]]; \[IndentingNewLine]indx\[IndentingNewLine]]; \
\[IndentingNewLine]\[IndentingNewLine]gencand[Lstate_, indx_] :=
Module[\[IndentingNewLine]{Ladd}, \[IndentingNewLine]Ladd = \
\[IndentingNewLine]If[
indx \[NotEqual] 0, \[IndentingNewLine]If[
indx \[LessEqual] r0, \[IndentingNewLine]Prepend[
Mdiv[\([indx]\)],
indx], \[IndentingNewLine]Mdiv[\([indx]\)]\
\[IndentingNewLine]], \[IndentingNewLine]Range[
r]\[IndentingNewLine]]; \[IndentingNewLine]Select[Ladd,
Lstate[\([#]\)] \[Equal]
0 &]\[IndentingNewLine]]; \[IndentingNewLine]\
\[IndentingNewLine]addcand[Lstate_, quo_, Ladd_] :=
Module[\[IndentingNewLine]{ans = {}, Lstate2, quo2, len,
i}, \[IndentingNewLine]len =
Length[Ladd]; \[IndentingNewLine]For[i = 1,
i \[LessEqual] len, \(i++\), \[IndentingNewLine]Lstate2 =
ReplacePart[Lstate, 1,
Ladd[\([i]\)]\ ]; \[IndentingNewLine]quo2 =
quo/\((Lp[\([Ladd[\([i]\)]\ ]\)] -
1)\); \[IndentingNewLine]\(\((Lstate2[\([Ladd[\([#]\)]\ \
]\)] = \(-1\))\) &\) /@ \[IndentingNewLine]Range[
i - 1]; \[IndentingNewLine]\(If[
Lstate2[\([#]\)] \[Equal]
0 && \[IndentingNewLine]Mod[quo2,
Lp[\([#]\)] - 1] \[NotEqual]
0, \[IndentingNewLine]Lstate2[\([#]\)] = \(-1\)] &\) \
/@ \[IndentingNewLine]Range[r]; \[IndentingNewLine]AppendTo[
ans, {Lstate2,
quo2}]\[IndentingNewLine]]; \[IndentingNewLine]ans\
\[IndentingNewLine]]; \[IndentingNewLine]\[IndentingNewLine]genans[L_] :=
Module[\[IndentingNewLine]{ans = {}, Lstate, quo, res, add2, i,
j}, \[IndentingNewLine]For[i = 1,
i \[LessEqual]
Length[L], \(i++\), \[IndentingNewLine]{Lstate, quo} =
L[\([i]\)]; \[IndentingNewLine]For[add2 = 0,
add2 \[LessEqual] 1, \(add2++\), \[IndentingNewLine]If[
add2 \[Equal] 1, \[IndentingNewLine]Lstate[\([1]\)] =
1\[IndentingNewLine]]; \[IndentingNewLine]For[j = 1,
j \[LessEqual]
s, \(j++\), \[IndentingNewLine]If[\((\((j \[LessEqual]
r0 && \[IndentingNewLine]Lstate[\([j]\)] \
\[NotEqual] 1)\) || \[IndentingNewLine]j > r0)\) && \[IndentingNewLine]Mod[
quo, Lq[\([j]\)]\ ] \[Equal]
0, \[IndentingNewLine]Break[]\[IndentingNewLine]]\
\[IndentingNewLine]]; \[IndentingNewLine]If[
j \[NotEqual]
s + 1, \[IndentingNewLine]Continue[]\[IndentingNewLine]]; \
\[IndentingNewLine]res =
Cases[\[IndentingNewLine]Transpose[{Lp,
Lstate}], \[IndentingNewLine]{x_, 1} \[Rule]
x]; \[IndentingNewLine]res =
m\ Times @@ res/
Times @@ \((res - 1)\); \[IndentingNewLine]ans = {ans,
res}\[IndentingNewLine]]\[IndentingNewLine]]; \
\[IndentingNewLine]ans\[IndentingNewLine]]; \[IndentingNewLine]\
\[IndentingNewLine]Switch[m, \[IndentingNewLine]0,
Return[{0}], \[IndentingNewLine]1,
Return[{1,
2}], \[IndentingNewLine]_?\((OddQ[#] || Negative[#] &)\),
Return[{}]\[IndentingNewLine]]; \[IndentingNewLine]init[]; \
\[IndentingNewLine]main[]\[IndentingNewLine]]\)], "Input"],
Cell["Here finding all 847 solutions takes only a few seconds.", "Text"],
Cell[CellGroupData[{
Cell[BoxData[
\(Timing[Length[invphi[\(2\^5\) \(3\^5\) 5\^5]\ ]\ ]\)], "Input"],
Cell[BoxData[
\({6.369999999999891`\ Second, 847}\)], "Output"]
}, Open ]],
Cell["\<\
In this case, factoring each of the resulting integers would take about 150 \
seconds, so checking the answer takes many times longer than finding it.\
\>", "Text"],
Cell[CellGroupData[{
Cell[BoxData[
\(Timing[
invphi[\(2\^3\) \(31\^3\) \(313\^3\) \(619\^3\) \(1511\^3\)
1733\^3]\ ]\)], "Input"],
Cell[BoxData[
\({6.480000000000018`\ Second, \
{62244230888910733677707758056446080408564,
46683173166683050258280818542334560306423,
93366346333366100516561637084669120612846}}\)], "Output"]
}, Open ]],
Cell["Copyright \[Copyright] 1999, Maxim Rytin", "Text"]
},
FrontEndVersion->"4.0 for Microsoft Windows",
ScreenRectangle->{{0, 1024}, {0, 723}},
WindowToolbars->{},
WindowSize->{1016, 696},
WindowMargins->{{0, Automatic}, {Automatic, 0}}
]
(***********************************************************************
Cached data follows. If you edit this Notebook file directly, not using
Mathematica, you must remove the line containing CacheID at the top of
the file. The cache data will then be recreated when you save this file
from within Mathematica.
***********************************************************************)
(*CellTagsOutline
CellTagsIndex->{}
*)
(*CellTagsIndex
CellTagsIndex->{}
*)
(*NotebookFileOutline
Notebook[{
Cell[1717, 49, 753, 25, 60, "Text"],
Cell[2473, 76, 1164, 31, 185, "Text"],
Cell[3640, 109, 1492, 38, 276, "Text"],
Cell[5135, 149, 2009, 55, 248, "Text"],
Cell[7147, 206, 3151, 94, 398, "Text"],
Cell[10301, 302, 870, 15, 354, "Input"],
Cell[11174, 319, 259, 5, 127, "Input"],
Cell[CellGroupData[{
Cell[11458, 328, 99, 2, 55, "Input"],
Cell[11560, 332, 58, 1, 32, "Output"],
Cell[11621, 335, 58, 1, 32, "Output"]
}, Open ]],
Cell[11694, 339, 804, 21, 83, "Text"],
Cell[12501, 362, 1696, 55, 130, "Text"],
Cell[14200, 419, 851, 25, 134, "Text"],
Cell[15054, 446, 3510, 113, 307, "Text"],
Cell[18567, 561, 591, 19, 60, "Text"],
Cell[19161, 582, 4228, 141, 269, "Text"],
Cell[23392, 725, 1762, 52, 152, "Text"],
Cell[25157, 779, 487, 18, 60, "Text"],
Cell[25647, 799, 7188, 125, 3137, "Input"],
Cell[32838, 926, 72, 0, 37, "Text"],
Cell[CellGroupData[{
Cell[32935, 930, 83, 1, 37, "Input"],
Cell[33021, 933, 67, 1, 32, "Output"]
}, Open ]],
Cell[33103, 937, 174, 3, 37, "Text"],
Cell[CellGroupData[{
Cell[33302, 944, 129, 3, 37, "Input"],
Cell[33434, 949, 213, 4, 55, "Output"]
}, Open ]],
Cell[33662, 956, 56, 0, 37, "Text"]
}
]
*)
(***********************************************************************
End of Mathematica Notebook file.
***********************************************************************)