４台クラスターで、高次方程式を解く I (重解なし)

前回、Raspberry Pi ４台のクラスターで単純ループの計算時間が1/4に短縮できる結果が得られました。この４台のクラスターで何か計算できないかと考えました。素数計算、multiprocessingなどと候補はあったのですが、高次方程式を解くことにしました。今まで高次方程式を解く時は、エクセルを使って散布図・折れ線グラフを描き、x軸付近をズームアップすることにより解を得ていました。Pythonに入門してプログラムで求めてみます。

方程式を解くならSciPyが最適だと知りました。NumPyとSciPyを追加インストールしました。ところがSciPyの使い方が分かりませんでした。結局、またPythonプログラムらしくない、Pythonプログラムを自分で作ることにしました。今回は１０次方程式を解きます。プログラム作成途上なので、重解の無い、因数分解の可能な式を使います。（有効数字と桁数の議論はあるようですが今回はふれません）下表に、今回作成した１０次方程式を解くプログラムの計算時間を示しました。

|機種|Raspberry Pi 4B|Raspberry Pi 4B|Raspberry Pi 4B|
|:--|:---:|:---:|:--:|
|台数|1台|1台|4台のクラスター|
|スレッド数|1 スレッド|4 スレッド|1 スレッド（各Pi）|
|処理時間|52秒|53秒|13秒|
|importモジュール|time|time　threading|time　socket|

最初に１台で１スレッドと４スレッドを比べましたが52秒と53秒と差はありません。４スレッドにすると少し複雑なプログラムになりますが、計算時間の短縮効果はありません。（前回は例外的に半分にループ計算時間が短縮されたのですが）今回は１スレッドのプログラムを作ることにしました。Raspberry Pi ４台のクラスターで、各１スレッドの１０次方程式を解くプログラムを作成しました。結果は、上表のように52秒→13秒と計算時間を1/4に短縮することができました。４台のクラスター効果がそのまま出ています。有効数字６桁の10次方程式の関数の計算ループが極端な律速なので、複数台のクラスターのLANターンオーバーによる負の影響はまだ見出せません。

#重解の無い３次方程式の例
10次方程式を解く前に、重解が無く、答えが分かっている３次方程式を解くプログラムを作ってみました。
y = 2.1 - 2.1425x - 2.3x^2 + x^3 という３次方程式です。

``` python
eq-3ji.py
import time

# 3次方程式の係数（重解なし）
a0 =  2.1
a1 = -2.1425
a2 = -2.3
a3 =  1.

answerszero = []

def y(R_f, R_w, X_init, X_end, A0, A1, A2, A3):
    for x0 in range(X_init, X_end):
        X = x0/R_f
        Y = A0 + A1*X**1 + A2*(X)**2 + A3*(X)**3
        if Y == 0.0:
            answerszero.append((X, Y))
        elif 0<Y and Y<R_w:
            answerszero.append((X, Y))
        elif 0>Y and Y>-R_w:
            answerszero.append((X, Y))

initial_time=time.time()
r_w=0.00002#候補少な過ぎor多過ぎ

x_init=500000
x_end=2000000
r_f=-1000000
y(r_f,r_w,x_init,x_end,a0,a1,a2,a3)

x_init=500000
x_end=3000000
r_f=1000000
y(r_f,r_w,x_init,x_end,a0,a1,a2,a3)

print(f"Time :")
print(str(time.time() - initial_time))
print(f"Answers : (x, y)")
print(answerszero)
time.sleep(10)
# Time :
# 3.85....
# Answers : (x, y)
# (-1.150175, 4.698....e-06), (-1.150176, -2.418....e-06)
# (0.652651, 1.443....e-06), (0.652652, -2.422....e-06)
# (2.797524, -2.428....e-06), (2.797525, 6.038....e-06)
```
解前後の候補をピックアップしています。
xが-1.150175から-1.150176でy値が正から負へと変わっています。xが0.652651から0.652652でy値が正から負へと変わっています。xが2.797524から2.797525でy値が負から正へと変わっています。重解が無いので、この３点で、３次曲線がx軸を横切っています。つまりy = 2.1 - 2.1425x - 2.3x^2 + x^3 = 0 という３次方程式の３個の解です。
--1.15017，0.65265，2.79752

