espressomd-users
[Top][All Lists]

## [ESPResSo] Fwd: analyze nbhood & Bond-angle Interactions

 From: Lorenzo Isella Subject: [ESPResSo] Fwd: analyze nbhood & Bond-angle Interactions Date: Fri, 19 Oct 2007 17:37:43 +0200

```Dear All,
Thanks for the very useful bits of advice I got on the list.
I am now running some experiments with really simple systems, made up
of 3 particles.
The idea is how suggested on this list: when particles get "close
enough", I stick them together with a potential which depends on the
angle (a rigid bond).
I have a few questions for the list.
(1) if I look at Eq. 5.16 of the Espresso manual (bond_angle_cosine) I
do not see how the length of the bond is defined. I mean: the
particles will have a preferred orientation but what is the range of
the interaction? Can I specify it?
(2) the command analyze nbhood  tells me how many particles I have
within a given distance from a certain particle. It returns a list of
particle numbers, but what is the meaning of the order in which the
numbers appear? I mean, say I am looking for the particles surrounding
particle 1 and I get the array 0 1 2. Does that tell me anything about
which particle is in the middle?
Particles 0 and 2 could be on the same side wrt particle 1, or
particle 1could be somewhere in between. This is important seen that I
need to know which is the central particle to set up the bond-angle
interaction.
(3)In the end it all boils down to the problems I am experiencing in
setting up a system with a few particles able to interact via a
bond-angle interaction.
In the code below, everything should work (in the end I check that my
particle 1 has the bonded interaction), but from VMD that does not
seem to be the case...
Cheers

Lorenzo

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

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

inter 0 angle 700.0 0.8
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
}

puts "[part ]"

setmd time_step 0.01; setmd skin 0.4
set temp 0.0005; set gamma 0.01
thermostat langevin \$temp \$gamma
set integ_steps 200
set vmd "yes"

if { \$vmd == "yes" } {
prepare_vmd_connection tutorial 3000
exec sleep 4
imd positions
}

set activate_fene 1
set fene_warning 1
set integ_step 20

for {set i 0} { \$i < 2000 } { 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

if { \$vmd == "yes" } { imd positions }
set acti [analyze nbhood 1 4.5]

set neighbor [analyze nbhood 1 4.]
#puts " the neighbors are, \$neighbor "
set mylen [llength \$neighbor]

if { \$mylen == 3 && \$activate_fene == 1 } {

part 1 bond 0 0 2
puts "activated angle bond!"
set activate_fene 0
}
}

puts "[part ]"

puts "So far so good"

```