[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
RE: [ESPResSo] Question about bond-angle interaction
From: |
Limbach, Hans Joerg, LAUSANNE, NRC-FS |
Subject: |
RE: [ESPResSo] Question about bond-angle interaction |
Date: |
Thu, 18 Oct 2007 08:50:32 +0200 |
Hello Lorenzo,
With two particles sticking yoour problem is tough to solve, since all
the interactions between two particles are only distance dependend
(included the bonded ones). There are also non spherical potentials in
Espresso like the Gay-Burne potential, but this would make it too
complicated in my opinion.
My suggestion would be the following:
Case 1: two particles collide:
FENE interaction:
inter 0 fene 7 2 (<k> <max dist>)
Identities: p1, p2
part $p1 bond 0 $p2 (Syntax: part <part identity 1> bond <bond type
number> <part identity 2>)
Case 2: three particles collide:
Here depending on the bond angle alpha between the three particles p1,
p2 and p3 (p2 being the particle in the middle) you have to define a new
bonded interaction:
set i 1
inter $i angle 1000 $alpha
incr i
part $p2 bond $i $p1 $p3
(But this would mean you end up with a huge number of potentials (code
will get slow) so it is probably better to define a finite number of
bond angle potentials e.g. one for every 10 degree step)
Case 3: four particles collide:
The same as in case 2, but with a dihedral potential.
This will freeze in the 3D structure in which the particles collide with
each other.
Hope this works,
Hanjo
> -----Original Message-----
> From: address@hidden
> [mailto:address@hidden On Behalf Of
> Lorenzo Isella
> Sent: jeudi, 18. octobre 2007 01:38
> To: address@hidden
> Subject: [ESPResSo] Question about bond-angle interaction
>
> 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
> # follow what your simulation is doing
> 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
>
>
>
>
> _______________________________________________
> ESPResSo mailing list
> address@hidden
> https://fias.uni-frankfurt.de/mailman/listinfo/espresso
>