今回は目視で最終的な解を得ました。この解を自動判定するアルゴリズムがNewton法などです。これらは、次回のIIで重解がある高次方程式を考えていますので、その時に一緒にやってみる予定です。今回は目視で。

#重解の無い10次方程式へのステップ
◇２次方程式 (解の公式あり)
(x-1)(x-2) =
x**2 - 3*x + 2 = 0 
◇３次方程式 (解の公式あり)
(x-1)(x-2)(x-3) = 
(x-1) * (x**2 - 5* x + 6) =
x**3 - 6*x**2 + 11* x - 6 = 0 
◇４次方程式 (解の公式あり)
(x-1)(x-2)(x-3)(x-4) = 
(x**2 - 3*x + 2) * (x**2 - 7* x + 12) =
x**4 - 10*x**3 + 35*x**2 - 50* x + 24 = 0
◇５次方程式
(x-1)(x-2)(x-3)(x-4)(x-5) = 
(x**4 - 10*x**3 + 35*x**2 - 50* x + 24)(x-5) =
x**5 - 15*x**4 + 85*x**3 - 225*x**2 + 274* x - 120 = 0
◇６次方程式 (6次式までエクセル近似あり)
(x-1)(x-2)(x-3)(x-4)(x-5)(x-6) = 
(x**5 - 15*x**4 + 85*x**3 - 225*x**2 + 274* x - 120)(x-6) =
x**6 - 21*x**5 + 175*x**4 - 735*x**3 +1624*x**2 - 1764* x + 720 = 0
◇７次方程式
(x-1)(x-2)(x-3)(x-4)(x-5)(x-6)(x-7) = 
(x**6 - 21*x**5 + 175*x**4 - 735*x**3 +1624*x**2 - 1764* x + 720)(x-7) =
x**7- 28*x**6 + 322*x**5 - 1960*x**4 + 6769*x**3 - 13132*x**2 + 13068*x - 5040 = 0
◇８次方程式
(x-1)(x-2)(x-3)(x-4)(x-5)(x-6)(x-7)(x-8) = 
(x**7- 28*x**6 + 322*x**5 - 1960*x**4 + 6769*x**3 - 13132*x**2 + 13068*x - 5040)(x-8) =
x**8 - 36*x**7+ 546*x**6 - 4536*x**5 + 22449*x**4 - 67284*x**3 + 118124*x**2 - 109584* x + 40320 = 0
◇９次方程式
(x-1)(x-2)(x-3)(x-4)(x-5)(x-6)(x-7)(x-8)(x-9) = 
(x**8 - 36*x**7+ 546*x**6 - 4536*x**5 + 22449*x**4 - 67284*x**3 + 118124*x**2 - 109584* x + 40320 )(x-9) =
x**9- 45*x**8 + 870*x**7 - 9450*x**6 + 63273*x**5 - 269325*x**4 + 723680*x**3 - 1172700*x**2 + 1026576*x - 362880 = 0
◇１０次方程式
(x-1)(x-2)(x-3)(x-4)(x-5)(x-6)(x-7)(x-8)(x-9)(x-10) = 
(x**9- 45*x**8 + 870*x**7 - 9450*x**6 + 63273*x**5 - 269325*x**4 + 723680*x**3 - 1172700*x**2 + 1026576*x - 362880)(x-10) =
x**10 - 55*x**9 + 1320*x**8 - 18150*x**7 + 157773*x**6 - 902055*x**5 + 3416930*x**4 - 8409500*x**3 + 12753576*x**2 - 10628640*x + 3628800 = 0

10次方程式までの道のりは長かったのですが、どうにか計算用の係数が得られました。そして10次方程式を解くプログラムを作りました。有効数字を小数点以下６桁→７桁に上げましたので、1台のRaspberry Pi 4Bで480秒です。約8分です。ヒートシンクが温かくなりました。

