espressomd-users
[Top][All Lists]

 From: Lorenzo Isella Subject: [ESPResSo] Question about Date: Thu, 18 Oct 2007 01:37:16 +0200 User-agent: Icedove 1.5.0.12 (X11/20070730)

```Dear All,
```
I have already posted on this list, mainly to ask whether I could used Espresso for simulating the dynamics of nanoparticles rather than molecules/polimers. My problem is to have these particles "sticking" when they collide. This is essential for me, since I want to study the complicated, absolutely non-spherical, aggregates they give rise to. I have been toying around both with Lennard-Jones and Morse potential. The problem (which I think is due to the potentials being a function of the distance between two particles only) is that an aggregate is indeed formed, but I always end up with a compact, almost spherical structure.
```To make it less an abstract discussion, here is the code I have used:

set n_part 20; set density 0.07
set box_l [expr pow(\$n_part/\$density,1./3.)]

setmd box_l \$box_l \$box_l \$box_l
setmd periodic 1 1 1

set q 0; set type 0

```
# I now set the particle charge to 0 and see if the code works the same of not
```
for {set i 0} { \$i < \$n_part } {incr i} {
set posx [expr \$box_l*[t_random]]
set posy [expr \$box_l*[t_random]]
set posz [expr \$box_l*[t_random]]
part \$i pos \$posx \$posy \$posz q \$q type \$type
}

setmd time_step 0.01; setmd skin 0.4
set temp 0.0005; set gamma 0.01
thermostat langevin \$temp \$gamma

#set sig 1.0; set cut [expr 1.12246*\$sig]
#set eps 1.0; set shift [expr 0.25*\$eps]
#inter 0 0 lennard-jones \$eps \$sig \$cut \$shift 0
inter 0 0  morse 1000. 0.01 1. 1.3

#inter 3 angle 1000.0 3.14

set integ_steps 200

for {set cap 20} {\$cap < 200} {incr cap 20} {
puts "t=[setmd time] E=[analyze energy total]"
inter ljforcecap \$cap; integrate \$integ_steps
set min [analyze mindist]

}

puts "Warmup finished. Minimal distance now \$min"
#uncap forces
inter ljforcecap 0

puts "end of equilibration"

set vmd "yes"

if { \$vmd == "yes" } {
# This calls a small tcl script which starts the program    #
# VMD and opens a socket connection between ESPResSo and    #
# VMD.                                                      #
prepare_vmd_connection tutorial 3000

# Just wait a moment until VMD has started.                 #
# The 'exec' command is quite useful since with that you can#
# call any other program from within your simulation script.#
exec sleep 4

# The additional command imd steers the socket connection   #
# to VMD, e.g. sending the actual coordinates               #
imd positions
}

#prepare the saving of the results

#set obs [open "rg_lorenzo.dat" "w"]

set integ_step 20

for {set i 0} { \$i < 1000 } { incr i} {
set temp [expr [analyze energy kinetic]/(1.5*\$n_part)]
puts "t=[setmd time] E=[analyze energy total], T=\$temp"
integrate \$integ_steps

puts "save configuration"

set f [open "config_\$i" "w"]
blockfile \$f write tclvariable {box_l density}
blockfile \$f write variable box_l
blockfile \$f write particles {id pos type}
close \$f
# if you have turned on the vmd option you can now
if { \$vmd == "yes" } { imd positions }

#Now I calculate the radius of gyration
# set rg [lindex [analyze rg] 0]
#puts \$obs "[setmd time] \$rg"

}

#close \$obs
puts "end of integration"

#plotObs "rg.dat" { 1:2 } labels { "time" "rg" } out "rg"
#exec gv rg.ps

set f [open "config_1" "r"]
while { [blockfile \$f read auto] != "eof" } {}
close \$f

puts "ok reading the block file"

set rdf [analyze rdf 0 0 0.9 [expr \$box_l/2] 100]
set rlist ""
set rdflist ""
foreach value [lindex \$rdf 1] {
lappend rlist [lindex \$value 0]
lappend rdflist [lindex \$value 1]
}

```
Any suggestions about a different interaction to use? In a previous post I was suggested a combination of Lennard-Jones and a stochastically added FENE bond, but this is still obscure to me. So, I was thinking about some interaction which depends on the angle determined by 3 particles (bond-angle interaction). However, setting the interaction in the previous code as:
```
inter 4 angle 1000.0 3.14
part 0  bond 4 0 0

```
Leads to some strange errors. Actually, the number after part should probably be the individual number stuck on each particles and not the type but then how can I say that I want all the particles to interact according to that potential? And then: "angle" stands for BOND_ANGLE_COSINE according to my myconfig.h file.
```
```
I am a bit confused here, I may be using funny numbers in the code. Any help is really appreciated.
```Many thanks

```
Lorenzo
```

```