2元１次不定方程式( indefinite equation?)＜4/4＞「2022年度 大学入学共通テスト 本試験｜数学Ⅰ第4問」(配点20点)中国剰余定理関数( 手計算,WolframAlpha,Sage,sympy )でやってみたい。

オリジナル

https://qiita.com/mrrclb48z/items/af08059157cfbce8f0fe

```
整数の性質
(1) 5**4*x-2**4*y-1
(3) 5**5*x-2**5*y-1
(4)11**5*x-2**5*y-1
```
# 解法1:手計算とwolframのPowerModで

手計算なのに？逆元の計算にwolframalphaのPowerModを使いました。

(例)2元1次不定方程式
X=ai(mod mi) i=1,2
X=0 (mod 5^4)
X=1 (mod 2^4)
x=X/5^4

```
(1)
a    m  (M/m)    (M/m)^-1  a*(M/m)*(M/m)^-1   
0  ＜以下省略＞
1  2^4    5^4        1     5^4*1

x= 5^4   /5^4    =   1
y=(5^4-1)/2^4    =  39
x=2^4*k+ 1
y=5^4*k+39

(2)省略

(3)
1  2^5       5^5    29    5^5*29

x= 5^5*29   /5^5 =  29
y=(5^5*29-1)/2^5 =2832
x=2^5*k+  29
y=5^5*k+2832

(4)
1    2^5   11^5      19  11^5*19

x= 11^5*19  /11^5=   19
y=(11^5*19-1)/2^5=95624
x= 2^5*k+   19
y=11^5*k+95624
```
逆元の計算


https://www.wolframalpha.com/input?i=powermod%28+5%5E4%2C-1%2C2%5E4%29%2Cpowermod%285%5E5%2C-1%2C2%5E5%29%2Cpowermod%2811%5E5%2C-1%2C2%5E5%29&lang=ja

結果
{1, 29, 19}


# 解法2: 解法1(手計算とwolframalphaのPowerModで)をvbaで

```vb
 Private Function mul_inv(a As Long, n As Long) As Variant
    If n < 0 Then n = -n
    If a < 0 Then a = n - ((-a) Mod n)
    Dim t As Long: t = 0
    Dim nt As Long: nt = 1
    Dim r As Long: r = n
    Dim nr As Long: nr = a
    Dim q As Long
    Do While nr <> 0
        q = r \ nr
        tmp = t
        t = nt
        nt = tmp - q * nt
        tmp = r
        r = nr
        nr = tmp - q * nr
    Loop
    If r > 1 Then
        mul_inv = "a is not invertible"
    Else
        If t < 0 Then t = t + n
        mul_inv = t
    End If
End Function
Sub myLine(ir, aa As Long, mm As Long, bb As Long)
    Cells(ir, 4) = mul_inv(bb, mm)
    Cells(ir, 5) = aa * bb * Cells(ir, 4)
    Cells(ir, 6) = (Cells(ir, 5) - 1) / mm
End Sub
Sub aaa_mi()
    ActiveSheet.Cells.Clear
    Cells(1, 1) = 1
    Cells(1, 2) = 2 ^ 4
    Cells(1, 3) = 5 ^ 4
    Cells(2, 1) = 1
    Cells(2, 2) = 2 ^ 5
    Cells(2, 3) = 5 ^ 5
    Cells(3, 1) = 1
    Cells(3, 2) = 2 ^ 5
    Cells(3, 3) = 11 ^ 5
    For i = 1 To 3
        Call myLine(i, Cells(i, 1), Cells(i, 2), Cells(i, 3))
    Next
End Sub
'1   16     625   1      625    39
'1   32    3125  29    90625  2832
'1   32  161051  19  3059969 95624


```

# 解法3:WolframAlphaのChineseRemainderで 

(1)ア～ウ

https://www.wolframalpha.com/input?i=ChineseRemainder%5B%7B0%2C1%7D%2C+%7B5%5E4%2C2%5E4%7D%5D-2%5E4y%3D1&lang=ja

y = 39

(1)エ～ク 



https://www.wolframalpha.com/input?i=ChineseRemainder%5B%7B0%2C1%7D%2C+%7B5%5E4%2C2%5E4%7D%2C10*5%5E4%5D-2%5E4y%3D1&lang=ja

