ORCA中如何中和DNA带电体系的电荷

Cadmus · 发表于 Post on 2023-10-13 17:58:08

各位老师好，我想使用ORCA对含有5个碱基对的DNA分子进行结构优化和单点能计算。在使用ORCA前，使用xtb进行了结构优化（隐式溶剂为水），发现在不设置--chrg参数前，优化后的结构会产生弯曲，如下图。
查找资料后发现DNA分子带局部电荷，设置--chrg -5后不再出现弯曲。
我想在xtb优化的基础上用ORCA使用水为隐式溶剂进行结构优化和单点能计算，请问各位老师我在ORCA中如何中和DNA分子的电荷，以避免分子弯曲？

sobereva · 发表于 Post on 2023-10-13 18:07:09

照样设净电荷就可以。要么你让PO4不处于去质子化状态，这时候都可以当真空计算

Cadmus · 发表于 Post on 2023-10-13 18:13:45

sobereva 发表于 2023-10-13 18:07
照样设净电荷就可以。要么你让PO4不处于去质子化状态，这时候都可以当真空计算

好的老师，我试试，谢谢您

wzkchem5 · 发表于 Post on 2023-10-13 20:31:34

严格的处理应该是加一些显式水分子，并用Na+中和电荷，显式水分子和Na+可以用MM描述。否则DNA结构尽管不弯曲，但其他方面未必可靠

Cadmus · 发表于 Post on 2023-10-16 07:35:45

wzkchem5 发表于 2023-10-13 20:31
严格的处理应该是加一些显式水分子，并用Na+中和电荷，显式水分子和Na+可以用MM描述。否则DNA结构尽管不弯 ...

好的老师，我试试这个方案，谢谢您

Cadmus · 发表于 Post on 2023-10-17 17:40:30

wzkchem5 发表于 2023-10-13 20:31
严格的处理应该是加一些显式水分子，并用Na+中和电荷，显式水分子和Na+可以用MM描述。否则DNA结构尽管不弯 ...

老师好，向您请教个问题。我按照您说的这个方案，首先用MD跑了一个轨迹，抽取了其中一帧，然后用GROMACS使用parmbsc1力场生成了一个top文件（包含水分子），再用parmed转换成了prmtop文件。然后用orca_ff转换成了.ORCAFF.prms文件。用orca做几何优化，但是报错ORCA finished by error termination in GSTEP，请您看看我是哪里出了问题，谢谢老师！
inp文件：（大概有10000个原子，服务器内存1024G）

! QMMM BLYP D3 GCP(DFT/SVP) def2-SVP def2/J opt noautostart miniprint nopop
%maxcore 10000
%pal nprocs 20 end
%qmmm
QMAtoms {0:314} end
ORCAFFFilename "clusters.pdb01.ORCAFF.prms"
end
* xyz 0 1
O 10.05900000 -3.83100000 -6.46100000
H 10.71900000 -4.23100000 -7.09100000
C 10.39900000 -2.48100000 -6.33100000
H 10.50900000 -2.05100000 -7.31100000
H 11.32900000 -2.48100000 -5.75100000
C 9.39900000 -1.81100000 -5.39100000
H 9.79900000 -0.76100000 -5.27100000
O 8.08900000 -1.92100000 -5.88100000
C 7.29900000 -2.31100000 -4.79100000
H 7.13900000 -1.36100000 -4.22100000
N 5.97900000 -2.89100000 -5.19100000
C 5.60900000 -4.18100000 -5.35100000
H 6.11900000 -5.07100000 -4.99100000
N 4.40900000 -4.32100000 -5.88100000
C 3.96900000 -3.03100000 -6.01100000
C 2.71900000 -2.50100000 -6.50100000
...
...
...
O 14.94900000 17.53900000 24.34900000
H 15.72900000 17.71900000 23.82900000
H 15.27900000 17.06900000 25.11900000
K -9.97100000 -4.60100000 14.69900000
K -10.95100000 -6.92100000 10.72900000
K 13.82900000 15.49900000 30.83900000
K -10.54100000 18.70900000 -0.68100000
K -6.64100000 16.50900000 4.23900000
K -14.46100000 -15.62100000 38.25900000
K -3.09100000 9.00900000 -1.87100000
K 2.27900000 -8.76100000 23.22900000
*

