Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F77993921
test_mpi.py
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sat, Aug 17, 17:45
Size
3 KB
Mime Type
text/x-python
Expires
Mon, Aug 19, 17:45 (2 d)
Engine
blob
Format
Raw Data
Handle
19963337
Attached To
rGEAR Gear
test_mpi.py
View Options
#!/usr/bin/env python
from
pNbody
import
*
from
pNbody
import
ic
from
pNbody
import
libutil
from
numpy
import
*
import
sys
import
time
import
copy
from
PyGear
import
gadget
nb
=
Nbody
(
"snap.dat"
,
ftype
=
"gadget"
)
###################################
# init PyGear
gadget
.
InitMPI
()
# init MPI
gadget
.
InitDefaultParameters
()
# init default parameters
gadget
.
Info
()
params
=
{}
params
[
'ErrTolTheta'
]
=
0.7
params
[
'DesNumNgb'
]
=
15
params
[
'MaxNumNgbDeviation'
]
=
3
params
[
'UnitLength_in_cm'
]
=
3.085e+21
params
[
'UnitMass_in_g'
]
=
4.435693e+44
params
[
'UnitVelocity_in_cm_per_s'
]
=
97824708.2699
params
[
'BoxSize'
]
=
1.0
params
[
'PartAllocFactor'
]
=
10
gadget
.
SetParameters
(
params
)
###################################
# load particles
gadget
.
LoadParticles
(
array
(
nb
.
npart
),
nb
.
pos
,
nb
.
vel
,
nb
.
mass
,
nb
.
num
,
nb
.
tpe
)
###################################
# compute Delaunay
t1
=
time
.
time
()
gadget
.
ConstructDelaunay
()
t2
=
time
.
time
()
print
"Delaunay"
,
"in"
,(
t2
-
t1
),
"s"
###################################
# compute Voronoi
t1
=
time
.
time
()
gadget
.
ComputeVoronoi
()
t2
=
time
.
time
()
print
"Voronoi"
,
"in"
,(
t2
-
t1
),
"s"
#sys.exit()
###################################
# retrieve ghost points
gpos
=
gadget
.
GetAllGhostPositions
()
gnb
=
Nbody
(
pos
=
gpos
,
ftype
=
'gadget'
)
gnb
.
u
=
zeros
(
gnb
.
nbody
)
gnb
.
rho
=
gadget
.
GetAllGhostvDensities
()
gnb
.
rsp
=
gadget
.
GetAllGhostvVolumes
()
gnb
.
init
()
# save
gnb
.
rename
(
'ghost.dat'
)
gnb
.
set_pio
(
"yes"
)
gnb
.
write
()
###################################
# retrieve points
nb
.
pos
=
gadget
.
GetAllPositions
()
nb
.
vel
=
gadget
.
GetAllVelocities
()
nb
.
mass
=
gadget
.
GetAllMasses
()
nb
.
num
=
gadget
.
GetAllID
()
nb
.
tpe
=
None
nb
.
u
=
gadget
.
GetAllEnergySpec
()
# not defined
nb
.
rho
=
gadget
.
GetAllvDensities
()
# note the v
nb
.
rsp
=
gadget
.
GetAllvVolumes
()
nb
.
rsp
=
gadget
.
GetAllDensities
()
# put sph densities in rsp
nb
.
init
()
# save
nb
.
rename
(
'snap2.dat'
)
nb
.
set_pio
(
"yes"
)
nb
.
write
()
# save all
#nbtot = nb + gnb
nbtot
=
copy
.
deepcopy
(
nb
)
# this is important to avoid to mix particles
nbtot
.
append
(
gnb
,
do_not_sort
=
True
)
# then they are incompatible with vpos below
nbtot
.
rename
(
'snap+ghost.dat'
)
nbtot
.
write
()
###################################
# get some values
TriangleList
=
gadget
.
GetAllDelaunayTriangles
()
rho
,
mn
,
mx
,
cd
=
libutil
.
set_ranges
(
nbtot
.
rho
,
scale
=
'log'
)
rho
=
rho
/
255.
#rho = clip(nbtot.rho,0,1)
###################################
# plot
###################################
#sys.exit()
import
Ptools
as
pt
import
plot
########################
# draw a box
plot
.
draw_box
(
x
=
array
([
0
,
1
,
1
,
0
]),
y
=
array
([
0
,
0
,
1
,
1
]))
########################
# draw the triangles
i
=
0
for
Triangle
in
TriangleList
:
P1
=
Triangle
[
'coord'
][
0
]
P2
=
Triangle
[
'coord'
][
1
]
P3
=
Triangle
[
'coord'
][
2
]
#plot.draw_triangle(P1,P2,P3,c='r')
#cm = 1/3.*(P1+P2+P3)
#pt.text(cm[0],cm[1],Triangle['id'],fontsize=12,horizontalalignment='center',verticalalignment='center')
########################
# draw the points
#pt.scatter( nb.pos[:,0], nb.pos[:,1],lw=0,s=50,color='r')
#pt.scatter(gnb.pos[:,0],gnb.pos[:,1],lw=0,s=50)
########################
# draw the Voronoi
cmap
=
pt
.
GetColormap
(
'rainbow4'
,
revesed
=
False
)
# single point with vpoints
for
i
in
xrange
(
nbtot
.
nbody
):
vpos
=
gadget
.
GetvPointsForOnePoint
(
i
)
#pt.scatter(nb.pos[i,0],nb.pos[i,1],lw=0,s=50)
plot
.
draw_cell
(
vpos
,
color
=
[
rho
[
i
],
rho
[
i
],
rho
[
i
]],
alpha
=
0.8
)
pt
.
axis
([
-
0.1
,
1.1
,
-
0.1
,
1.1
])
pt
.
show
()
Event Timeline
Log In to Comment