``` python
eq-10ji-option-jit.py
import time
#from numba import njit

answersX  = []
answersY  = []

#@njit(cache=True)
def calc_2s(thread_NO,R_F,R_W,X_init,X_end,A0,A1,A2,A3,A4,A5,A6,A7,A8,A9,A10):
    totalsX=[6]
    totalsY=[6]
    for x in range(X_init, X_end):
        X = x/R_F
        Y =A0+A1*(X)**1+A2*(X)**2+A3*(X)**3+A4*(X)**4+A5*(X)**5+A6*(X)**6+A7*(X)**7+A8*(X)**8+A9*(X)**9+A10*(X)**10
        if Y == 0.0:
            totalsX.append(X)
            totalsY.append(Y)
        elif 0 < Y and Y <  R_W:
            totalsX.append(X)
            totalsY.append(Y)
        elif 0 > Y and Y > -R_W:
            totalsX.append(X)
            totalsY.append(Y)
    return totalsX, totalsY

n=0
print(f"pre-cache")
x_init=0
x_end=5
r_f=1
r_w=0.000000001
a0=1.0
a1=-3.0
a2=2.0
a3=0.
a4=0.0
a5=0.0
a6=0.0
a7=0.0
a8=10.0
a9=0.0
a10=0.0
answersX, answersY = calc_2s(n,r_f,r_w,x_init,x_end,a0,a1,a2,a3,a4,a5,a6,a7,a8,a9,a10)
answersX.clear()
answersY.clear()
print(f"Calculation Start")
x_init=0
x_end=120000001
r_f=10000000
r_w=0.011
a0=3628800.0
a1=-10628640.0
a2=12753576.0
a3=-8409500.0
a4=3416930.0
a5=-902055.0
a6=157773.0
a7=-18150.0
a8=1320.0
a9=-55.0
a10=1.0
initial_time=time.time()
answersX, answersY = calc_2s(n,r_f,r_w,x_init,x_end,a0,a1,a2,a3,a4,a5,a6,a7,a8,a9,a10)
print(f"Calculation Finished")
answersX.remove(answersX[0])
answersY.remove(answersY[0])
answersXY=[]
for i in range(len(answersX)):
    answersXY.append((answersX[i], answersY[i]))
print(f"Time :")
print(str(time.time() - initial_time))
print(f"Answers :")
print(answersXY)
time.sleep(10)
```

##基本オプション @njit
プログラム中の'#'を外すと、numbaの@njitを使ってRaspberry Piでプログラムが動きます。3.9秒です。計算時間は@njit無と比較すると、1/120に短縮されます。@njitは、Raspberry Pi OS に最初から入っている関数の一部に対してのみの対応ですが、今後さらに対応が拡大されるようです。
Raspberry Pi 4Bの各モジュールのVersion (2020DEC24)
OS : Raspberry Pi OS with desktop : 2020-12-02-raspios-buster-armhf.img
Numpy : 1.19.4
LLvmlite : 0.27.0
Numba : 0.39.0
Cython : 0.29.21
(Scipy : 1.5.4)
* 追加モジュールのRaspberryPiへのインストールは下記のようにupgradeを繰り返しました。結果として上記のようなVersionになりました。
sudo apt-get install libatlas-base-dev    (y)
sudo apt install python3-numpy
sudo pip3 install scipy==1.5.4   (これでscipy-1.5.4)
pip install matplotlib 　　 (Failed出るが次へ)
pip3 install --upgrade matplotlib  (Failed出るが次へ)
sudo apt-get install build-essential
sudo apt-get install llvm  (y)　(Failed出るが次へ)
pip3 install llvmpy    　(Failed出るが次へ)
pip3 install llvmlite      (Failed出るが次へ)
sudo apt install python3-llvmlite　(Failed出るが次へ)
pip3 install cython　　　(Failed出るが次へ)
pip3 install numba　　　(Failed出るが次へ)
sudo pip3 install --upgrade numpy  (●numpy-1.19.4)
sudo pip3 install --upgrade llvmlite  (●Failed出るが次へ)
sudo pip3 install --upgrade numba　(Failed出るが次へ)
(Thonny warning msg : require 0.3.9)
pip3 install numba==0.39.0  (●完了numba-0.39.0)
pip3 install llvmlite==0.27.0

