计算化学公社

 找回密码 Forget password
 注册 Register
Views: 983|回复 Reply: 1
打印 Print 上一主题 Last thread 下一主题 Next thread

[VMD] 想要统计离子与膜之间的静电作用,但是运行下面的代码会发生“vecsub: two vectors...

[复制链接 Copy URL]

1

帖子

0

威望

23

eV
积分
24

Level 1 能力者

想要统计离子与膜之间的静电作用,但是运行下面的代码会发生“vecsub: two vectors don't have the same size”报错,请问原因是?


# Load the necessary files
mol load psf XXXXX.psf
mol addfile XXXXXXX.dcd first 1000 last 4000 waitfor all


# Add ions to the system
set ion [atomselect top "resname XXX"]
set num_ions [$ion num]


set wat [atomselect top "resname XXX"]
set num_wats [$wat num]




set sel [atomselect top "all"]
$sel set beta 0.0
set num_atoms [expr {[molinfo top get numatoms] - $num_ions - $num_wats}]
$sel set beta 1.0


# Set up output file for energy vs. time data
set outfile [open "energy_vs_time.dat" w]
puts $outfile "Time (ps) Energy (kJ/mol)"


# Calculate the Coulombic interaction between ions and XXX
set sel_ion [atomselect top "resname XXX and name XXX"]
set sel_cof [atomselect top "segname XXXXX"]
set ion_charge 1.0
set dielectric 80.0
set rij_cutoff 12.0
set k_e 1389.3548
set start_frame 0
set end_frame [molinfo top get numframes]
for {set frame $start_frame} {$frame < $end_frame} {incr frame} {
    animate goto $frame
    set time [expr {[molinfo top get timesteps] * $frame / 1000.0}]
    set total_energy 0.0
    for {set i 0} {$i < $num_ions} {incr i} {
        set ion_pos [$sel_ion get {x y z}]
        set ion_q [expr {$ion_charge * $k_e}]
        for {set j 0} {$j < $num_atoms} {incr j} {
            set cof_pos [$sel_cof get {x y z}]
set rij [vecdist $ion_pos $cof_pos]
            if {$rij < $rij_cutoff} {
                set q_cof [$sel_cof get charge]
                set q_ion $ion_q
                set energy [expr {$k_e * $q_ion * $q_cof / ($dielectric * $rij)}]
                set total_energy [expr {$total_energy + $energy}]
            }
        }
    }
    puts $outfile "$time $total_energy"
}
close $outfile


























































































6万

帖子

99

威望

5万

eV
积分
120179

管理员

公社社长

2#
发表于 Post on 2023-6-7 10:01:26 | 只看该作者 Only view this author
一点点调试,弄清楚vecsub是哪一步涉及的,然后检查传进去的两个矢量
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

本版积分规则 Credits rule

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2025-8-17 01:46 , Processed in 0.205619 second(s), 26 queries , Gzip On.

快速回复 返回顶部 返回列表 Return to list