python中国余数定理_残数系统

间隙

背景知识：我会承认，在几个月前提出这个问题时，我还没有一种方法可以解决这个问题的难处：确定要使用的质数的正确数目。这个网站上有很多非常聪明的人，我真的希望有人会找到一种相当快的方法。但是，由于这没有发生，因此我甚至不确定是否确实有可能解决此问题。因此，我不得不花时间设计一种方法。我相信我所做的一切不会违反这项挑战中的任何规则，当然，我希望对此进行事实检查。

我也对选择代码高尔夫球感到遗憾，因为解决方案的深度比通常适合标记格式的深度略多。话虽如此，为遵循站点规则，本文底部有我解决方案的“ golf”版本。

码

### The first 100 primes;

primes := Primes{[1..100]};

### In many of the functions below, the 'string' variable is a string of digits

###

### Returns the 'index' digit of 'string' as an integer

GetValueAsInt := function(string, index)

return IntChar(string[index]) - 48;

end;

### Used in the 'modulus' function. See that comment for more information.

### Calculates the contribution to the modulus of a digit 'digit' in the 10^power place.

### 'integer' is the modulus

digit_contribution := function(digit, integer, power)

local result, i;

result := 1;

for i in [0..power-1] do

result := ( result * (10 mod integer) ) mod integer;

od;

result := (result * (digit mod integer) ) mod integer;

return result;

end;

### This modulus function is used to calculate the modulus of large numbers without storing them

##### as large numbers.

### It does so by breaking them into digits, and calculating the contribution of each digit.

### e.g. 1234 mod 5 = (1000 mod 5)(1 mod 5) + (200 mod 5)(2 mod 5) + (10 mod 5)(3 mod 5) + (4 mod 5)

### It actually mods after every calculation to ensure that we never get a number larger

##### than the modulus ('integer') squared, which will never be even close to 10^64-1

modulus := function(string, integer)

local i, result, digit, len;

len := Length(string);

result := 0;

for i in [1..len] do

digit := IntChar(string[i]) -48;

result := ( result + digit_contribution(digit, integer, len-i) ) mod integer;

od;

return result;

end;

### This returns the product of the first i-1 primes (mod j). It must not (and does not)

##### ever store an integer larger than 2^64-1

phi_i := function(i,j)

local index, result;

result := 1;

for index in [1..i-1] do

result := ( result * primes[index] ) mod primes[j];

od;

return result;

end;

### Calculates the first residues of 'string' mod the first 100 primes

get_residues := function(string)

local p, result;

result := [];

for p in primes do

Add( result, modulus(string, p) );

od;

return result;

end;

### Gets the ith element in the partial_chinese array, given the previous elements

### See the explanation section and partial_chinese function for more info

get_partial_i := function( i, residues, previous_array )

local index, result;

result := residues[i];

for index in [1..Length(previous_array)] do

result := ( result - previous_array[index]*phi_i(index,i) ) mod primes[i];

od;

result := ( result / phi_i(i,i) ) mod primes[i];

return result;

end;

### returns an array such that the sum of prod_primes(i)*array[i] is equal to the integer value

##### that is represented by the residues. (It basically just does the CRT without

##### actually summing everything.) prod_primes(i) is the product of the first i-1 primes

### See the explanation for a bit more info

### This is what allows us to determine the minimal number of primes to represent a RNS number

partial_chinese := function( string )

local array, i, residues;

residues := get_residues(string);

array := [];

for i in [1 .. Length(primes)] do

Add( array, get_partial_i( i, residues, array ) );

od;

return array;

end;

### Same as partial_chinese but takes input in a different form.

partial_chinese_from_residues := function(residues)

local array, i;

array := [];

for i in [1 .. Length(primes)] do

Add( array, get_partial_i( i, residues, array ) );

od;

return array;

end;

### gives you the number of primes needed to represent an integer. Basically asks how

##### many trailing zeros there are in the chinese array.

get_size := function(string)

local array, i, len, result;

array := partial_chinese(string);

len := Length(array);

for i in [0..len-1] do

if not (array[len-i] = 0) then

return len -i;

fi;

od;

Print("ERROR: get_size().\n");

return 0;

end;

### Same as above but with different input format

get_size_from_residues := function(residues)

local array, i, len, result;

array := partial_chinese_from_residues(residues);

len := Length(array);

for i in [0..len-1] do

if not (array[len-i] = 0) then

return len -i;

fi;

od;

Print("ERROR: get_size().\n");

return 0;

end;

### the actual function. inputs are all strings

f := function(in1, in2, opperation)

local residues_1, residues_2, residues_result, i;

residues_1 := get_residues(in1);

residues_2 := get_residues(in2);

residues_result := [];

if opperation = "+" then

for i in [1..Length(primes)] do

Add( residues_result, ( residues_1[i] + residues_2[i] ) mod primes[i]);

od;

elif opperation = "*" then

for i in [1..Length(primes)] do

Add( residues_result, ( residues_1[i] * residues_2[i] ) mod primes[i]);

od;

elif opperation = "-" then

for i in [1..Length(primes)] do

Add( residues_result, ( residues_1[i] - residues_2[i] ) mod primes[i]);

od;

fi;

