[ESPResSo-users] WARNING: pressure calculated, but MMM2D pressure not implemented

[padidela uday]

Dear all,

I am trying to implement charged wall constraints in one of my scripts along with the particles between the constraints. But however I land up in few troubles with it.

After equilibrating the system a warning is displayed on screen during sampling

calculating pressure for electrostatics method that does have it implemented

WARNING: pressure calculated, but MMM2D pressure not implemented

Is there something wrong with the code could you help me to handle this?

Thank you for your time

set n_solute 200
set n_solvent 200
set top_wall 100
set bottom_wall 100
set n_wall 200
set n_part [expr ($n_solute + $n_solvent + $top_wall + $bottom_wall)]
set box_l 200
set sig1 2.0
set box_lxy [expr sqrt(2) * $box_l]
set box_lz [expr 0.5*$box_l]

setmd box_l $box_lxy $box_lxy [expr $box_lz +  $sig1]
setmd periodic 1 1 0 

for {set i 0} { $i < 200 } {incr i} {
    set posx [expr $box_lxy*[t_random]]
    set posy [expr $box_lxy*[t_random]]
    set posz [expr $box_lz*[t_random] + $sig1]
    set vx [gauss_random]
    set vy [gauss_random]
    set vz [gauss_random]
    part $i pos $posx $posy $posz type 0 v $vx $vy $vz q -2 mass 1.0
}

for {set i 200} { $i < 400 } {incr i} {
    set posx [expr $box_lxy*[t_random]]
    set posy [expr $box_lxy*[t_random]]
    set posz [expr $box_lz*[t_random] + $sig1]
    set vx [gauss_random]
    set vy [gauss_random]
    set vz [gauss_random]
    part $i pos $posx $posy $posz type 1 v $vx $vy $vz q 1 mass 1.0
}

#wall at top
constraint wall normal 0 0  1 dist $sig1 type 2

for {set i 400} { $i < 500 } {incr i} {
    set posx [expr $box_lxy*[t_random]]
    set posy [expr $box_lxy*[t_random]]
    set posz [expr ($box_lz + $sig1)]

    part  [expr $i + $one_wall] pos $posx $posy $posz type 2  q 1
}

constraint wall normal 0 0 -1 dist [expr -$box_lz - $box_lz] type 2

#wall at bottom

for {set i 500} { $i < 600 } {incr i} {
    set posx [expr $box_lxy*[t_random]]
    set posy [expr $box_lxy*[t_random]]
    set posz [expr ($box_lz - $box_lz)*[t_random]]
    part  [expr $i + $two_wall] pos $posx $posy $posz type 2  q 1
}

set sigma [expr -0.25*$n_wall/($box_lxy*$box_lxy)]
constraint plate height 0 sigma $sigma
constraint plate height [expr $box_lz + $sig1] sigma $sigma

setmd time_step 0.00002624
setmd skin 0.3
set temp 4.9576
set gamma 1
thermostat langevin $temp $gamma

#  POTENTIAL
set sig1 2.0
set sig2 1.15
set sig3 0.3
set eps1 1.0
set eps2 0.7
set eps3 0.5
set cut1 [expr 1.12246*$sig1]
set cut2 [expr 1.12246*$sig2]
set cut3 [expr 1.12246*$sig3]
set shift1 [expr 0.25*$eps3]
set shift2 [expr 0.25*$eps3]
set shift3 [expr 0.25*$eps3]
inter 0 0 lennard-jones $eps3 $sig1 $cut1 $shift1  0 
inter 1 0 lennard-jones $eps3 $sig2 $cut2 $shift2  0 
inter 1 1 lennard-jones $eps3 $sig3 $cut3 $shift3  0 
inter 2 0 lennard-jones $eps3 $sig1 $cut1 $shift1  0
inter 2 1 lennard-jones $eps3 $sig2 $cut2 $shift2  0
#inter 2 2 lennard-jones $eps3 $sig3 $cut3 $shift3  0 

cellsystem layered 3
inter coulomb 1.0 mmm2d 1e-4

if { [regexp "ROTATION" [code_info]] } {
      set deg_free 6
} {
      set deg_free 3
}

# prepare vmd connection
if { $vmd_output=="yes" } {
   prepare_vmd_connection "vmd"
    imd listen 10000
}

inter ljforcecap individual
set inter_steps 10000
for {set cap 10} {$cap < 100} {incr  10} {
       set rad [expr 1.0 - 0.5*$cap/10.0]
       set lb [expr 1.0 * $cap / 10.0]
       inter 0 0 lennard-jones $eps3 $sig1 $cut1 $shift1  0  $rad
       inter 1 0 lennard-jones $eps3 $sig2 $cut2 $shift2  0  $rad
       inter 1 1 lennard-jones $eps3 $sig3 $cut3 $shift3  0  $rad
       inter coloumb $lb mmm2d 1e-4
       set temp [expr [analyze energy kinetic]/(($deg_free/2.0)*$n_part)]
    puts "t=[setmd time] E=[analyze energy total] T=$temp"
         if {$vmd_output=="yes"} {imd positions}
    inter ljforcecap $cap;
    integrate $inter_steps
}
inter ljforcecap 0
inter coulomb 1.0 mmm2d 1e-4

set h [open "production.data" "w"]
set g [open "analysis.data" "w"]
set n_part [expr ($n_solute + $n_solvent + $top_wall + $bottom_wall)]
for {set i 0} { $i <4000 } {incr i} {
    set temp [expr [analyze energy kinetic]/(($deg_free/2.0)*$n_part)]
        puts $h " $i [setmd time] [analyze energy kinetic] $temp"
        integrate 200

    set f [open "config_$i" "w"]
    blockfile $f write tclvariable {box_lxy box_lz}
    blockfile $f write variable box_l
    blockfile $f write particles {id pos type v f}
    set temp [expr [analyze energy kinetic]/(($deg_free/2.0)*$n_part)]
    puts $f " \ { energy [analyze energy total] $temp}"
    puts $g " [analyze energy kinetic] [analyze pressure total] [analyze energy total] $temp"

    close $f
}
close $g

With Kind Regards

Uday