复制代码

out文件：

...
...
...
Active atoms ... All atoms
Constrain bonds and angle of MM waters ... NO
*****************************
* Geometry Optimization Run *
*****************************
Geometry optimization settings:
Update method Update .... BFGS
Choice of coordinates CoordSys .... Z-matrix Internals
Initial Hessian InHess .... Almoef's Model
Convergence Tolerances:
Energy Change TolE .... 5.0000e-06 Eh
Max. Gradient TolMAXG .... 3.0000e-04 Eh/bohr
RMS Gradient TolRMSG .... 1.0000e-04 Eh/bohr
Max. Displacement TolMAXD .... 4.0000e-03 bohr
RMS Displacement TolRMSD .... 2.0000e-03 bohr
Strict Convergence .... False
------------------------------------------------------------------------------
ORCA OPTIMIZATION COORDINATE SETUP
------------------------------------------------------------------------------
The optimization will be done in new redundant internal coordinates
Making redundant internal coordinates ... (new redundants)
ORCA finished by error termination in GSTEP
Calling Command: /public/home/cxc/orca504/orca_gstep clusters.ginp.tmp
[file orca_tools/qcmsg.cpp, line 465]:
.... aborting the run

复制代码

wzkchem5 · 发表于 Post on 2023-10-17 17:47:48

Cadmus 发表于 2023-10-17 10:40
老师好，向您请教个问题。我按照您说的这个方案，首先用MD跑了一个轨迹，抽取了其中一帧，然后用GROMACS ...

体系太大了，默认的内坐标优化器跑不动，应该把Opt改成L-Opt，用L-BFGS算法在直角坐标下优化

Cadmus · 发表于 Post on 2023-10-17 19:11:09

wzkchem5 发表于 2023-10-17 17:47
体系太大了，默认的内坐标优化器跑不动，应该把Opt改成L-Opt，用L-BFGS算法在直角坐标下优化

好的，谢谢您！

Cadmus · 发表于 Post on 2023-10-18 15:21:37

本帖最后由 Cadmus 于 2023-10-18 17:19 编辑

wzkchem5 发表于 2023-10-17 17:47
体系太大了，默认的内坐标优化器跑不动，应该把Opt改成L-Opt，用L-BFGS算法在直角坐标下优化

老师您好，我改成L-Opt后几何优化开始运行了，但是运行了大概100多步后ORCA崩溃了，也没有输出什么结构可以人工检查。想请问您这样应该怎么处理？下面是截取的一些log文件信息，输入的文件中将QM的部分改成GFN2-XTB来计算了，完整的log文件上传附件。