Print(residues_1{[1..get_size(in1)]}, " ", opperation, " ", residues_2{[1..get_size(in2)]}, " = ", residues_result{[1..get_size_from_residues(residues_result)]} );

end;

说明：

首先，我们计算两个输入的所有100个残基。我们使用modulus代码中的函数执行此操作。我尝试要小心，以便我们mod在每一步之后都使用内置函数。这样可以确保我们永远不会有大于的数字，该数字540^2比第100个素数平方小1。

拥有所有残差后，我们可以执行给定的操作并mod再次执行每个条目。现在，我们为结果指定了一个唯一的标识符，但是我们需要确定代表结果和每个输入的最小条目数。

到目前为止，弄清我们实际上需要多少残留物是最困难的部分。为了确定这一点，我们执行了中国剩余定理(CRT)的大多数步骤。但是，显然我们必须进行修改，以免最终数字不会太大。

设prod(i)第一个i-1素数之和。例如，

prod(1) = 1

prod(2) = 2

prod(3) = 6

prod(4) = 30

etc

让X是一个整数。让{r_i}是的残基X，即

r_i = X mod p_i

哪里p_i是i个素数。这是针对1

现在，我们将使用CRT找到一个序列{u_i}使得在总和i的prod(i) * u_i等于X。请注意，从u_i技术上讲，每个也是残基，例如u_i < p_i。此外，如果X < prod(i)那么u_i = 0。这是至关重要的。这意味着通过检查尾随零，我们可以确定r_i实际需要X在RNS中表示多少个残基。

如果您希望检查的某些序列u_i，该partial_chinese函数将返回该u_i序列。

通过弄乱CRT，我能够找到这些u_i值的递归公式，从而解决了确定需要多少残基的问题。

公式为：

u_i = [ r_i - SUM ] / prod(i) (mod p_i)

凡SUM超过总和j in [1,i)的u_j * prod(i)。

当然，prod(i)在某些情况下，由于它太大而无法实际计算。为此，我使用了该phi_i功能。此函数返回prod(j) (mod p_i)。它mod步步为营，因此我们永远不会真正计算太大的东西。

如果您好奇该公式的来源，我建议您使用几个CRT示例，可以在Wikipedia页面上找到它们。

最后，对于每个输入以及我们的输出，我们计算u_i序列，然后确定尾随零。然后我们r_i从残基序列的末尾扔掉那么多。

“ Golfed”代码，2621字节

primes:=Primes{[1..100]};GetValueAsInt:=function(string,index)return IntChar(string[index])-48;end;digit_contribution := function(digit, integer, power)local result, i;result:=1;for i in [0..power-1] do result := ( result * (10 mod integer) ) mod integer;od;result:=(result*(digit mod integer) ) mod integer;return result;end;modulus:=function(string, integer)local i,result,digit,len;len:=Length(string);result:=0;for i in [1..len] do digit:= IntChar(string[i])-48;result:=(result+digit_contribution(digit,integer,len-i)) mod integer;od;return result;end;phi_i:=function(i,j)local index,result;result:=1;for index in [1..i-1] do result:=(result*primes[index] ) mod primes[j];od;return result;end;get_residues:=function(string) local p,result;result:=[];for p in primes do Add(result,modulus(string,p));od;return result;end;get_partial_i:=function(i,residues,previous_array)local index,result;result:=residues[i];for index in [1..Length(previous_array)] do result:=(result-previous_array[index]*phi_i(index,i) ) mod primes[i];od;result:=(result/phi_i(i,i)) mod primes[i];return result;end;partial_chinese:=function(string)local array,i,residues;residues:=get_residues(string);array:=[];for i in [1 .. Length(primes)] do Add(array,get_partial_i(i,residues,array));od;return array;end;partial_chinese_from_residues:=function(residues)local array,i;array:=[];for i in [1..Length(primes)] do Add(array,get_partial_i(i,residues,array));od;return array;end;get_size:=function(string)local array,i,len,result;array:=partial_chinese(string);len:=Length(array);for i in [0..len-1] do if not (array[len-i] = 0) then return len-i;fi;od;Print("ERROR: get_size().\n");return 0;end;get_size_from_residues:=function(residues)local array,i,len,result;array:=partial_chinese_from_residues(residues);len:=Length(array);for i in [0..len-1] do if not (array[len-i]=0) then return len-i;fi;od;Print("ERROR: get_size().\n");return 0;end;f:=function(in1,in2,opperation)local residues_1,residues_2,residues_result,i;residues_1:=get_residues(in1);residues_2:=get_residues(in2);residues_result:=[];if opperation = "+" then for i in [1..Length(primes)] do Add(residues_result,(residues_1[i]+residues_2[i] ) mod primes[i]);od;elif opperation = "*" then for i in [1..Length(primes)] do Add(residues_result,(residues_1[i]*residues_2[i])mod primes[i]);od;elif opperation = "-" then for i in [1..Length(primes)] do Add(residues_result,(residues_1[i]-residues_2[i]) mod primes[i]);od;fi;Print(residues_1{[1..get_size(in1)]}, " ", opperation, " ", residues_2{[1..get_size(in2)]}, " = ", residues_result{[1..get_size_from_residues(residues_result)]} );end;