[akitao]

東工大の北尾です。お世話になっております。昨日うまく投稿できなかったようなので再投稿させていただきます。2重投稿になってしまっていた場合は申し訳ありません。

GENESIS2.0 beta (Nov)版 spdynでNPT MDを実行しました。

AMBER pmemdで等方的NPTで平衡化後のあと、GENESISIでVVER, NPE, NHC, ISOの条件で更に平衡化しました。元々BOXの大きさはX=Y>Zです。

既にNPで平衡化しているので、全体の体積はあまり変化しないのですがGENESISに移行後XYが減少し、Zが増加します。長く計算するとBOXサイズは一定になりますが、初期設定からBOXのXYとZの比率が大きく変わってしまいます。ISOと設定しているはずなのですが、実際にはSEMI-ISOになっているように見えます。ログ上ではISOと認識されているように見えます。ちなみに別の系でも同じようなことが起こりました。またSEMI-ISOと設定した時にはSEMI-ISOになっているようです。

設定の組み合わせに問題があるのではないかと予想していますが、1週間ほど解決できていないのでご相談しました。何かお気づきの点があればアドバイスいただければ幸いです。よろしくお願いいたします。

[INPUT]
prmtopfile = lpr.prmtop   # topology file
ambcrdfile = equ102.rst7  # protein structure file
ambreffile = lpr.crd      # protein structure file

[OUTPUT]
dcdfile = equ103.dcd  # DCD trajectory file
rstfile = equ103.rst  # restart file

[ENERGY]
forcefield    = AMBER # AMBER force field
electrostatic = PME   # use Particle mesh Ewald method
switchdist    = 9.0   # switch distance
cutoffdist    = 9.0   # cutoff distance
pairlistdist  = 10.5  # pair-list distance
pme_nspline   = 4     # order of B-spline in [PME]
pme_max_spacing = 1.2 # max grid spacing allowed

[DYNAMICS]
integrator    = VVER      # [LEAP,VVER]
nsteps        = 2500000   # number of MD steps 1ns
timestep      = 0.002     # timestep (ps)
eneout_period = 500       # energy output period
crdout_period = 0          # coordinates output period
rstout_period = 2500000   # restart output period


[CONSTRAINTS]
rigid_bond  = YES # Turn on SHAKE/RATTLE
fast_water  = YES # Turn on SETTLE
water_model = WAT # residue name of the rigid water

[ENSEMBLE]
ensemble    = NPT      # Isothermal-isobaric ensemble
tpcontrol   = NHC      # NHC thermostat and barostat
temperature = 300.0    # target temperature (K)
pressure    = 1.0      # target pressure (atm)
isotropy    = ISO      # isotropic

[BOUNDARY]
type          = PBC       # [PBC,NOBC]
box_size_x    = 279.8376514  # box size (x) in [PBC]
box_size_y    = 279.8376514  # box size (y) in [PBC]
box_size_z    = 235.7621046  # box size (z) in [PBC]

[SELECTION]
group1        = an:CA  # index of restraint group 1

[RESTRAINTS]
nfunctions    = 1     # number of functions
function1     = POSI  # restraint function type
constant1     = 0.1  # force constant
select_index1 = 1     # restrained group

Attachments:

You must be logged in to view attached files.

[jwjung]

Dear Kitao sensei,

This is Jaewoon Jung in genesis developer team. I’m sorry to reply in English.

I checked the attached output files, and it seems that the ratio does not change.
At 0-th step, the box size in x, y, and z dimensions are 279.8376, 279.8376, and 235.7621. The ratio of box size between x and z dimensions is 279.8376/235.7621=1.18694904737. At 2000-th step, the box size in each dimension is 279.0337, 279.0337, and 235.0873. The ratio is 279.0337/235.0873=1.18693651252. The ratio is slightly changed, but it seems to be just due to insufficient number of digits in outputs.