最後にwarning msgに従って、各モジュールのVersionを合わせたらRaspberry Pi 4Bで＠jitが偶然に動き出したという状況です。またnumbaはrequired : numpy, llvmliteと表示されてますが、Windows 10(x64)ではVersionを合わせても動きませんでした。Windows 10(x64)は2020DEC25現在、Python 3.8＋AnacondaのJupiter Notebookなら@njitがすぐ動きだしたのですが。参考にRaspberry Pi OS with desktop (2020DEC02) のRaspberry Pi 4Bでのpip3 listを添付します。

``` python
pi@raspberrypi:～ $ pip3 list 
Package　　　Version
- - - - -    - - - - - 
asnicrypto　　　0.24.0 
astroid　　　　2.1.0 
asttokens 　　　1.1.13 
automationhat 　　0.2.0 
beautifulsoup4 　　4.7.1 
blinker 　　　　1.4 
blinkt 	　　　　0.1.2 
buttonshim 　　　0.0.2
cairocffi　　　　1.2.0
capixxx 　　　　0.1.3 
certifi 　　　　2018.8.24 
cffi           1.14.4
chardet 　　　　3.0.4 
Click 	　　　　7.0 
colorama 　　　0.3.7
colorzero 　　　1.1
cookies 　　　　2.2.1 
cryptography 　　2.6.1 
cupshelpers 　　1.0 
cycler 　　　　0.10.0 
cython 　　　　0.29.21 
decorator 　　　　4.3.0 
docutils　　　　0.14
drumhat　　　0.1.0 
entrypoints　　　0.3 
envirophat　　　1.0.0 
Explorer HAT 　　　0.4.2 
Flask　　　　　1.0.2 
fourletterphat 　　0.1.0 
gpiozero　　　　1.5.1 
html5lib　　　　1.0.1 
idna　　　　2.6 
isort　　　　　4.3.4 
itsdangerous　　　0.24 
jedi　　　　0.13.2 
Jinja2　　　　2.10 
keyring　　　　17.1.1 
keyrings.alt 　　　　3.1.1 
kiwisolver　　　　1.3.1 
lazy-object-proxy 　　1.3.1 
llvmlite　　　　0.27.0 
logilab-common 　　　1.4.2 
lxml　　　　　4.3.2 
MarkupSafe　　　　1.1.0 
matplolib         3.3.3
mccabe　　　　　　　0.6.1 
microdotphat 　　　0.2.1 
mote　　　　0.0.4
motephat 　　　0.0.3 
mypy 　　　　0.670 
mypy-extensions 　0.4.1 
numba 　　　0.39.0 
numpy 　　　1.19.4 
oauthlib 　　　2.1.0 
olefile 　　　0.46 
pantilthat 　　0.0.7 
parso 　　　0.3.1 
pexpect 　　4.6.0 
pgzero 　　　1.2 
phatbeat 　　0.1.1 
pianohat 　　0.1.0 
picamera 　　1.13 
piglow 　　　1.2.5 
pigpio 　　　1.44 
Pillow 　　　8.0.1 
pip 　　　　18.1 
psutil 　　　5.5.1 
pycairo 　　1.16.2
pycparser    2.20
pycrypto 　　2.6.1 
pycups　　　1.9.73 
pygame　　　1.9.4.post1 
Pygments　　　2.3.1 
PyGObject　　　3.30.4 
pyinotify　　　0.9.6 
PyJWT　　　1.7.0 
pylint　　　2.2.2 
pyOpenSSL　　　19.0.0 
pyparsing　　　2.4.7 
pyserial　　　3.4 
pysmbc　　　　1.0.15.6 
python-apt　　　1.8.4.1 
python-dateutil 　　2.8.1 
pyxdg　　　　0.25 
rainbowhat　　0.1.0 
report lab　　　　3.5.13 
requests　　　　2.21.0 
requests-oauthlib 　　1.0.0 
responses　　　　0.9.0 
roman　　　　2.0.0 
RPi.GPIO　　　　0.7.0 
RTIMULib　　　　7.2.1 
scipy　　　　1.5.4 
scrollphat　　　0.0.7
scrollphathd 　1.2.1 
Secret Storage 　2.3.1 
Send2Trash　　1.5.0 
sense-hat　　2.2.0 
setuptools　　40.8.0 
simplejson　　3.16.0 
six　　　　1.12.0 
skywriter　　0.0.7 
sn3218　　　1.2.7 
soupsieve　　1.8 
spidev　　　3.4 
ssh-import.id 　5.7 
thonny　　　3.3.0 
touchphat　　0.0.1 
twython　　　3.7.0 
typed-ast　　1.3.1 
unicornhathd 　0.0.4 
urllib3　　　1.24.1 
webencodings 　　0.5.1 
Werkzeug　　　0.14.1 
wheel　　　0.32.3 
wrapt　　　1.10.11 
pi@raspberrypi:～ $
```

