# Created by charmm2lammps v1.9.1 on Mon 13 Dec 2021 12:44:15 AM +0330
# Command: charmm2lammps.pl all27_prot_lipid water1_autopsf

echo            both
units           real
boundary        p p p
neigh_modify    delay 2 every 1
timestep        1
atom_style      full

bond_style      harmonic
angle_style     harmonic

read_data       config.data

bond_style harmonic	
bond_coeff 2 469.0  1.0     #H-O
bond_coeff 1 450.0  1.418      #C-C

angle_style harmonic
angle_coeff 1  85.0  109.47    #H-O-H
#angle_coeff 2  55.0  109.47    #C-C-C


pair_style       lj/charmm/coul/long 8 10 

pair_coeff     1  1     0.0     0.0     #H-H
pair_coeff     2  2    0.1553    3.166     #O-O
pair_coeff     3  3     0.07    3.9848     #C-C

pair_coeff     1  2    0.0     0.0     #H-O
pair_coeff     1  3      0.0     1.9924     #H-C
pair_coeff     2  3      0.104264     3.5754     #O-C

kspace_style    pppm 1e-4

neighbor	      2 bin
neigh_modify     exclude type 3 3 
neigh_modify     delay 0 every 1 check yes

group water_molec type 1 2

velocity water_molec create 300.0 12345678 dist uniform 

fix NVT water_molec nvt temp 300.0 300.0 100   

fix rigid water_molec shake 1.0e-6 500 0 t 1 2 a 1 b 2

group water_box_1 id  <> 1 11763
group piston id <> 11764  12771 
group cnt id <> 13780  14419
group water_box_2 id <> 16436 28198
group sheet2 id <> 12772 13779
group sheet34 id <> 14420 16435

group g_all union cnt sheet2 sheet34 

velocity piston   set 0 0 0 units box
velocity g_all     set 0 0 0 units box



fix NOFORCE1 piston setforce 0 0 0
fix NOFORCE2 g_all setforce 0 0 0

region   cyl_1 cylinder z 25 25 6.77 59 101 units box
region   cyl_2 cylinder z 25 25 6.77 159 161 units box

group cyl_1 region cyl_1
group cyl_2 region cyl_2

delete_atoms group cyl_1
delete_atoms group cyl_2

fix push piston nve
fix force piston aveforce 0 0 0.364593


#region box1 block EDGE EDGE EDGE EDGE 0 60 units box
#region box2 block EDGE EDGE EDGE EDGE 100 160 units box
region  cyl_out cylinder z 25 25 6.77 61 99 units box

group O type 2

################# Pressure Tensor calculation ###############

variable pxy equal pxy
variable pxz equal pxz
variable pyz equal pyz
variable pxx equal pxx  
variable pyy equal pyy  
variable pzz equal pzz 

group mine dynamic O region cyl_out every 5

fix P_tensor1 mine ave/time 1 10 30 v_pxx v_pyy v_pzz v_pxy v_pxz v_pyz file ptensor1.profile #ptensor1

fix P_tensor2 water_box_1 ave/time 1 10 30 v_pxx v_pyy v_pzz v_pxy v_pxz v_pyz file ptensor2.profile #ptensor2

thermo 100
thermo_style custom step temp press v_pxx v_pyy v_pzz v_pxy v_pxz v_pyz
################# End Pressure Tensor calculation ###############



dump 2 all custom 1000 dump2.lammpstrj id type x y z
run 1000000 #1ns 1000000