Defined regions:
all: 11912 atoms
qm: 315 atoms
mm: 11597 atoms
active: 11912 atoms
Will write restart file to "clusters_final.mdrestart" after each step.
Will write energy statistics to "clusters_final-minimize-ener.csv" after each step.
The system contains 35736 active degrees of freedom (11912 active atoms * 3).
This is a 35736-dimensional optimization problem.
Allocating 14.49 MiB of memory for L-BFGS working space...
| t_SCF | t_Grad | Temp | E_Pot | Delta E | RMS(Grad) | Max(Grad) | Max(Step)
Step | [s] | [s] | [K] | [Hartree] | [Hartree] | [kJ/(mol*A)] | [kJ/(mol*A)] | [Angstrom]
------------|--------|--------|----------|-----------------|--------------|--------------|--------------|--------------
1 17.2 0.2 6724.02477224
2 7.0 0.1 6648.47954071 -75.54523152 5706.0325 313823.5542 0.158983
3 8.7 0.1 6572.96014517 -75.51939554 5674.6626 312191.0077 0.237112
4 8.4 0.2 6498.18941686 -74.77072831 5643.9818 310595.9461 0.175621
5 8.3 0.1 6423.84604817 -74.34336869 5612.5361 309000.6667 0.172805
6 8.1 0.1 6349.90012622 -73.94592195 5582.2098 307405.4352 0.166918
7 7.9 0.1 6276.33335404 -73.56677218 5553.2103 305810.2773 0.166985
8 8.5 0.1 6203.14139928 -73.19195476 5523.5595 304215.0378 0.173139
9 8.1 0.2 6130.35503062 -72.78636866 5497.6284 302619.8276 0.225309
10 8.1 0.2 6057.91608822 -72.43894240 5466.2542 301024.5725 0.200905
11 8.4 0.1 5985.89056611 -72.02552211 5438.0600 299429.3217 0.181453
12 7.8 0.2 5914.26172225 -71.62884386 5408.2646 297834.0512 0.169440
...
...
...
120 9.0 0.2 969.38198049 38.20202116 2766.7444 146514.3079 0.274969
121 8.2 0.2 915.18657051 -54.19540998 2725.7507 147562.0515 0.177773
122 8.3 0.2 904.54227616 -10.64429436 2605.5680 145992.3962 0.274969
123 7.8 0.1 897.14643187 -7.39584429 2618.0783 146643.0328 0.114131
124 7.8 0.1 850.56683880 -46.57959307 2746.3052 145039.0807 0.274969
125 8.5 0.1 901.16216312 50.59532431 2608.0345 143466.1624 0.274969
126 8.7 0.1 849.45424282 -51.70792029 2651.2302 144827.2925 0.237842
127 8.2 0.1 812.63597629 -36.81826653 9167.3952 679989.1396 0.274969
*** SIGSEGV caught: Segmentation Fault ***
ORCA_MD apparently has crashed :-( Sorry.
Please send the log file ("clusters_final.md.log") to the developers, enabling them to fix this error.
Decoding line numbers... ("addr2line <addr> -f -i -C -e /public/home/orca504/orca_md")
# 20 - [ 0x5a8143] - ?? - ??:0
# 19 - [ 0x5a857d] - ?? - ??:0
# 18 - [ 0x5a86b2] - ?? - ??:0
# 16 - [ 0x8c8f59] - ?? - ??:0
# 15 - [ 0x8eb3d9] - ?? - ??:0
# 14 - [ 0x5bb871] - ?? - ??:0
# 13 - [ 0x530522] - ?? - ??:0
# 12 - [ 0x513cfa] - ?? - ??:0
# 11 - [ 0x5129e9] - ?? - ??:0
# 10 - [ 0x550462] - ?? - ??:0
# 9 - [ 0x54ef15] - ?? - ??:0
# 8 - [ 0x56908f] - ?? - ??:0
# 7 - [ 0x568c27] - ?? - ??:0
# 6 - [ 0x5687ef] - ?? - ??:0
# 5 - [ 0x5685fc] - ?? - ??:0
# 4 - [ 0x58ba0a] - ?? - ??:0
# 1 - [ 0x509d76] - ?? - ??:0
Delivering control to calling process.
ERROR : MD Program returns an error
cannot continue the md simulation
COMMAND : /public/home/orca504/orca_md clusters_final orca
RETURN CODE : 134

复制代码

wzkchem5 · 发表于 Post on 2023-10-18 15:53:29

Cadmus 发表于 2023-10-18 08:21
老师您好，我改成L-Opt后几何优化开始运行了，但是运行了大概100多步后ORCA崩溃了，也没有输出什么结构可 ...

检查是不是初始结构不合理，例如有原子靠得太近，某些力场参数设置错误，等等。正常结构优化，哪怕是你这么大的体系，能量不太可能下降几千Hartree。

Cadmus · 发表于 Post on 2023-10-18 17:18:59

wzkchem5 发表于 2023-10-18 15:53
检查是不是初始结构不合理，例如有原子靠得太近，某些力场参数设置错误，等等。正常结构优化，哪怕是你这 ...

好的，谢谢您

		自动登录 Automatic login	找回密码 Forget password
密码 Password			注册 Register

[ORCA] ORCA中如何中和DNA带电体系的电荷