##基本オプション multiprocessing
``` python
eq-10ji-multiprocessing-option-jit.py
import time
from multiprocessing import Pool
#from numba import njit

#@njit(cache=True)
def calc_2s(n0,R_f,R_w,X_init,X_end,A0,A1,A2,A3,A4,A5,A6,A7,A8,A9,A10):
    totalsXY=[6]
    for x0 in range(int(X_init), int(X_end)):
        X = x0/R_f
        Y = A0+A1*(X)**1+A2*(X)**2+A3*(X)**3+A4*(X)**4+A5*(X)**5+A6*(X)**6+A7*(X)**7+A8*(X)**8+A9*(X)**9+A10*(X)**10
        if Y == 0.0:
            totalsXY.append(X)
            totalsXY.append(Y)
        elif 0 < Y and Y <  R_w:
            totalsXY.append(X)
            totalsXY.append(Y)
        elif 0 > Y and Y > -R_w:
            totalsXY.append(X)
            totalsXY.append(Y)
    return totalsXY


def wrapper_func(args):
    return calc_2s(*args)

def multi_processing_pool(Args_Lists):
    p = Pool(4)
    process_XY = p.map(wrapper_func, Args_Lists)
    p.close()
    return process_XY

if __name__=="__main__":


    print(f"pre-cache")
    n=0
    x_Range_to   = 12
    r_f = 10
    r_w = 0.1
    f_range = x_Range_to * r_f
    a0=1.
    a1=0.
    a2=0.
    a3=0.
    a4=0.
    a5=0.
    a6=0.
    a7=0.
    a8=0.
    a9=0.
    a10=0.
    ArgsLists= [(n,r_f,r_w,i*f_range/4,(i+1)*f_range/4,a0,a1,a2,a3,a4,a5,a6,a7,a8,a9,a10) for i in range(4)]
    processXY = multi_processing_pool(ArgsLists)
    processXY.clear()
    
    print(f"Calculation Start")
    x_Range_to   = 12
    r_f = 10000000
    r_w = 0.011
    f_range = x_Range_to * r_f
    a0 = 3628800.0
    a1 =-10628640.0
    a2 = 12753576.0
    a3 =-8409500.0
    a4 = 3416930.0
    a5 =-902055.0
    a6 = 157773.0
    a7 =-18150.0
    a8 = 1320.0
    a9 =-55.0
    a10= 1.0
    initial_time=time.time()
    ArgsLists= [(n,r_f,r_w,i*f_range/4,(i+1)*f_range/4,a0,a1,a2,a3,a4,a5,a6,a7,a8,a9,a10) for i in range(4)]
    processXY = multi_processing_pool(ArgsLists)
    print(f"Calculation Finished")
   
    processXYstr = str(processXY)
    processXYstr = processXYstr.replace('[[', '')
    processXYstr = processXYstr.replace(']]', '')
    processesXYstr=processXYstr.split('], [')
    process0XYstr=processesXYstr[0]
    process1XYstr=processesXYstr[1]
    process2XYstr=processesXYstr[2]
    process3XYstr=processesXYstr[3]
    process0sXYstr=process0XYstr.split(', ')
    process0sXYstr.remove(process0sXYstr[0])
    process1sXYstr=process1XYstr.split(', ')
    process1sXYstr.remove(process1sXYstr[0])
    process2sXYstr=process2XYstr.split(', ')
    process2sXYstr.remove(process2sXYstr[0])
    process3sXYstr=process3XYstr.split(', ')
    process3sXYstr.remove(process3sXYstr[0])
    answersXY=[]
    for i in range(int(len(process0sXYstr)/2)):
        answersXY.append((float(process0sXYstr[2*i]), float(process0sXYstr[2*i+1])))
    for i in range(int(len(process1sXYstr)/2)):
        answersXY.append((float(process1sXYstr[2*i]), float(process1sXYstr[2*i+1])))
    for i in range(int(len(process2sXYstr)/2)):
        answersXY.append((float(process2sXYstr[2*i]), float(process2sXYstr[2*i+1])))
    for i in range(int(len(process3sXYstr)/2)):
        answersXY.append((float(process3sXYstr[2*i]), float(process3sXYstr[2*i+1])))

    print(f"Time :")
    print(str(time.time() - initial_time))
    print(f"Answers :")
    print(answersXY)

```