解
y = 664

(3)サ～ツ

https://www.wolframalpha.com/input?i=ChineseRemainder%5B%7B0%2C1%7D%2C+%7B5%5E5%2C2%5E5%7D%2C5%5E5*10%5E2%5D-2%5E5y%3D1&lang=ja

解
y = 12207

(4)テ～ノ

https://www.wolframalpha.com/input?i=ChineseRemainder%5B%7B0%2C1%7D%2C+%7B11%5E5%2C2%5E5%7D%5D-2%5E5y%3D1&lang=ja

y = 95624

# 解法4:Sageのcrtで 
```python
var('x,y,k,kk')
def mySage_Solve2(myA,myB,myC,myXmin):
    mymy=mySage_Solve(myA,myB,myC)
    kk=solve(abs(myB)*k+mymy[0]>=myXmin,k)[0][0].rhs()
    kk=math.ceil(float(kk))
    xo=abs(myB)*kk+mymy[0]
    yo=solve(myA*xo+myB*y+myC,y)[0].rhs()
    return xo,yo
def mySage_Solve(myA,myB,myC):
    xo=crt([0,-myC],[myA,-myB])/myA
    yo=solve(myA*xo+myB*y+myC,y)[0].rhs()
    return xo,yo
print("#(1)ア～ウ",mySage_Solve ( 5**4,-2**4,-1)    )
print("#(1)エ～ク",mySage_Solve2( 5**4,-2**4,-1, 10))
print("#(3)サ～ツ",mySage_Solve2( 5**5,-2**5,-1,100))
print("#(4)テ～ノ",mySage_Solve (11**5,-2**5,-1)    )

m =mySage_Solve (5**4,-2**4,-1)[1]
KE=solve(625**2-5**x,x)[0].rhs()
KO=solve(625**2-(2**KE*m**2+2**x*m+1),x)[0].rhs()
print("#(2)ケ～コ",KE,KO)
#(1)ア～ウ (1, 39)
#(1)エ～ク (17, 664)
#(3)サ～ツ (125, 12207)
#(4)テ～ノ (19, 95624)
#(2)ケ～コ 8 5
```

# 解法5:sympyのsolve_congruenceで 
```python
from sympy import *
from sympy.ntheory.modular import solve_congruence
import math
var('x,y,k,kk')
def mySym_Solve2(myA,myB,myC,myXmin):
    mymy=mySym_Solve(myA,myB,myC)
    kk=solve(abs(myB)*k+mymy[0]-myXmin,k)[0]
    kk=math.ceil(float(kk))
    xo=abs(myB)*kk+mymy[0]
    yo=solve(myA*xo+myB*y+myC,y)[0]
    return xo,yo
def mySym_Solve(myA,myB,myC):
    xo=solve_congruence((0, myA),(-myC,-myB))[0]//myA
    yo=solve(myA*xo+myB*y+myC,y)[0]
    return xo,yo
print("#(1)ア～ウ",mySym_Solve ( 5**4,-2**4,-1)    )
print("#(1)エ～ク",mySym_Solve2( 5**4,-2**4,-1, 10))
print("#(3)サ～ツ",mySym_Solve2( 5**5,-2**5,-1,100))
print("#(4)テ～ノ",mySym_Solve (11**5,-2**5,-1)    )

m =mySym_Solve (5**4,-2**4,-1)[1]
KE=solve(625**2-5**x,x)[0]
KO=solve(625**2-(2**KE*m**2+2**x*m+1),x)[0]
print("#(2)ケ～コ",KE,KO)
#(1)ア～ウ (1, 39)
#(1)エ～ク (17, 664)
#(3)サ～ツ (125, 12207)
#(4)テ～ノ (19, 95624)
#(2)ケ～コ 8 5
```

# 参考

https://qiita.com/mrrclb48z/items/ff978d6cd87e05aa486a


(20220511)
ChineseRemainder「整数合同式...を満たすx>=dの最小のxを与える。」ありました。
**原点の移動ができました。「原点の移動」は、世界基準？かもしれません。**

https://reference.wolfram.com/language/ref/ChineseRemainder.html


(20220512)
solve_congruence

https://docs.sympy.org/latest/modules/ntheory.html#sympy.ntheory.modular.solve_congruence