-
Notifications
You must be signed in to change notification settings - Fork 0
/
randomWalk.py
190 lines (164 loc) · 7.08 KB
/
randomWalk.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
# Python code for 2D random walk.
# from cmath import sqrt
from cmath import sqrt
# from logging import _Level
import logging
import math
import numpy as np
import random
import pandas
import plotly.graph_objects as go
from config import *
from numba import jit
from octTreeCPU import *
from util import particlesToDF
# -----------------------------------------functions for walk algorithms: --------------------------------------------------------------------------------------------------------------------------------------------------------------------
def getInitialSphere(
particlesNumber=particlesNumber,
porosityFraction=porosityFraction,
sphereRadius=sphereRadius,
):
return particlesToDF(getInitialSphereNumpy())
def getInitialSphereNumpy(
particlesNumber=particlesNumber,
porosityFraction=porosityFraction,
sphereRadius=sphereRadius,
):
# child variables
vacancies = round(particlesNumber * porosityFraction)
totalPositions = particlesNumber + vacancies
# creating full array
x = 0
y = 0
z = 0
initialSphere = []
# From the panda library, a Data frame is just a c++ maps that are in the form of arrays
# This while loop is getting all the x, y, z integer values insde of 3D sphere (in all 8 3D quadrants). Each of these combination of x,y,z values,
# of which there are ~840, will be the vector components of a vector whose magnitude is an integer, conforming to a radius of a smaller sphere than the original one
# set in the limit
# initialSphere.index returns total rows in the initialSphere dataFrame
for x in range(-sphereRadius, sphereRadius):
for y in range(-sphereRadius, sphereRadius):
for z in range(-sphereRadius, sphereRadius):
if math.sqrt(x**2 + y**2 + z**2) <= sphereRadius:
initialSphere.append([np.int32(x), np.int32(y), np.int32(z)])
totalPositions = len(initialSphere)
vacancies = round(totalPositions * porosityFraction)
# print("Initial sphere complete")
# this randomizes the areas where the particles can and cannot go in the sphere
for i in range(0, vacancies):
val = random.randint(0, len(initialSphere) - 1)
initialSphere.pop(val)
return np.array(initialSphere, dtype=np.int32)
def randomWalkCPU(
initialSphere,
capillaryRadius,
n=n,
maxTries=maxTries,
sphereRadius=sphereRadius,
):
# Constraints for cell movement
squaredRadius = sphereRadius**2
squaredCapillaryRadius = capillaryRadius**2
# creating two array for containing x and y coordinate
# of size equals to the number of size and filled up with 0's
particles = [
initialSphere
] # particles is now a list containing the first timestep result
# random walking for n timesteps and add to each timestep result to particles
for i in range(1, n + 1):
particles.append(particles[-1].copy())
particleN = -1
for particle in particles[-1]:
particleN += 1
for j in range(maxTries):
walkedParticle = particle.copy()
val = random.randint(1, 6)
if val == 1:
walkedParticle[0] += 1
elif val == 2:
walkedParticle[0] -= 1
elif val == 3:
walkedParticle[1] += 1
elif val == 4:
walkedParticle[1] -= 1
elif val == 5:
walkedParticle[2] += 1
else:
walkedParticle[2] -= 1
x_2 = walkedParticle[0] ** 2
y_2 = walkedParticle[1] ** 2
z_2 = walkedParticle[2] ** 2
# comparing this values to the previous dataFrame (particles[i-1]) means we don't want the particle to move to a past position, nor do we want it to move to the x,y,z coordinate of a current position, last part is we want to squared distance to be within squared capillary radius
if (
not any(
np.equal(
particles[i],
[
walkedParticle[0],
walkedParticle[1],
walkedParticle[2],
],
).all(1)
)
) and ( # No particle in this timestep shares the same position
(x_2 + y_2 + z_2) < squaredRadius
or (x_2 + z_2) < squaredCapillaryRadius
or (y_2 + z_2) < squaredCapillaryRadius
): # If inside the capillary radius
particles[i][particleN] = walkedParticle
break
(i % 100 == 0) and print("Time steps elapsed: " + str(i))
return particles
def randomWalkCPUOctTree(
initialSphere,
capillaryRadius,
n=n,
maxTries=maxTries,
sphereRadius=sphereRadius,
):
# Constraints for cell movement
squaredRadius = sphereRadius**2
squaredCapillaryRadius = capillaryRadius**2
# creating two array for containing x and y coordinate
# of size equals to the number of size and filled up with 0's
particles = [
initialSphere
] # particles is now a list containing the first timestep result
# random walking for n timesteps and add to each timestep result to particles
for i in range(1, n + 1):
boundRange = (sphereRadius + 1 + i) * 2
tree = buildTreeCPU(boundRange=boundRange)
particles.append(particles[i - 1].copy())
# Now walk the particles (the insert function returns if successful or not)
particleN = -1
for particle in particles[i]:
particleN += 1
for j in range(maxTries):
walkedParticle = particle.copy()
val = random.randint(1, 6)
if val == 1:
walkedParticle[0] += 1
elif val == 2:
walkedParticle[0] -= 1
elif val == 3:
walkedParticle[1] += 1
elif val == 4:
walkedParticle[1] -= 1
elif val == 5:
walkedParticle[2] += 1
else:
walkedParticle[2] -= 1
x_2 = walkedParticle[0] ** 2
y_2 = walkedParticle[1] ** 2
z_2 = walkedParticle[2] ** 2
# If able to insert into tree (unique value for this timestep) and inside the capillary
if tree.insert(walkedParticle) and (
(x_2 + y_2 + z_2) < squaredRadius
or (x_2 + z_2) < squaredCapillaryRadius
or (y_2 + z_2) < squaredCapillaryRadius
):
particles[i][particleN] = walkedParticle
break
(i % 100 == 0) and print("Time steps elapsed: " + str(i))
return particles