4coreのRaspberry Pi 4Bで動くマルチプロセスのプログラムです。4coreすべて使って稼働します。計算時間は1台のRaspberry Pi 4Bで120秒です。シングルプロセスが1台のRaspberry Pi 4Bで1coreのみ稼働で480秒でしたから、1/4に計算時間が短縮されます。さらにプログラム中の'#'を外すと、numbaの@njitを使ってRaspberry Piでプログラムが動きます。Raspberry Pi 4Bで1.4秒です。*argsなどpythonでは見慣れないC言語由来らしき呪文のような２行がありますが、関数名calc_2s, wrapper_func以外はcopyして置いておきます。呪文ですから(もともとwrapper関数はmapで複数引数が使えるようになるものなのですがjitの実体がよくわからないので)。これでマルチプロセスが@njitの影響を受けなくなります。計算時間は、シングルプロセスの@njit有3.9秒と比較すると、マルチプロセスの@njit有1.4秒は、約1/3に短縮されます。4coreなので1/4短縮を予想していましたが、マルチプロセス処理部分は@njit化に対応してないので仕方ありません。それでもマルチプロセス化で約1/3に短縮されました。全く高速化をしてない@njit無のシングルプロセスの480秒と比較すると、1/340に計算時間が短縮されます。

次に4台クラスターでの並列処理のプログラムです。ほとんど前回のプログラムを移植したので、前回のように「解決不可能なエラー」というものはありませんでした。（次回のII 重解の有る場合は少し難しそうですが）今回作成したプログラムを次に示します(まだ、numbaの@jitでも、multiprocessingでも、動きません)。

#作成した重解の無い10次方程式用の４クラスターのpythonプログラム
cluster2s-21.pyをPiフォルダーに置いて、python3 で実行して待機。他の３台はcluster2s-21.pyをコピペし、それぞれ 22, 23, 24 へ変更する。同様にPiフォルダーに置いて、python3 で実行して待機。今回も、pythonプログラムらしくない、pythonプログラムなのですが、公開します。

