Python Script for creating STAR Particles of different Shape.
from scipy import *
import sys,os,glob
import pylab as plt
cols=['dimgray','darkgray','k','indianred','darkred','salmon','orangered','chocolate','peru','burlywood','wheat','goldenrod','khaki','y','k','k']
# 1. Generate a R-radius circle at (0,0), put as many constituent particles as possible with cutoff radius .25
# 2. Add N-number L-length "wings" by adding diagonally additional particles
# 3. Displace and Rotate Randomly in box coordinates BOXx, BOXy
# 4. Write data.star file
R=.5
r=0.22
N_stars=50
Nwings=7
L=1.0
BOXx=15
BOXy=15
# Make 1000 trials to place randomly a disk in the circle, reject particles that are too close.
from sklearn.metrics import pairwise_distances_argmin_min as pd_argmin
star_ps=[]
for i in range(10000):
xc=2.*(random.random()-0.5)*R
yc=2.*(random.random()-0.5)*R
# Find d, compare with R
d=sqrt(xc**2+yc**2)
if d <= R: #candidate for acceptance
if star_ps!=[]:
XY=array([xc,yc])
i_close, d_close = pd_argmin(XY, array(star_ps))
if d_close >= 2 * r: # don't allow for any overlap
star_ps.append([xc,yc])
else:
star_ps.append([xc,yc])
Ncore=len(star_ps)
print Ncore, '--in core'
# Add N-number L-length wings by defining angle, and Number of particles in each wing, then extend it perpendicularly one by one
theta=2*pi/Nwings
n_in_wing=int(L/(2*r) + 1)
print n_in_wing, 'n_inwing'
for i in range(Nwings):
rx = (R + r) * cos(theta*i)
ry = (R + r) * sin(theta*i)
# add n_in_wing more particles along line
star_ps.append([rx,ry])
for j in range(1,n_in_wing):
rx_add=rx + j * 2 * r *cos(theta*i)
ry_add=ry + j * 2 * r *sin(theta*i)
star_ps.append([rx_add,ry_add])
print len(star_ps) - Ncore, '-- in wings'
stars=[]
# Displace particle N times, Calculate pbcs
Nsmall_Total=0
for i in range(N_stars):
#Center randomly in box
Xc=(random.random()-0.5)* BOXx*2
Yc=(random.random()-0.5)* BOXy*2
#Displace
Xnew=array(star_ps)[:,0]+Xc
Ynew=array(star_ps)[:,1]+Yc
Bx=rint(Xnew/(2*BOXx))
By=rint(Ynew/(2*BOXy))
Xnew -= Bx * 2*BOXx
Ynew -= By * 2*BOXy
stars.append([Xnew,Ynew,Bx,By])
Nsmall_Total += len(Xnew)
fig=plt.figure()
ax=fig.add_subplot(111)
i=0
for star in stars:
ax.plot(star[0],star[1],'o',c=cols[i%12],markersize=8,alpha=0.5)
i+=1
ax.set_xlim((-BOXx,BOXx))
ax.set_ylim((-BOXy,BOXy))
print str(Nsmall_Total)+' total atoms...'
V=BOXx * BOXy
Vs= float(Nsmall_Total) * pi * r**2
phi=str(round(Vs/V,2)).ljust(4,'0')
print phi+'--packing fraction'
fig.savefig('NewConfig_N'+str(N_stars)+'_L'+str(L)+'_Nw'+str(Nwings)+'_phi'+phi+'.png')
star_ps=array(star_ps)
Nsmall_in_star=len(star_ps)
fig=plt.figure()
ax=fig.add_subplot(111)
ax.plot(star_ps[:Ncore,0],star_ps[:Ncore,1],'o',c='red',markersize=8,alpha=0.5)
ax.plot(star_ps[Ncore:,0],star_ps[Ncore:,1],'o',c='red',markersize=8,alpha=0.5)
fig.savefig('new_particle.png')
# Make data_new.star file by copying and modifying data.star
f=open('data_new_N'+str(N_stars)+'_L'+str(L)+'_Nw'+str(Nwings)+'_phi'+phi+'.star','w')
f.writelines('LAMMPS data file for Nanoparticles \n')
f.writelines(str(Nsmall_Total)+' atoms \n')
f.writelines('2 atom types \n')
f.writelines(str(-BOXx)+' '+str(BOXx)+' xlo xhi \n')
f.writelines(str(-BOXy)+' '+str(BOXy)+' ylo yhi \n')
f.writelines('-0.5 0.5 zlo zhi \n')
f.writelines(' \n')
f.writelines('Atoms \n')
f.writelines(' \n')
cnt=0
for star in stars:
Xs=star[0]
Ys=star[1]
Bx=star[2]
By=star[3]
for x,y,bx,by in zip(Xs,Ys,Bx,By):
x0=round(x,4)
y0=round(y,4)
bx0=int(bx)
by0=int(by)
cnt+=1
f.writelines(str(cnt)+' 1 1 1 '+str(x0)+' '+str(y0)+' 0 '+str(bx0)+' '+str(by0)+' 0 \n')
f.writelines(' \n')
f.writelines('Molecules \n')
f.writelines(' \n')
cnt2=0
cnt=0
for star in stars:
cnt2+=1
Xs=star[0]
Ys=star[1]
Bx=star[2]
By=star[3]
for x,y,bx,by in zip(Xs,Ys,Bx,By):
cnt+=1
f.writelines(str(cnt)+' '+str(cnt2)+' \n')
f.close()
from scipy import *
import sys,os,glob
import pylab as plt
cols=['dimgray','darkgray','k','indianred','darkred','salmon','orangered','chocolate','peru','burlywood','wheat','goldenrod','khaki','y','k','k']
# 1. Generate a R-radius circle at (0,0), put as many constituent particles as possible with cutoff radius .25
# 2. Add N-number L-length "wings" by adding diagonally additional particles
# 3. Displace and Rotate Randomly in box coordinates BOXx, BOXy
# 4. Write data.star file
R=.5
r=0.22
N_stars=50
Nwings=7
L=1.0
BOXx=15
BOXy=15
# Make 1000 trials to place randomly a disk in the circle, reject particles that are too close.
from sklearn.metrics import pairwise_distances_argmin_min as pd_argmin
star_ps=[]
for i in range(10000):
xc=2.*(random.random()-0.5)*R
yc=2.*(random.random()-0.5)*R
# Find d, compare with R
d=sqrt(xc**2+yc**2)
if d <= R: #candidate for acceptance
if star_ps!=[]:
XY=array([xc,yc])
i_close, d_close = pd_argmin(XY, array(star_ps))
if d_close >= 2 * r: # don't allow for any overlap
star_ps.append([xc,yc])
else:
star_ps.append([xc,yc])
Ncore=len(star_ps)
print Ncore, '--in core'
# Add N-number L-length wings by defining angle, and Number of particles in each wing, then extend it perpendicularly one by one
theta=2*pi/Nwings
n_in_wing=int(L/(2*r) + 1)
print n_in_wing, 'n_inwing'
for i in range(Nwings):
rx = (R + r) * cos(theta*i)
ry = (R + r) * sin(theta*i)
# add n_in_wing more particles along line
star_ps.append([rx,ry])
for j in range(1,n_in_wing):
rx_add=rx + j * 2 * r *cos(theta*i)
ry_add=ry + j * 2 * r *sin(theta*i)
star_ps.append([rx_add,ry_add])
print len(star_ps) - Ncore, '-- in wings'
stars=[]
# Displace particle N times, Calculate pbcs
Nsmall_Total=0
for i in range(N_stars):
#Center randomly in box
Xc=(random.random()-0.5)* BOXx*2
Yc=(random.random()-0.5)* BOXy*2
#Displace
Xnew=array(star_ps)[:,0]+Xc
Ynew=array(star_ps)[:,1]+Yc
Bx=rint(Xnew/(2*BOXx))
By=rint(Ynew/(2*BOXy))
Xnew -= Bx * 2*BOXx
Ynew -= By * 2*BOXy
stars.append([Xnew,Ynew,Bx,By])
Nsmall_Total += len(Xnew)
fig=plt.figure()
ax=fig.add_subplot(111)
i=0
for star in stars:
ax.plot(star[0],star[1],'o',c=cols[i%12],markersize=8,alpha=0.5)
i+=1
ax.set_xlim((-BOXx,BOXx))
ax.set_ylim((-BOXy,BOXy))
print str(Nsmall_Total)+' total atoms...'
V=BOXx * BOXy
Vs= float(Nsmall_Total) * pi * r**2
phi=str(round(Vs/V,2)).ljust(4,'0')
print phi+'--packing fraction'
fig.savefig('NewConfig_N'+str(N_stars)+'_L'+str(L)+'_Nw'+str(Nwings)+'_phi'+phi+'.png')
star_ps=array(star_ps)
Nsmall_in_star=len(star_ps)
fig=plt.figure()
ax=fig.add_subplot(111)
ax.plot(star_ps[:Ncore,0],star_ps[:Ncore,1],'o',c='red',markersize=8,alpha=0.5)
ax.plot(star_ps[Ncore:,0],star_ps[Ncore:,1],'o',c='red',markersize=8,alpha=0.5)
fig.savefig('new_particle.png')
# Make data_new.star file by copying and modifying data.star
f=open('data_new_N'+str(N_stars)+'_L'+str(L)+'_Nw'+str(Nwings)+'_phi'+phi+'.star','w')
f.writelines('LAMMPS data file for Nanoparticles \n')
f.writelines(str(Nsmall_Total)+' atoms \n')
f.writelines('2 atom types \n')
f.writelines(str(-BOXx)+' '+str(BOXx)+' xlo xhi \n')
f.writelines(str(-BOXy)+' '+str(BOXy)+' ylo yhi \n')
f.writelines('-0.5 0.5 zlo zhi \n')
f.writelines(' \n')
f.writelines('Atoms \n')
f.writelines(' \n')
cnt=0
for star in stars:
Xs=star[0]
Ys=star[1]
Bx=star[2]
By=star[3]
for x,y,bx,by in zip(Xs,Ys,Bx,By):
x0=round(x,4)
y0=round(y,4)
bx0=int(bx)
by0=int(by)
cnt+=1
f.writelines(str(cnt)+' 1 1 1 '+str(x0)+' '+str(y0)+' 0 '+str(bx0)+' '+str(by0)+' 0 \n')
f.writelines(' \n')
f.writelines('Molecules \n')
f.writelines(' \n')
cnt2=0
cnt=0
for star in stars:
cnt2+=1
Xs=star[0]
Ys=star[1]
Bx=star[2]
By=star[3]
for x,y,bx,by in zip(Xs,Ys,Bx,By):
cnt+=1
f.writelines(str(cnt)+' '+str(cnt2)+' \n')
f.close()
No comments:
Post a Comment