| ＜Python (Object Oriented Programming)＞
[Main] rndPilates.py
import matplotlib.pyplot as plt
import subpilates as subp
import time
w = subp.setVoid()
ww = subp.vScan(0,10,0,6,w.voidx,w.voidy,w.Nvoid)
www = subp.vNeighbor(ww.vposnx, ww.vposny, ww.Nvposn)
wwww = subp.MonteCarlo(10000,100,www.pos2D,www.Neibor,www.Nvp)
ax = plt.axes()
ax.set_ylim(0.0, 0.2)
plt.plot(wwww.x, wwww.y)
plt.show()
w5 = subp.MonteCarloDetail(0,10,0,6,Nanim,ww.vposnx, ww.vposny,
www.pos2D,www.Neibor,www.Nvp)
[Module] subpilates.py
import numpy as np
import random as ra
import numpy.random as rnd
import matplotlib.pyplot as plt
import time
separation = 4
tanimn = 1
class setVoid:
def __init__(self):
voidx = np.array([5, 0, 10])
voidy = np.array([0, 3, 3])
num = len(voidx)
self.voidx = voidx
self.voidy = voidy
self.Nvoid = num
def info(self):
print("'setVoid' gives")
print(self.voidx)
print(self.voidy)
print(self.Nvoid," items of void")
class vScan:
def __init__(self,ix0,ix1,iy0,iy1,voidx,voidy,Nvoid):
def isColln(ix,iy):
q = False
for i in range(Nvoid):
if (voidx[i] == ix) and (voidy[i] == iy):
print('Bumping on (',ix,',',iy,')')
q = True
break
return q
vposnx=[]
vposny=[]
m = 0
for iy in range(iy0,iy1+1):
for ix in range(ix0,ix1+1):
if (ix + iy) % 2 == 1:
if not isColln(ix,iy):
vposnx.append(ix)
vposny.append(iy)
m += 1
self.vposnx = vposnx
self.vposny = vposny
self.Nvposn = m
class vNeighbor:
def __init__(self,vpx,vpy,Nvp):
def isNear(xm,xmm,ym,ymm):
d = (xm - xmm)**2 + (ym - ymm)**2
return d <= separation
Ns = 0
pos2D = np.ndarray(shape=(Nvp,8), dtype=int)
Neibor = []
for m in range(Nvp):
k = 0
s = []
for mm in range(Nvp):
if (m!=mm)and isNear(vpx[m],vpx[mm],vpy[m],vpy[mm]):
s.append(mm+1)
k += 1
Neibor.append(k)
Ns += k
for kk in range(k):
pos2D[m,kk] = s[kk]
self.Neibor = Neibor
self.pos2D = pos2D
self.Ns = Ns
self.Nvp = Nvp
self.prob = float(self.Ns)/float(Nvp*(Nvp-1))
def result(self):
print('Probability is calculated as ',f'{self.prob:,.4f}')
def isNeighbor(ip,jp,Neibor,pos2D):
p = False
for m in range(Neibor[ip]):
if (jp+1) == pos2D[ip,m]:
p = True
break
return p
class MonteCarlo:
def __init__(self,Nloop,Nstep,pos2D,Neibor,Nvp):
Nhit = 0
x = []
y = []
Nsize = 0
for loop in range(Nloop):
m1 = int(Nvp*ra.random())
m2 = int((Nvp-1)*ra.random())
if m1 == m2:
m2 = Nvp - 1
if isNeighbor(m1,m2,Neibor,pos2D):
Nhit += 1
if (loop+1) % Nstep == 0:
y.append(float(Nhit)/float(loop+1))
x.append(loop+1)
print(loop+1 ,': ', y[Nsize])
Nsize += 1
self.x = x
self.y = y
self.Nloop = Nloop
self.Nsize = Nsize
class MonteCarloDetail:
････････
|
＜Modula-2 (Structured Programming)＞
[MAIN] MODULE randPilates.MOD
MODULE randPilates;
FROM SubPilates IMPORT PilatesStudioFast, array2D, ･･･;
VAR
Npostn : ARRAY [0.._positn-1] OF INTEGER;
xypostn: ARRAY [0.._positn-1] OF INTEGER;
postn : array2D;
Nvoid : INTEGER;
Nvposn : INTEGER;
void : ARRAY [0..10] OF vect;
vposn : ARRAY [0.._positn-1] OF vect;
･･･;
BEGIN
vScan(vposn, Nvposn, void, Nvoid, ix0,ix1, iy0,iy1);
vNeighbor(vposn, Nvposn, Nvneighbor, pos2D);
PilatesStudioFast(Nloop,Nvposn,x,y,n,Npostn,postn);
･･･;
END randPilates.
[MODULE] SubPilates
DEFINITION MODULE SubPilates;
CONST
_trials = 1000;
_positn = 42; (* array size *)
_site = 8; (* max number of neighbors *)
TYPE
array2D = ARRAY [0.._positn - 1] OF ARRAY [0.._site - 1] OF INTEGER;
vect = RECORD ix, iy: INTEGER; END;
IMPLEMENTATION MODULE SubPilates;
PROCEDURE vScan(VAR vposn: ARRAY OF vect;
VAR Nvposn: INTEGER;
VAR void: ARRAY OF vect;
Nvoid: INTEGER;
ix0,ix1: INTEGER;
iy0,iy1: INTEGER);
VAR
i,j,m: INTEGER;
p: vect;
BEGIN
m:= 0;
FOR j:=iy0 TO iy1 DO
FOR i:= ix0 TO ix1 DO
IF (i + j) MOD 2 = 1 THEN
WITH p DO
ix:= i; iy:= j;
IF NOT isCollision(p, void, Nvoid) THEN
vposn[m]:= p;
INC(m);
END;
END;
END;
END;
END;
Nvposn:= m;
END vScan;
PROCEDURE vNeighbor(VAR vposn: ARRAY OF vect;
Nvposn: INTEGER;
VAR Nneighbor: ARRAY OF INTEGER;
VAR pos2D: array2D);
VAR
k,kk,m,mm,ns: INTEGER;
s: ARRAY [0..9] OF INTEGER;
BEGIN
ns:= 0;
FOR m:= 0 TO Nvposn-1 DO
k:= 0;
FOR mm:= 0 TO Nvposn-1 DO
IF (m<>mm) AND isNear(vposn[m], vposn[mm]) THEN
s[k]:= mm+1;
INC(k);
END;
END;
Nneighbor[m]:= k;
INC(ns,k);
FOR kk:=0 TO k-1 DO
pos2D[m,kk]:= s[kk];
END;
END;
END vNeighbor;
PROCEDURE isNeighbor(ipos, jpos: INTEGER;
VAR Npostn: ARRAY OF INTEGER;
VAR postn: array2D): BOOLEAN;
VAR
p: BOOLEAN;
m: INTEGER;
BEGIN
p:= FALSE;
m:= 0;
LOOP
IF jpos+1 = postn[ipos,m] THEN
p:= TRUE;
EXIT;
END;
INC(m);
IF m = Npostn[ipos] THEN EXIT END;
END;
RETURN p;
END isNeighbor;
PROCEDURE PilatesStudioFast(Nloop: INTEGER;
N: INTEGER;
VAR x, y: ARRAY OF REAL;
VAR n: INTEGER;
VAR Npostn: ARRAY OF INTEGER;
VAR postn: array2D);
VAR
ipos,jpos: INTEGER;
loop : INTEGER;
Neighbor : INTEGER;
BEGIN
Neighbor:= 0; n:= 0;
FOR loop:= 1 TO Nloop DO
ipos:= INT(FLOAT(N)*Lib.RAND());
jpos:= INT(FLOAT(N - 1)*Lib.RAND());
IF ipos=jpos THEN jpos:= N - 1; END;
IF isNeighbor(ipos, jpos, Npostn, postn) THEN
INC(Neighbor);
END;
y[n]:= FLOAT(Neighbor)/FLOAT(loop);
x[n]:= FLOAT(loop);
INC(n);
END;
END PilatesStudioFast;
PROCEDURE PilatesStudio(Nloop: INTEGER;
VAR Neighbor: INTEGER;
npostn: INTEGER;
VAR Npostn: ARRAY OF INTEGER;
VAR xypostn: ARRAY OF INTEGER;
VAR postn: array2D);
･･････
END SubPilates.
|