``` python
cluster2s-21.py
import time
import socket

######  Change manually 2 Points to program NO  ########

answerszero = []

def calc_2s(n0,R_f,R_w,X_init,X_end,A0,A1,A2,A3,A4,A5,A6,A7,A8,A9,A10):
    print(f"START_CPU")
    for x0 in range(X_init, X_end):
        X = x0/R_f
        Y = A0+A1*(X)**1+A2*(X)**2+A3*(X)**3+A4*(X)**4+A5*(X)**5+A6*(X)**6+A7*(X)**7+A8*(X)**8+A9*(X)**9+A10*(X)**10
        if Y == 0.0:
            answerszero.append((X,Y))
        elif 0 < Y and Y <  R_w:
            answerszero.append((X,Y))
        elif 0 > Y and Y > -R_w:
            answerszero.append((X,Y))
    print(f"END_OF_CPU")

socket00 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
socket00.bind(('192.168.2.21', 60021))######## Point 1 : Change IP address and port to program NO ########
socket00.listen(5)
clientsocket, address = socket00.accept()
liststrEncode = b''
liststrEncode = clientsocket.recv(1016)
liststr=liststrEncode.decode("utf-8")
print(liststr)
# [[3628800.0, -10628640.0, 12753576.0, -8409500.0, 3416930.0, -902055.0, 157773.0, -18150.0, 1320.0, -55.0, 1.0],\
# [2, 1000000, 0.01, 0, 3000001]]
liststr = liststr.replace('[[', '')
liststr = liststr.replace(']]', '')
listsstr=[]
listsstr=liststr.split('], [')
list00str=listsstr[0]
listCLUSTERstr=listsstr[1]

lists00str=[]
lists00str=list00str.split(', ')
a0 =float(lists00str[0])
a1 =float(lists00str[1])
a2 =float(lists00str[2])
a3 =float(lists00str[3])
a4 =float(lists00str[4])
a5 =float(lists00str[5])
a6 =float(lists00str[6])
a7 =float(lists00str[7])
a8 =float(lists00str[8])
a9 =float(lists00str[9])
a10=float(lists00str[10])

listsCLUSTERstr=[]
listsCLUSTERstr=listCLUSTERstr.split(', ')
calc_Number=int(listsCLUSTERstr[0])
r_f        =int(listsCLUSTERstr[1])
r_w      =float(listsCLUSTERstr[2])
x_init     =int(listsCLUSTERstr[3])
x_end      =int(listsCLUSTERstr[4])

print(lists00str)
print(listsCLUSTERstr)

clientsocket.send(bytes("OK  Start ....", 'utf-8'))
clientsocket.close()
socket00.close()

if calc_Number == 2:
    initial_time = time.time()
    print(f"CPU is calculating now !")
    n=0
    calc_2s(n,r_f,r_w,x_init,x_end,a0,a1,a2,a3,a4,a5,a6,a7,a8,a9,a10)
    print(f"TIME :")
    print(str(time.time() - initial_time))
    print(f"ANSWER :")
    print(answerszero)

    socket00 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
    socket00.connect(('192.168.2.20', 60021))######## Point 2 : Change only port to program NO ########
    Answerszerostr = str(answerszero)
    socket00.send(bytes(answerszerostr, 'utf-8'))
    msg = b''
    msg = socket00.recv(16)
    socket00.close()
    print(msg.decode("utf-8"))

time.sleep(10)
```

