Dear Dr.Chaban, I tried to get the same results as you got for density and total energy starting from the 8 chains is the box. I used the revised mdp file and since PME option takes much time I used shift instead for a moment to verify if I can get the density you obtained. Also I used constraints for bonds. density I get after 1800ps is 850 kg/m3 while what you got is a slightly different (after 1000ps density varies around 880!). Also the total energy you have is bout 2300 whereas I get somethong around 12000!. I dont know why they are so different..
In the second trial I used "constraints = none" just the "shift" for electrostatics is not the same as the mdp file you suggested. again results are not the same for total energy (20000 Kj/mol). Below is the average values...also simulation crashes at 1550 ps and the last density I am getting is 834! still less than what you got. dont really know why results are not the same. The reason simulation crashes could be due to dt=0.002 while removing constraints? I have two more inquiries: 1- What I see is that the system is getting compressed to a size of around 3 3 3 nm is the first 20ps abruptly and in the remainder of time system reaches equilibrium... This is very different approach than what I was trying doing before. I used berendsen scheme for pressure coupling and box reduced is size slowly but I had difficulty reaching the size I want and most of the simulations crashed. Does this mean berenden is not suitable for compressing system for all systems? 2- You said I dont need NVT at all. But if I want to study the system I get from NPT giving the density I want, can I not do NVT for a fixed system size and calcualte energies or other quantities I am after? Thanks for your attention Moeed Energy Average RMSD Fluct. Drift Tot-Drift ------------------------------------------------------------------------------- Bond 3885.23 141.365 139.725 0.0478427 74.3668 Angle 5864.23 149.886 149.885 -0.00114897 -1.78597 Ryckaert-Bell. 2243.62 144.068 134.737 -0.113664 -176.679 LJ-14 1167.5 30.8742 30.6402 -0.00845667 -13.1451 Coulomb-14 -553.954 31.3413 30.7865 -0.0130831 -20.3364 LJ (SR) -4898.73 405.434 386.584 -0.272304 -423.27 Coulomb (SR) 1250.22 56.2361 56.1252 0.0078672 12.2288 Potential 8958.11 379.038 344.366 -0.352947 -548.621 Kinetic En. 10834.2 156.085 156.066 -0.00535835 -8.32903 Total Energy 19792.3 414.016 381.523 -0.358305 -556.95 Temperature 300.069 4.32299 4.32248 -0.000148407 -0.230684 Pressure (bar) 29.5896 1263.17 1263.16 0.0112876 17.5454 Box-X 3.67408 1.72247 1.70558 -0.000536305 -0.833634 Box-Y 3.44445 1.61482 1.59898 -0.000502786 -0.781532 Box-Z 2.45126 1.66735 1.64995 -0.000535423 -0.832263 Density (SI) 815.504 72.8491 70.7123 0.0390324 60.6721 pV 54.9014 2111.1 2111.07 -0.0208174 -32.3586 Vir-XX 3606.98 1466.95 1466.92 0.0213665 33.2121 Box-Vel-XX -0.0182685 0.234428 0.232397 6.8633e-05 0.106683 Mu-Y 0.102393 0.966426 0.966077 -5.79406e-05 -0.0900631 Heat Capacity Cv: 12.4757 J/mol K (factor = 0.000207551) ; Run control integrator = md ; type of dynamics algorithm. Here md uses a leap-frog algorithm for integrating Newtons's eq of motion dt = 0.002 ; in ps ! nsteps = 900000 ; length of simulation= nsteps*dt nstcomm = 1 ; frequency for center of mass motion removal ; Output control nstenergy = 100 ; frequency to write energies to energy file. i.e., energies and other statistical data are stored every 10 steps nstxout = 100 ; frequency to write coordinates/velocity/force to output trajectory file. how often snapshots are collected= nstxout*dt nstvout = 0 nstfout = 0 nstlog = 1000 ; frequency to write energies to log file nstxtcout = 0 ; frequency to write coordinates to xtc trajectory ; Neighbor searching nstlist = 10 ; frequency to update neighbor list. Neighborlist will be updated at least every 10 steps. Manual p80 ns_type = grid ; make a grid in the box and only check atoms in neighboring grid cells when constructing a new neighbor list every nstlist steps ; Electrostatics/VdW coulombtype = Shift ;PME rev ; tells gromacs how to model electrostatics. Shift: Coulomb/LJ potential is decreased over the whole range and forces decay smoothly to zero between vdw-type = Shift ; rcoulomb-switch/rvw-switch & rcoulomb/rvdw rcoulomb-switch = 0.9 ;0 ; where to start switching the Coulomb potential rvdw-switch = 0.9 ;0 ; where to start switching the LJ potential ; Cut-offs rlist = 1.1 ; in nm. Cut-off distance for short-range neighbor list rcoulomb = 1.1 ;1.0 ; distance for coulomb cut-off rvdw = 1.0 ; distance for coulomb cut-off ; Temperature coupling Tcoupl = v-rescale ;berendsen tc-grps = System ;HEX ; groups to couple to thermostat; Berendsen temperature coupling is on in these groups tau_t = 0.1 ;0.1 ; time constant for T coupling ref_t = 300 ;300 ; reference T for coupling. When you alter the T, don't forget to change the gen_temp for velocity generation ; Pressure coupling Pcoupl = Parrinello-Rahman ;berendsen Pcoupltype = semiisotropic ;isotropic ; isotropic: means the box expands and contracts in all directions(x,y,z)in order to maintain the proper pressure; semiisotropic: isotropic in x & y directions tau_p = 1 1 ;0.5 ; time constant for coupling in ps compressibility = 4.5e-5 4.5e-5 ; compressibility of solvent used in simulation in 1/bar ref_p = 30.0 30.0 ; reference P for coupling in bar ; Velocity generation Generate velocites is on at 300 K. Manual p155 gen_vel = yes ; generate velocites according to Maxwell distribution at T: gen_temp with random gen seed gen_seed gen_temp = 300.0 ; T for Maxwell distribution gen_seed = 173529 ; used to initialize random generator for random velocities ; Bonds constraints = none ;all-bonds ; sets the LINCS constraint for all bonds constraint-algorithm = lincs
-- gmx-users mailing list gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
