계산 phi(k)1
https://stackoverflow.com/questions/1024640
italiano
english
français
española
中国
日本の
العربية
Deutsch
한국어
Português
Russian문제
주어진 큰 N,내가 필요한 반복을 통해 모든 phi(k)는 1 < k < N:
- 시간 복잡성야
O(N logN) - 메모리 복잡해야 하위
O(N)(이후의 값 N 될 것입니다 주위에 1012)
이것이 가능합니까?그렇다면,어떻게 합니까?
해결책
이것은 아래 코드 예제에 구현 된 바와 같이, 메모리 복잡성 O (SQRT (N)) 및 CPU 복잡성 O (N * LOG (LOG (N))를 최적화 된 윈도우의 Eratosthenes를 사용하여 수행 할 수 있습니다.
언어가 지정되지 않았고 Python을 알지 못하면 vb.net에서 구현했지만 필요한 경우 C#으로 변환 할 수 있습니다.
Imports System.Math
Public Class TotientSerialCalculator
'Implements an extremely efficient Serial Totient(phi) calculator '
' This implements an optimized windowed Sieve of Eratosthenes. The'
' window size is set at Sqrt(N) both to optimize collecting and '
' applying all of the Primes below Sqrt(N), and to minimize '
' window-turning overhead. '
' '
' CPU complexity is O( N * Log(Log(N)) ), which is virtually linear.'
' '
' MEM Complexity is O( Sqrt(N) ). '
' '
' This is probalby the ideal combination, as any attempt to further '
'reduce memory will almost certainly result in disproportionate increases'
'in CPU complexity, and vice-versa. '
Structure NumberFactors
Dim UnFactored As Long 'the part of the number that still needs to be factored'
Dim Phi As Long 'the totient value progressively calculated'
' (equals total numbers less than N that are CoPrime to N)'
'MEM = 8 bytes each'
End Structure
Private ReportInterval As Long
Private PrevLast As Long 'the last value in the previous window'
Private FirstValue As Long 'the first value in this windows range'
Private WindowSize As Long
Private LastValue As Long 'the last value in this windows range'
Private NextFirst As Long 'the first value in the next window'
'Array that stores all of the NumberFactors in the current window.'
' this is the primary memory consumption for the class and it'
' is 16 * Sqrt(N) Bytes, which is O(Sqrt(N)).'
Public Numbers() As NumberFactors
' For N=10^12 (1 trilion), this will be 16MB, which should be bearable anywhere.'
'(note that the Primes() array is a secondary memory consumer'
' at O(pi(Sqrt(N)), which will be within 10x of O(Sqrt(N)))'
Public Event EmitTotientPair(ByVal k As Long, ByVal Phi As Long)
'===== The Routine To Call: ========================'
Public Sub EmitTotientPairsToN(ByVal N As Long)
'Routine to Emit Totient pairs {k, Phi(k)} for k = 1 to N'
' 2009-07-14, RBarryYoung, Created.'
Dim i As Long
Dim k As Long 'the current number being factored'
Dim p As Long 'the current prime factor'
'Establish the Window frame:'
' note: WindowSize is the critical value that controls both memory'
' usage and CPU consumption and must be SQRT(N) for it to work optimally.'
WindowSize = Ceiling(Sqrt(CDbl(N)))
ReDim Numbers(0 To WindowSize - 1)
'Initialize the first window:'
MapWindow(1)
Dim IsFirstWindow As Boolean = True
'adjust this to control how often results are show'
ReportInterval = N / 100
'Allocate the primes array to hold the primes list:'
' Only primes <= SQRT(N) are needed for factoring'
' PiMax(X) is a Max estimate of the number of primes <= X'
Dim Primes() As Long, PrimeIndex As Long, NextPrime As Long
'init the primes list and its pointers'
ReDim Primes(0 To PiMax(WindowSize) - 1)
Primes(0) = 2 '"prime" the primes list with the first prime'
NextPrime = 1
'Map (and Remap) the window with Sqrt(N) numbers, Sqrt(N) times to'
' sequentially map all of the numbers <= N.'
Do
'Sieve the primes across the current window'
PrimeIndex = 0
'note: cant use enumerator for the loop below because NextPrime'
' changes during the first window as new primes <= SQRT(N) are accumulated'
Do While PrimeIndex < NextPrime
'get the next prime in the list'
p = Primes(PrimeIndex)
'find the first multiple of (p) in the current window range'
k = PrevLast + p - (PrevLast Mod p)
Do
With Numbers(k - FirstValue)
.UnFactored = .UnFactored \ p 'always works the first time'
.Phi = .Phi * (p - 1) 'Phi = PRODUCT( (Pi-1)*Pi^(Ei-1) )'
'The loop test that follows is probably the central CPU overhead'
' I believe that it is O(N*Log(Log(N)), which is virtually O(N)'
' ( for instance at N = 10^12, Log(Log(N)) = 3.3 )'
Do While (.UnFactored Mod p) = 0
.UnFactored = .UnFactored \ p
.Phi = .Phi * p
Loop
End With
'skip ahead to the next multiple of p: '
'(this is what makes it so fast, never have to try prime factors that dont apply)'
k += p
'repeat until we step out of the current window:'
Loop While k < NextFirst
'if this is the first window, then scan ahead for primes'
If IsFirstWindow Then
For i = Primes(NextPrime - 1) + 1 To p ^ 2 - 1 'the range of possible new primes'
'Dont go beyond the first window'
If i >= WindowSize Then Exit For
If Numbers(i - FirstValue).UnFactored = i Then
'this is a prime less than SQRT(N), so add it to the list.'
Primes(NextPrime) = i
NextPrime += 1
End If
Next
End If
PrimeIndex += 1 'move to the next prime'
Loop
'Now Finish & Emit each one'
For k = FirstValue To LastValue
With Numbers(k - FirstValue)
'Primes larger than Sqrt(N) will not be finished: '
If .UnFactored > 1 Then
'Not done factoring, must be an large prime factor remaining: '
.Phi = .Phi * (.UnFactored - 1)
.UnFactored = 1
End If
'Emit the value pair: (k, Phi(k)) '
EmitPhi(k, .Phi)
End With
Next
're-Map to the next window '
IsFirstWindow = False
MapWindow(NextFirst)
Loop While FirstValue <= N
End Sub
Sub EmitPhi(ByVal k As Long, ByVal Phi As Long)
'just a placeholder for now, that raises an event to the display form'
' periodically for reporting purposes. Change this to do the actual'
' emitting.'
If (k Mod ReportInterval) = 0 Then
RaiseEvent EmitTotientPair(k, Phi)
End If
End Sub
Public Sub MapWindow(ByVal FirstVal As Long)
'Efficiently reset the window so that we do not have to re-allocate it.'
'init all of the boundary values'
FirstValue = FirstVal
PrevLast = FirstValue - 1
NextFirst = FirstValue + WindowSize
LastValue = NextFirst - 1
'Initialize the Numbers prime factor arrays'
Dim i As Long
For i = 0 To WindowSize - 1
With Numbers(i)
.UnFactored = i + FirstValue 'initially equal to the number itself'
.Phi = 1 'starts at mulplicative identity(1)'
End With
Next
End Sub
Function PiMax(ByVal x As Long) As Long
'estimate of pi(n) == {primes <= (n)} that is never less'
' than the actual number of primes. (from P. Dusart, 1999)'
Return (x / Log(x)) * (1.0 + 1.2762 / Log(x))
End Function
End Class
O (n * log (log (n))) 에서이 루틴은 평균적으로 O (log (log (n)))에서 각 숫자를 고려하고 있으며 여기에 답장 중 일부입니다. 사실, n = 10^12에서 그것은 2400 더 빨리!
2GHz Intel Core 2 노트북 에서이 루틴을 테스트했으며 초당 3,000,000 개 이상의 PHI () 값을 계산합니다. 이 속도에서 10^12 값을 계산하는 데 약 4 일이 걸립니다. 또한 오류없이 최대 100,000,000의 정확성을 테스트했습니다. 64 비트 정수에 기반을두고 있으므로 최대 2^63 (10^19)의 모든 것이 정확해야합니다 (누구도 너무 느리지 만).
또한 실행/테스트를위한 Visual Studio WinForm (VB.NET)도 있습니다. 원하는 경우 제공 할 수 있습니다.
궁금한 점이 있으면 알려주세요.
주석에 요청 된대로 C# 버전의 코드를 추가했습니다. 그러나 현재 다른 프로젝트의 중간에 있기 때문에 직접 변환 할 시간이 없으므로 온라인 VB에서 C# 변환 사이트 중 하나를 사용했습니다 (http://www.carlosag.net/tools/codetranslator/). 따라서 이것이 자동으로 변환되었으며 아직 테스트하거나 확인할 시간이 없었습니다.
using System.Math;
public class TotientSerialCalculator {
// Implements an extremely efficient Serial Totient(phi) calculator '
// This implements an optimized windowed Sieve of Eratosthenes. The'
// window size is set at Sqrt(N) both to optimize collecting and '
// applying all of the Primes below Sqrt(N), and to minimize '
// window-turning overhead. '
// '
// CPU complexity is O( N * Log(Log(N)) ), which is virtually linear.'
// '
// MEM Complexity is O( Sqrt(N) ). '
// '
// This is probalby the ideal combination, as any attempt to further '
// reduce memory will almost certainly result in disproportionate increases'
// in CPU complexity, and vice-versa. '
struct NumberFactors {
private long UnFactored; // the part of the number that still needs to be factored'
private long Phi;
}
private long ReportInterval;
private long PrevLast; // the last value in the previous window'
private long FirstValue; // the first value in this windows range'
private long WindowSize;
private long LastValue; // the last value in this windows range'
private long NextFirst; // the first value in the next window'
// Array that stores all of the NumberFactors in the current window.'
// this is the primary memory consumption for the class and it'
// is 16 * Sqrt(N) Bytes, which is O(Sqrt(N)).'
public NumberFactors[] Numbers;
// For N=10^12 (1 trilion), this will be 16MB, which should be bearable anywhere.'
// (note that the Primes() array is a secondary memory consumer'
// at O(pi(Sqrt(N)), which will be within 10x of O(Sqrt(N)))'
//NOTE: this part looks like it did not convert correctly
public event EventHandler EmitTotientPair;
private long k;
private long Phi;
// ===== The Routine To Call: ========================'
public void EmitTotientPairsToN(long N) {
// Routine to Emit Totient pairs {k, Phi(k)} for k = 1 to N'
// 2009-07-14, RBarryYoung, Created.'
long i;
long k;
// the current number being factored'
long p;
// the current prime factor'
// Establish the Window frame:'
// note: WindowSize is the critical value that controls both memory'
// usage and CPU consumption and must be SQRT(N) for it to work optimally.'
WindowSize = Ceiling(Sqrt(double.Parse(N)));
object Numbers;
this.MapWindow(1);
bool IsFirstWindow = true;
ReportInterval = (N / 100);
// Allocate the primes array to hold the primes list:'
// Only primes <= SQRT(N) are needed for factoring'
// PiMax(X) is a Max estimate of the number of primes <= X'
long[] Primes;
long PrimeIndex;
long NextPrime;
// init the primes list and its pointers'
object Primes;
-1;
Primes[0] = 2;
// "prime" the primes list with the first prime'
NextPrime = 1;
// Map (and Remap) the window with Sqrt(N) numbers, Sqrt(N) times to'
// sequentially map all of the numbers <= N.'
for (
; (FirstValue <= N);
) {
PrimeIndex = 0;
// note: cant use enumerator for the loop below because NextPrime'
// changes during the first window as new primes <= SQRT(N) are accumulated'
while ((PrimeIndex < NextPrime)) {
// get the next prime in the list'
p = Primes[PrimeIndex];
// find the first multiple of (p) in the current window range'
k = (PrevLast
+ (p
- (PrevLast % p)));
for (
; (k < NextFirst);
) {
// With...
UnFactored;
p;
// always works the first time'
(Phi
* (p - 1));
while (// TODO: Warning!!!! NULL EXPRESSION DETECTED...
) {
(UnFactored % p);
UnFactored;
(Phi * p);
}
// skip ahead to the next multiple of p: '
// (this is what makes it so fast, never have to try prime factors that dont apply)'
k = (k + p);
// repeat until we step out of the current window:'
}
// if this is the first window, then scan ahead for primes'
if (IsFirstWindow) {
for (i = (Primes[(NextPrime - 1)] + 1); (i
<= (p | (2 - 1))); i++) {
// the range of possible new primes'
// TODO: Warning!!! The operator should be an XOR ^ instead of an OR, but not available in CodeDOM
// Dont go beyond the first window'
if ((i >= WindowSize)) {
break;
}
if ((Numbers[(i - FirstValue)].UnFactored == i)) {
// this is a prime less than SQRT(N), so add it to the list.'
Primes[NextPrime] = i;
NextPrime++;
}
}
}
PrimeIndex++;
// move to the next prime'
}
// Now Finish & Emit each one'
for (k = FirstValue; (k <= LastValue); k++) {
// With...
// Primes larger than Sqrt(N) will not be finished: '
if ((Numbers[(k - FirstValue)].UnFactored > 1)) {
// Not done factoring, must be an large prime factor remaining: '
(Numbers[(k - FirstValue)].Phi * (Numbers[(k - FirstValue)].UnFactored - 1).UnFactored) = 1;
Numbers[(k - FirstValue)].Phi = 1;
}
// Emit the value pair: (k, Phi(k)) '
this.EmitPhi(k, Numbers[(k - FirstValue)].Phi);
}
// re-Map to the next window '
IsFirstWindow = false;
this.MapWindow(NextFirst);
}
}
void EmitPhi(long k, long Phi) {
// just a placeholder for now, that raises an event to the display form'
// periodically for reporting purposes. Change this to do the actual'
// emitting.'
if (((k % ReportInterval)
== 0)) {
EmitTotientPair(k, Phi);
}
}
public void MapWindow(long FirstVal) {
// Efficiently reset the window so that we do not have to re-allocate it.'
// init all of the boundary values'
FirstValue = FirstVal;
PrevLast = (FirstValue - 1);
NextFirst = (FirstValue + WindowSize);
LastValue = (NextFirst - 1);
// Initialize the Numbers prime factor arrays'
long i;
for (i = 0; (i
<= (WindowSize - 1)); i++) {
// With...
// initially equal to the number itself'
Phi = 1;
// starts at mulplicative identity(1)'
}
}
long PiMax(long x) {
// estimate of pi(n) == {primes <= (n)} that is never less'
// than the actual number of primes. (from P. Dusart, 1999)'
return ((x / Log(x)) * (1 + (1.2762 / Log(x))));
}
}
다른 팁
아무도 Phi (k)를 계산하는 더 빠른 방법을 찾지 못했습니다 (일명, Euler의 totient 기능) 먼저 k의 주요 요인을 찾는 것보다. 세계 최고의 수학자들은이 문제를 해결하는 더 빠른 방법을 찾는 것이 약점을 창출하기 때문에 1977 년부터 문제에 많은 CPU주기를 던졌습니다. RSA 공개 키 알고리즘. (RSA의 대중과 개인 키는 모두 phi (n)를 기준으로 계산되며, 여기서 n은 두 개의 큰 프라임의 산물입니다.)
PHI (k)의 계산은 K의 주요 인자화를 사용하여 수행되어야하며, 이는 유일한 현명한 방법입니다. 그것에 새로 고침이 필요한 경우 Wikipedia는 공식을 운반합니다.
이제 1에서 큰 n 사이의 모든 숫자의 모든 프라임 디바이저를 계산해야한다면, 결과를보기 전에 노년기로 죽을 것입니다. 가능한 주요 요인, 즉 모든 프라임은 N보다 작거나 동일합니다.
따라서 귀하의 문제는 숫자의 모든 디바이저를 계산하는 것과 비슷할 것입니다. 만 인수 분해에서 특정 프라임을 사전에 찾을 수있는 최대 횟수가 무엇인지 모릅니다. Python 목록에서 Tim Peters가 원래 작성한 반복자 조정 (Something Something 나는 블로그를 작성했다...) totient 함수를 포함시키기 위해, K, phi (k) 쌍을 생성하는 Python에서 가능한 구현은 다음과 같습니다.
def composites(factors, N) :
"""
Generates all number-totient pairs below N, unordered, from the prime factors.
"""
ps = sorted(set(factors))
omega = len(ps)
def rec_gen(n = 0) :
if n == omega :
yield (1,1)
else :
pows = [(1,1)]
val = ps[n]
while val <= N :
pows += [(val, val - pows[-1][0])]
val *= ps[n]
for q, phi_q in rec_gen(n + 1) :
for p, phi_p in pows :
if p * q > N :
break
else :
yield p * q, phi_p * phi_q
for p in rec_gen() :
yield p
n 아래의 모든 주요 요소를 계산하는 데 도움이 필요하면 그것에 대해 블로그를 작성했습니다... 명심하십시오.12 그 자체로는 꽤 놀라운 업적입니다 ...
이것은 Project Euler 245에서 나온 것입니까? 나는 그 질문을 기억하고 그것을 포기했다.
내가 Totient를 계산하는 가장 빠른 방법은 K에 반복되는 요소가 없다는 점을 감안할 때 주요 요인 (P-1)을 함께 곱하는 것이 었습니다 (문제를 정확하게 기억하면 결코 그렇지 않았습니다).
따라서 요인을 계산하려면 아마도 사용하는 것이 가장 좋습니다. gmpy.next_prime 또는 pyecm (타원 곡선 인수).
Jaime이 제안한대로 주요 요인을 체질 할 수도 있습니다. 최대 10의 숫자12, 최대 원수는 백만 미만이므로 합리적입니다.
요인화를 메모하면 PHI 기능 속도가 빨라질 수 있습니다.
이러한 종류의 문제 내가 사용하여 반복기는 반환에 대한 각 정수 m < N 목록 소수 < sqrt(N 는)분할 m.을 구현한 반복기를 나는 배열을 사용하여 A 의 길이 R 가 R >sqrt(N).각 지점에서 배열 A 담고 목록을의 소수가 그를 나눌 m..m+R-1 입니다.I.e. A[m % R]소수를 포함 분할 m.각각의 프라임 p 에 정확히 하나의 목록,즉에 A[m % R 에 대한 가장 작은 범위의 정수 m .. m+R-1 는 의해 나눌 p.생성할 때 다음의 요소를 반복기 단순히 목록 A[m % R]이 반환됩니다.다음 목록의 소수에서 제거됩 A[m % R]그리고 각각의 프라임 p 가 추가됩 A[(m+p)%R].
의 목록과 함께 소수 < sqrt(N)나누어 m 그것은 쉽게 찾을 수 있의 인수 분해 m, 가 있기 때문에,대부분에서 하나의 주요보다 큰 sqrt(N).
이 방법은 복잡도 O(N 로그인(log(N)))가정에서는 모두 작업 목록을 포함하여 작업을 O(1).메모리 요구 사항 O(sqrt(N)).
불행하게,일부한 지속적인 오버헤드,여기에서 찾고 있었는 더 우아한 방법을 생성하는 값피(n)지만,그래서 나는 성공적이었다.
효율적인 파이썬 생성기가 있습니다. 경고는 결과를 순서대로 생성하지 않는다는 것입니다. 그것은 기반입니다 https://stackoverflow.com/a/10110008/412529 .
메모리 복잡성은 O (log (n))는 한 번에 단일 숫자에 대한 주요 요인 목록 만 저장하면됩니다.
CPU 복잡성은 O (n log log n)와 같은 간신히 초고선형입니다.
def totientsbelow(N):
allprimes = primesbelow(N+1)
def rec(n, partialtot=1, min_p = 0):
for p in allprimes:
if p > n:
break
# avoid double solutions such as (6, [2,3]), and (6, [3,2])
if p < min_p: continue
yield (p, p-1, [p])
for t, tot2, r in rec(n//p, partialtot, min_p = p): # uses integer division
yield (t*p, tot2 * p if p == r[0] else tot2 * (p-1), [p] + r)
for n, t, factors in rec(N):
yield (n, t)
나는 당신이 다른 길로 갈 수 있다고 생각합니다. PHI (k)를 얻기 위해 정수 K를 고려하는 대신, 프라임에서 1에서 N에서 N에서 모든 정수를 생성하고 Phi (k)를 빨리 얻을 수 있습니다. 예를 들어, pN N보다 작은 최대 프라임입니다.
피1 나 1 * p2 나 2 * ... * pN 나 N 각각의 위치제이 0에서 [log (n) / log (p까지 실행제이)] ([]는 바닥 함수입니다).
이렇게하면 PHI (K)를 빠르게 비싼 주요 요인화로 얻을 수 있으며 여전히 1과 N 사이의 모든 k를 반복 할 수 있습니다 (순서대로는 아니지만 주문에 관심이 없다고 생각합니다).
totients를 체로 체로 체결하십시오 N:
(define (totients n)
(let ((tots (make-vector (+ n 1))))
(do ((i 0 (+ i 1))) ((< n i))
(vector-set! tots i i))
(do ((i 2 (+ i 1))) ((< n i) tots)
(when (= i (vector-ref tots i))
(vector-set! tots i (- i 1))
(do ((j (+ i i) (+ i j))) ((< n j))
(vector-set! tots j
(* (vector-ref tots j) (- 1 (/ i)))))))))
이것은 p & q가 프라임 인 n = pq를 고려합니다.
Elixir 또는 Erlang에서 잘 작동합니다.
의사 랜덤 시퀀스에 대해 다른 발전기를 시도 할 수 있습니다. x*x + 1 일반적으로 사용됩니다.
이 라인 : defp f0(x, n), do: rem((x * x) + 1, n)
다른 가능한 개선 사항 : 더 나은 또는 대안 GCD, rem 그리고 ABS 기능
defmodule Factorizer do
def factorize(n) do
t = System.system_time
x = pollard(n, 2_000_000, 2_000_000)
y = div(n, x)
p = min(x, y)
q = max(x, y)
t = System.system_time - t
IO.puts "
Factorized #{n}: into [#{p} , #{q}] in #{t} μs
"
{p, q}
end
defp gcd(a,0), do: a
defp gcd(a,b), do: gcd(b,rem(a,b))
defp pollard(n, a, b) do
a = f0(a, n)
b = f0(f0(b, n), n)
p = gcd(abs(b - a), n)
case p > 1 do
true -> p
false -> pollard(n, a, b)
end
end
defp f0(x, n), do: rem((x * x) + 1, n)
end