Because we don’t have the same inputs, we tried to write outputs of barostat momentum given the BOXX=BOXY>BOXZ condition from our inputs, and found that the barostat momentum and the following box size scaling factors are same in all dimensions.

Just from the output, box size in dimension is also decreasing like those in x and y dimensions. If you observe the unexpected phenomena in long time MD (not from 2000 steps in the attached file), could you show us the log (part of the log is also okay)?

Beside the box size ratio issue written above, I have a few suggestions from the input.
1. In [ENSEMBLE] section, you can obtain better performance by writing group_tp=yes. The reference related to this can be found from https://aip.scitation.org/doi/full/10.1063/5.0027873 .
2. If you want to run NPT MD with positional restraint, how about including the position restraint virial in pressure? You can do it by writing pressure_position=yes in [RESTRAINTS] section.

Please let us know if our answer is not sufficient or you have additional questions/comments.

Best regards,
Jaewoon Jung

[akitao]

Dear Jung-san,

Thank you  very much for quickly responding to my question. I am also grateful for useful suggestions for the input. I will certainly follow your suggestions.

Actually I have done almost 100 ns MD before noticing this phenomenon. I will attach long-term plots of box lengths and volume as a PDF file. I will also attach the part of log file in the later run. I can share the full log but seems to be too big.

Also, is it possible that not including position restraint virial causes this kind of problem? I will try to run GENESIS after improving this setup.

Thank you very much for your kind support,

Akio Kitao

Attachments:

You must be logged in to view attached files.

[jwjung]

Dear Kitao-sensei,

Many thanks for the outputs.

First of all, position restraint virial is not related to the problem. It is just my recommendation.

One thing I’m not sure is the difference of box size between the figure and output files. For example, the initial box size is (276.0087,276.0087,240.6425) in the output file, but it seems to be about (279,279,235) in the attached figure. There are also differences between the output and figure in the box size at 50 ns, too. If you’re okay, could you check it?

Regardless of the difference between the output and figure, it would be a problem if the box size change is different in each dimension although isotropic condition is assigned. If you don’t mind, could you share the log up to 5 ns run? I think the file size could be too big to be attached in this webpage. I guess you can send it with another ways, i.e. via mail or box. My e-mail address is jung@riken.jp.

If you don’t mind, I’d like to ask you to run with tpcontrol=BUSSI up to 10 ns to understand if the same phenomena happens.

We’re sorry for the trouble and asking, and consider the reason of the problem and how to fix it.

Best regards,
Jaewoon Jung

[akitao]

Dear Jung-san

Thank you again for your reply

jwjung wrote:

One thing I’m not sure is the difference of box size between the figure and output files. For example, the initial box size is (276.0087,276.0087,240.6425) in the output file, but it seems to be about (279,279,235) in the attached figure. There are also differences between the output and figure in the box size at 50 ns, too. If you’re okay, could you check it?

Sorry for causing confusion. I sent the log of the first run in my first message (0-5 ns), but in my second message, I sent the log of the second run (5-10 ns). This must be the reason of these differences.

I have confirmed that the box size of the very initial run was correctly set, which means that the box sizes in the AMBER formatted restart file, those in box_size_x, y, z IN [BOUNDARY] and those of STEP 0 in log file agree exactly except for digits not shown in the log.

I will send the entire log by e-mail. I can also share all the input/output if you need.

I am now running the equivalent MD but with BUSSI. I will let you know the results once it finishes.

Thank you very much,
Akio Kitao

[akitao]

Dear Jung-san,

Thank you very much for identifying the cause of the problem.

I recompiled genesis as you suggested, ./configure –enable-mixed –host=Fugaku

I have confirmed that with mixed precision, pressure control realized ISO condition correctly.

Thank you very much for your kind support.

Best regards,

Akio Kitao

[jwjung]

Dear Kitao-sensei,

Thank you for the feedback.
It’s good to see that the problem is solved by using the mixed precision model.

Best regards,
Jaewoon Jung