``` python
cluster2s-20total.py
import time
import threading
import socket

## 共通 （10次方程式の係数）
list00 = [3628800.0, -10628640.0, 12753576.0, -8409500.0, 3416930.0, -902055.0, 157773.0, -18150.0, 1320.0, -55.0, 1.0]
## ４台のRaspberry Pi 4B
list21 = [2, 1000000, 0.01, 0, 3000001]
list22 = [2, 1000000, 0.01, 3000001, 6000001]
list23 = [2, 1000000, 0.01, 6000001, 9000001]
list24 = [2, 1000000, 0.01, 9000001, 12000001]

lists21 = []
lists21.append(list00)
lists21.append(list21)
lists22 = []
lists22.append(list00)
lists22.append(list22)
lists23 = []
lists23.append(list00)
lists23.append(list23)
lists24 = []
lists24.append(list00)
lists24.append(list24)

total = 0

## 各Piのtotals
totals21 = []
totals22 = []
totals23 = []
totals24 = []

def connectto_ipports(n):
    if n == 60021:
        print(n)
        socket21 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket21.connect(('192.168.2.21', 60021))
        list21str=''
        list21str=str(lists21)
        socket21.send(bytes(list21str, 'utf-8'))
        msg21 = b''
        msg21 = socket21.recv(16)
        socket21.close()
        print(msg21.decode("utf-8"))
    elif n == 60022:
        print(n)
        socket22 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket22.connect(('192.168.2.22', 60022))
        list22str=''
        list22str=str(lists22)
        socket22.send(bytes(list22str, 'utf-8'))
        msg22 = b''
        msg22 = socket22.recv(16)
        socket22.close()
        print(msg22.decode("utf-8"))
    elif n == 60023:
        print(n)
        socket23 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket23.connect(('192.168.2.23', 60023))
        list23str=''
        list23str=str(lists23)
        socket23.send(bytes(list23str, 'utf-8'))
        msg23 = b''
        msg23 = socket23.recv(16)
        socket23.close()
        print(msg23.decode("utf-8"))
    elif n == 60024:
        print(n)
        socket24 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket24.connect(('192.168.2.24', 60024))
        list24str=''
        list24str=str(lists24)
        socket24.send(bytes(list24str, 'utf-8'))
        msg24 = b''
        msg24 = socket24.recv(16)
        socket24.close()
        print(msg24.decode("utf-8"))

def accept_ports(n):
    if n == 60021:
        print(n)
        socket21 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket21.bind(('192.168.2.20', 60021))
        socket21.listen(5)
        clientsocket21, address21 = socket21.accept()
        total21en = b''
        total21en = clientsocket21.recv(1016)
        total21str = total21en.decode("utf-8")
        totals21.append(total21str)
        print(total21str)
        clientsocket21.send(bytes("Thank you ......", 'utf-8'))
        clientsocket21.close()
        socket21.close()
    elif n == 60022:
        print(n)
        socket22 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket22.bind(('192.168.2.20', 60022))
        socket22.listen(5)
        clientsocket22, address22 = socket22.accept()
        total22en = b''
        total22en = clientsocket22.recv(1016)
        total22str = total22en.decode("utf-8")
        totals22.append(total22str)
        print(total22str)
        clientsocket22.send(bytes("Thank you ......", 'utf-8'))
        clientsocket22.close()
        socket22.close()
    elif n == 60023:
        print(n)
        socket23 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket23.bind(('192.168.2.20', 60023))
        socket23.listen(5)
        clientsocket23, address23 = socket23.accept()
        total23en = b''
        total23en = clientsocket23.recv(1016)
        total23str = total23en.decode("utf-8")
        totals23.append(total23str)
        print(total23str)
        clientsocket23.send(bytes("Thank you ......", 'utf-8'))
        clientsocket23.close()
        socket23.close()
    elif n == 60024:
        print(n)
        socket24 = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
        socket24.bind(('192.168.2.20', 60024))
        socket24.listen(5)
        clientsocket24, address24 = socket24.accept()
        total24en = b''
        total24en = clientsocket24.recv(1016)
        total24str = total24en.decode("utf-8")
        totals24.append(total24str)
        print(total24str)
        clientsocket24.send(bytes("Thank you ......", 'utf-8'))
        clientsocket24.close()
        socket24.close()

threads = []

# Thread Client Start
for n in range(60021, 60025):
    thread = threading.Thread(target=connectto_ipports, args=(n,))
    thread.start()
    threads.append(thread)
# Wait Threads
for thread in threads:
    thread.join()
# Result Display
print("All Clusters Calculation Start")

threads.clear()

# Accept Thread Server Start
for n in range(60021, 60025):
    thread = threading.Thread(target=accept_ports, args=(n,))
    thread.start()
    threads.append(thread)
# Wait Threads
for thread in threads:
    thread.join()

# Result Display
print("Calculation Result:")
print(totals21)
print(totals22)
print(totals23)
print(totals24)

time.sleep(10) # for display

# cluster2s-20total.py 固定IPアドレス 192.168.2.20 (pyを最後に実行、計算開始)
# cluster2s-21.py      固定IPアドレス 192.168.2.21 (pyを先に実行して待機)
# cluster2s-22.py      固定IPアドレス 192.168.2.22 (pyを先に実行して待機)
# cluster2s-23.py      固定IPアドレス 192.168.2.23 (pyを先に実行して待機)
# cluster2s-24.py      固定IPアドレス 192.168.2.24 (pyを先に実行して待機)
## OSは全てRaspberry Pi OS with desktop : 2020-12-02-raspios-buster-armhf.img
```

最後までご覧いただき、ありがとうございました。

次回は重解ありの高次方程式を解く予定です。またQiitaの主旨に同意して投稿していますので、公開したプログラムはコピー、改変などでのご使用は自由です。著作権に関する問題も発生しません。

ただし、Raspberry Pi 4Bを使う場合にはCPUに特に大きめのヒートシンクが必須です。LAN通信の頻度の少ない今回のプログラムのような場合、LANチップは熱くなりません。しかし計算時間が継続するとCPUの激しい発熱で分かるように相当電力を使っています。電源USB Cタイプの後ろにある黒い小さなチップも熱くなりますので、ファンでの風流も必要です。