CN113593652A - 一种约束原子速度的方法及装置 - Google Patents

一种约束原子速度的方法及装置 Download PDF

Info

Publication number
CN113593652A
CN113593652A CN202110955870.4A CN202110955870A CN113593652A CN 113593652 A CN113593652 A CN 113593652A CN 202110955870 A CN202110955870 A CN 202110955870A CN 113593652 A CN113593652 A CN 113593652A
Authority
CN
China
Prior art keywords
constraint
atom
atoms
atomic
constrained
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110955870.4A
Other languages
English (en)
Other versions
CN113593652B (zh
Inventor
郭振江
郭璟
柳宇驰
郭晨
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Si Lang Science And Technology Co ltd
Original Assignee
Beijing Si Lang Science And Technology Co ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Si Lang Science And Technology Co ltd filed Critical Beijing Si Lang Science And Technology Co ltd
Priority to CN202110955870.4A priority Critical patent/CN113593652B/zh
Publication of CN113593652A publication Critical patent/CN113593652A/zh
Application granted granted Critical
Publication of CN113593652B publication Critical patent/CN113593652B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16CCOMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C10/00Computational theoretical chemistry, i.e. ICT specially adapted for theoretical aspects of quantum chemistry, molecular mechanics, molecular dynamics or the like
    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16CCOMPUTATIONAL CHEMISTRY; CHEMOINFORMATICS; COMPUTATIONAL MATERIALS SCIENCE
    • G16C20/00Chemoinformatics, i.e. ICT specially adapted for the handling of physicochemical or structural data of chemical particles, elements, compounds or mixtures
    • G16C20/50Molecular design, e.g. of drugs

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Computing Systems (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Health & Medical Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Chemical & Material Sciences (AREA)
  • Medicinal Chemistry (AREA)
  • Pharmacology & Pharmacy (AREA)
  • Crystallography & Structural Chemistry (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明提供一种约束原子速度的方法及装置,通过确定刚性模型中的原子构型,基于原子构型确定约束原子对并进行分组得到约束原子对组,基于每个组内约束条件构建N×N矩阵和N阶矢量;基于N×N矩阵和N阶矢量计算每个约束原子对对应的约束力;在每个约束原子对下,基于约束力、单步迭代时间和待约束原子的原子质量进行计算,得到约束精度;判断约束精度是否小于预设约束精度,若否,则重新构建N阶矢量,并执行迭代计算,直至当前得到的原子速度小于所述预设约束精度,结束迭代。在本方案中,通过快速收敛每个时刻刚性整体内原子速度,并将计算得到的速度用于下一时刻的原子运动,从而实现提高刚性模型约束原子速度的计算效率的目的。

Description

一种约束原子速度的方法及装置
技术领域
本发明涉及分子模拟技术领域,尤其涉及一种约束原子速度的方法及装置。
背景技术
在分子模拟计算中,存在不同类型的分子模型,如振动模型、刚性模型等。其中,刚性模型采用特殊的方法约束其原子间距离及速度。
目前最被认可的约束速度的方法为使用RATTLE方法。使用RATTLE方法约束刚性模型原子间速度,可以准确地获得刚性分子内原子在约束条件下的速度。由于使用RATTLE方法刚性模型每次迭代只考虑单一刚性约束原子对的速度约束,并没有考虑其它相互关联耦合的刚性原子的影响,为了使所有刚性模型内的原子都满足给定的约束精度,需要很多次迭代计算。而迭代计算次数过多则导致产生刚性模型约束原子速度的计算效率大大降低的问题。
发明内容
有鉴于此,本发明实施例提供一种约束原子速度的方法及装置,迭代计算中同时考虑相互关联耦合的刚性原子的影响,以实现在不降低计算精度的前提下减少迭代次数,从而提高刚性模型约束原子速度的计算效率的目的。
为实现上述目的,本发明实施例提供如下技术方案:
本发明实施例第一方面公开了一种约束原子速度的方法,所述方法包括:
确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成;
按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列;
利用所述每个约束原子对组基于相关约束条件构建N×N矩阵,以及构建N阶矢量,N的取值为约束原子对组内约束原子对的个数;
基于所述N×N矩阵和N阶矢量计算每个约束原子对内两个原子对应的约束力;
在每个约束原子下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到约束精度;
判断所述约束精度是否小于预设约束精度,若否,则返回执行构建N阶矢量这一步骤进行迭代计算,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。
可选的,利用所述约束原子对组基于相关约束条件构建N×N矩阵,以及构建N阶矢量,包括:
获取每个约束原子对对应的约束序号k,k∈[1,N];
获取每个约束原子对中两个关联耦合的原子k1和k2的相对矢量位置
Figure BDA0003220248380000026
,以及原子速度
Figure BDA0003220248380000021
针对每一约束原子对、基于单步迭代时间Δt、原子k1和k2的相对矢量位置
Figure BDA0003220248380000022
和原子速度
Figure BDA0003220248380000025
构建N阶矢量c,其中,第k阶的元素为
Figure BDA0003220248380000023
基于每一约束原子对组对应的N阶矢量和约束条件构建N×N矩阵,得到矩阵A
其中,所述矩阵A内的元素表示为:
Figure BDA0003220248380000024
其中,k和k′为约束原子对组内两个约束序号,k1和k2为约束序号k约束内的两个原子,k1′和k2′为k′约束内的两个原子,δ为狄拉克函数。
可选的,基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内两个原子对应的约束力,包括:
确定每个约束原子对组的N阶矢量c;
针对每个约束原子对组,基于对应的N阶矢量c和所述N×N矩阵的逆矩阵的乘积确定组内所有约束原子对对应的N阶矢量约束力集合B,其中,B=A -1c,A -1为所述N×N矩阵的逆矩阵。
可选的,所述在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,包括:
在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到对应的原子速度修正量;
基于所述原子速度修正量修正本次计算中的待约束原子的原子速度,得到本次计算后的待约束原子的原子速度和约束精度,在迭代未结束前,使本次计算后的待约束原子的原子速度参与下一次计算;
其中,构成约束原子对的两个原子k1和k2各自对应的原子速度修正量为:
Figure BDA0003220248380000031
Figure BDA0003220248380000032
其中,所述约束原子对的约束力
Figure BDA0003220248380000033
Figure BDA0003220248380000034
为矩阵A逆矩阵的第k行,Bk表示约束序号为k的约束原子对间的约束力,Δt为单步迭代的时间,mk1为原子k1的原子质量,mk2为原子k2的原子质量;
本次迭代后约束原子的速度为:
Figure BDA0003220248380000035
Figure BDA0003220248380000036
Figure BDA0003220248380000037
为本次计算后约束原子对k的约束精度。
本发明实施例第二方面公开了一种约束原子速度的装置,所述装置包括:
确定模块,用于确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成;
分组模块,用于按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列;
构建模块,用于基于所述每个约束原子对组内约束原子对的约束条件构建N×N矩阵,以及构建N阶矢量,N的取值为约束原子对组内约束原子对的个数;
约束力计算模块,用于基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内的两个原子对应的约束力;
迭代处理模块,用于在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到约束精度;判断所述约束精度是否小于预设约束精度,若否,则返回构建模块,使所述构建模块构建N阶矢量,并进行迭代计算,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。
可选的,所述构建模块,具体用于:
获取每个约束原子对对应的约束序号k,其中,k∈[1,N];
获取每个约束原子对中两个关联耦合的原子k1和k2的相对矢量位置
Figure BDA0003220248380000041
和矢量速度
Figure BDA0003220248380000042
针对每一约束原子对、基于单步迭代时间Δt、原子k1和k2的相对矢量位置
Figure BDA0003220248380000043
和矢量速度
Figure BDA0003220248380000046
构建N阶矢量c,基于每一约束原子对组对应的N阶矢量和约束条件构建N×N矩阵,得到矩阵A
其中,矢量c的第k阶元素
Figure BDA0003220248380000044
所述矩阵A内的元素表示为:
Figure BDA0003220248380000045
其中,k和k′为约束原子对组内的两个约束序号,k1和k2为约束序号k约束内的两个原子,k1′和k2′为k′约束内的两个原子,δ为狄拉克函数。
可选的,所述约束力计算模块,具体用于:
确定每个约束原子对组的N阶矢量c;针对每个约束原子对组,基于对应的N阶矢量c和所述N×N矩阵的逆矩阵的乘积确定组内所有约束原子对对应的N阶矢量约束力集合B,其中,B=A -1c,A -1为所述N×N矩阵的逆矩阵。
可选的,所述在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算的处理模块,具体用于:
在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到对应的原子速度修正量;基于所述原子速度修正量修正本次计算中的待约束原子的原子速度,得到本次计算后的待约束原子的原子速度和约束精度,在迭代未结束前,使本次计算后的待约束原子的原子速度参与下一次计算;
其中,构成约束原子对的两个原子k1和k2各自对应的原子速度修正量为:
Figure BDA0003220248380000051
Figure BDA0003220248380000052
其中,所述约束原子对的约束力
Figure BDA0003220248380000053
Figure BDA0003220248380000054
为矩阵A逆矩阵的第k行,Bk表示约束序号为k的约束原子对间的约束力,Δt为单步迭代的时间,mk1为原子k1的原子质量,mk2为原子k2的原子质量;
本次迭代后约束原子的速度为:
Figure BDA0003220248380000055
Figure BDA0003220248380000056
Figure BDA0003220248380000057
为本次计算后约束原子对k的约束精度。
基于上述本发明实施例提供的一种约束原子速度的方法及装置,通过确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成;按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列;基于所述组内约束原子对构建N×N矩阵,N的取值为约束原子对的个数;基于所述N×N矩阵和N阶矢量计算每个约束原子对内原子对应的约束力;在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行迭代计算;判断每一次迭代计算得到的约束精度是否小于预设约束精度,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。
在本方案中,由于N×N矩阵各元素的值在迭代过程中不会改变,因此不需要重新计算矩阵求逆,同时,在迭代计算过程中同时考虑相互有关联的约束原子对的相互影响,统一在矩阵中计算相互关联约束原子对内的原子的作用,从而可以大大降低约束迭代次数,且通过快速收敛每个时刻刚性整体内原子速度,并将计算得到的速度用于下一时刻的原子运动,使得两两原子之间运动的相对速度与二者的相对位置正交,从而实现提高刚性模型约束原子速度的计算效率的目的。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据提供的附图获得其他的附图。
图1为本发明实施例公开的一种原子相互关联刚性结构示意图;
图2为本发明实施例公开的一种约束原子速度的方法的流程图;
图3为本发明实施例公开的一种约束原子速度的装置的结构示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
在本申请中,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
在分子动力学模拟计算中,需要计算原子的速度及位置随着时间的演化情况。在一个时刻下,需要知道原子在该时刻的位置及速度,依据每个原子的受力求解牛顿方程,获得原子下一时刻的位置及速度,以此周而复始。
对于含氢原子的分子,氢原子与其它原子连接的共价键的键势非常大,所以一般都认为其键长(氢原子与相连原子间的距离)是固定不变的。那么在分子模拟计算中,就要将氢原子和与其相连的原子认为是一个整体来考虑,而且,经常一个原子会与多个氢原子相连接,则需要将其全部看成一个刚性整体。一般如水分子(包含一个氧原子和两个氢原子相连),则需要看成3个原子的整体来考虑,有机物上的甲基(一个碳原子和3个氢原子相连),则需要看成4个原子的整体来考虑。当几个原子被看成一个整体,那么在原子运动过程中,其速度需要满足一定的条件:两两原子之间运动的相对速度必须与二者的相对位置正交。
如图1所示,为原子相互关联刚性结构示意图。其中,黑色细线所指示的键为刚性约束键,刚性约束键两头的原子为刚性原子对,黑色虚线圈内表示一个相互关联的刚性结构。
由背景技术可知,使用RATTLE方法刚性模型每次迭代只考虑一个刚性原子对的速度约束,并没有考虑其它相互关联耦合的其它刚性原子的影响,为了使得所有刚性模型内的原子都满足给定的约束精度,需要很多次迭代计算。而迭代次数过多则导致产生刚性模型约束原子速度的计算效率大大降低的问题。
因此,本发明实施例公开了一种约束原子速度的方法及装置,通过快速收敛每个时刻刚性整体内原子速度,并将计算得到的速度用于下一时刻的原子运动,使得两两原子之间运动的相对速度与二者的相对位置正交,从而实现提高刚性模型约束原子速度的计算效率的目的。
如图2所示,为本发明实施例公开的一种约束原子速度的方法的流程示意图,所述方法包括如下步骤:
S201:确定刚性模型中的原子构型,基于所述原子构型确定约束原子对。
在S201中,每个约束原子对由两个关联耦合的原子构成。
S202:按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组。
在S202中,每一约束原子对组包含不同原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列。
在具体实现S202的过程中,将具有相同原子的约束原子对划分至一个约束原子对组内。也就是说,不同的约束原子对组内所包含的原子不同,同一个原子不会出现在两个或两个以上的约束原子对组内。
例如,约束原子对包括:1-2、1-3、2-4、5-7、5-8、7-8、9-10和9-11。
根据每个约束原子对所包含的原子,将具有相同原子的约束原子对划分至一个约束原子对组内。具体划分后可以得到三个约束原子对组。
其中,第一组:1-2、1-3和2-4;约束原子对数为3,序号为1,2和3。
第二组:5-7、5-8和7-8;约束原子对数为3,序号为1,2和3。
第三组:9-10、9-11;约束原子对数为2,序号为1和2。
S203:基于所述每个约束原子对组基于相关约束条件构建N×N矩阵,以及构建N阶矢量。
在S203中,N的取值为约束原子对的个数。
在具体执行S203的过程中,首先,获取每个约束原子对对应的约束序号k。
其中,k的取值范围为[1,N]。
其次,获取每个约束原子对中两个关联耦合的原子k1和k2的相对矢量位置
Figure BDA0003220248380000091
以及原子速度
Figure BDA0003220248380000092
其次,针对组内每一约束原子对、基于单步迭代时间Δt、原子k1和k2的相对矢量位置
Figure BDA0003220248380000099
和原子速度
Figure BDA0003220248380000093
构建约束原子对组的N阶矢量c。
其中,公式(1)示出了c的第k阶元素。
Figure BDA0003220248380000094
最后,基于组内每一约束原子对约束条件构建N×N矩阵,得到矩阵A
其中,矩阵A为N×N矩阵,所述矩阵A内的元素如公式(2)所示。
Figure BDA0003220248380000095
其中,k和k′为耦合约束内的两个约束序号。例如,矩阵A的第2行第3列的元素为A23,则k=2,k′=3。
k1和k2为约束序号k约束内的两个原子,k1′和k2′为k′约束内的两个原子,
Figure BDA0003220248380000096
为原子k1和k2的相对矢量位置,
Figure BDA0003220248380000097
为原子k1′和k2′的相对矢量位置,δ为狄拉克函数。
该狄拉克函数如公式(3)所示。
Figure BDA0003220248380000098
S204:基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内两个原子对应的约束力。
在具体执行S204的过程中,根据S203确定的矩阵A和N阶矢量c进行计算,具体计算过程为:
首先,确定每个约束原子对组的N阶矢量c。
然后,基于公式(4),针对每个约束原子对组,基于对应的N阶矢量c和所述N×N矩阵的逆矩阵的乘积确定组内每个约束原子对对应的约束力集合,即N阶矢量约束力B。
B=A -1c (4)
其中,A -1为所述N×N矩阵的逆矩阵。
S205:在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行迭代计算。
在具体执行S205的过程中,具体包括如下步骤:
首先,在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到对应的原子速度修正量。
其中,构成约束原子对的两个原子k1和k2各自对应的原子速度修正量分别为公式(5)和公式(6)所示:
Figure BDA0003220248380000101
Figure BDA0003220248380000102
其中,所述约束原子对的约束力Bk如公式(7)所示:
Figure BDA0003220248380000103
其中,
Figure BDA0003220248380000104
为矩阵A逆矩阵的第k行,Bk表示约束序号为k的约束原子对间的约束力,也是所有约束原子对对应的N阶矢量约束力集合B的第k个元素,Δt为单步迭代的时间,mk1为原子k1的原子质量,mk2为原子k2的原子质量。
然后,基于所述原子速度修正量修正本次计算中的待约束原子的原子速度,得到本次计算后的待约束原子的原子速度和约束精度,在迭代未结束前,使本次计算后的待约束原子的原子速度参与下一次计算;
若待约束原子为k1和k2,也就是说,基于
Figure BDA0003220248380000105
修正本次计算中的待约束原子k1的原子速度,基于
Figure BDA0003220248380000106
修正本次计算中的待约束原子k2的原子速度。
由此得到本次计算后的待约束原子的原子速度为:
Figure BDA0003220248380000111
Figure BDA0003220248380000112
以待约束原子k1和k2为例,本次计算后的待约束原子的约束精度为:
Figure BDA0003220248380000113
S206:判断所述约束精度是否小于预设约束精度,若否,返回执行S203重新构建N阶矢量,若是,则结束迭代,并执行S207。
在具体执行S206的过程中,针对每次计算,将得到的约束精度
Figure BDA0003220248380000114
与预设约束精度Δ进行比较,若
Figure BDA0003220248380000115
则返回执行S203重新构建N阶矢量,但不重新构建N×N矩阵,然后执行迭代运算,直至
Figure BDA0003220248380000116
则结束迭代。
S207:确定当前得到的原子速度为原子最终速度。
在S207中,确定当前得到的原子速度小于所述预设约束精度,结束迭代计算,使所述待约束原子基于所述原子最终速度进行运动。
基于本发明实施例公开的一种约束原子速度的方法,通过确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成;按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列;基于所述组内约束原子对构建N×N矩阵,N的取值为约束原子对的个数;基于所述N×N矩阵和N阶矢量计算每个约束原子对内原子对应的约束力;在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行迭代计算;判断每一次迭代计算得到的约束精度是否小于预设约束精度,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。在本方案中,由于N×N矩阵各元素的值在迭代过程中不会改变,因此不需要重新计算矩阵求逆,同时,在迭代计算过程中同时考虑相互有关联的约束原子对的相互影响,统一在矩阵中计算相互关联约束原子对内原子的作用,从而可以大大降低约束迭代次数,且通过快速收敛每个时刻刚性整体内原子速度,并将计算得到的速度用于下一时刻的原子运动,使得两两原子之间运动的相对速度与二者的相对位置正交,从而实现提高刚性模型约束原子速度的计算效率的目的。
基于上述本发明实施例公开的一种约束原子速度的方法,以水分子举例说明。
水分子结构包括3个原子,3个原子构型通常为三角约束结构,原子序号分别为1、2和3。
因此,由两两原子构成的约束原子对分别为:1-2、1-3和2-3。
基于公式(2)可知,公式AB=c的形式为:
Figure BDA0003220248380000121
根据公式(4),可以确定该刚性模型中的所有约束原子对下的约束力B,然后在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行迭代计算,并判断每一次迭代计算得到的约束精度是否小于预设约束精度,如果约束精度不小于预设约束精度,则继续迭代,如果约束精度小于预设精度,则结束迭代。
在该过程中,由于矩阵A不会改变,因此不需要重新计算矩阵求逆。通过本发明实施例公开的约束原子速度的方法,可以大大降低约束迭代次数,实现提高刚性模型约束原子速度的计算效率的目的。
基于上述本发明实施例公开的一种约束原子速度的方法,本发明实施例还公开了一种约束原子速度的装置,如图3所示,该装置包括:确定模块301,分组模块302,构建模块303,约束力计算模块304和处理模块305。
确定模块301,用于确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成。
分组模块302,用于按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列。
构建模块303,用于基于所述每个约束原子对组内约束原子对的约束条件构建N×N矩阵,以及构建N阶矢量,N的取值为约束原子对组内约束原子对的个数。
约束力计算模块304,用于基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内的两个原子对应的约束力。
处理模块305,用于在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到约束精度;判断所述约束精度是否小于预设约束精度,若否,则返回构建模块,使所述构建模块构建N阶矢量,并进行迭代计算,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。
所述构建模块303,具体用于:
获取每个约束原子对对应的约束序号k,其中,k∈[1,N]。
获取每个约束原子对中两个关联耦合的原子k1和k2的相对矢量位置r和矢量速度
Figure BDA0003220248380000131
针对每一约束原子对、基于单步迭代时间Δt、原子k1和k2的相对矢量位置
Figure BDA0003220248380000141
和矢量速度
Figure BDA0003220248380000142
构建N阶矢量c,基于每一约束原子对组内每一约束原子对的约束条件构建N×N矩阵,得到矩阵A
矢量c的第k阶的元素
Figure BDA0003220248380000143
所述矩阵A内的元素表示为:
Figure BDA0003220248380000144
其中,k和k′为约束原子对组内的两个约束序号,k1和k2为约束序号k约束内的两个原子,k1′和k2′为k′约束内的两个原子,δ为狄拉克函数。
所述约束力计算模块304,具体用于:
确定每个约束原子对组的N阶矢量c;针对每个约束原子对组,基于对应的N阶矢量c和所述N×N矩阵的逆矩阵的乘积确定组内所有约束原子对对应的N阶矢量约束力集合B,其中,B=A -1c,A -1为所述N×N矩阵的逆矩阵。
所述在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算的处理模块305,具体用于:
在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到对应的原子速度修正量;基于所述原子速度修正量修正本次计算中的待约束原子的原子速度,得到本次计算后的待约束原子的原子速度和约束精度,在迭代未结束前,使本次计算后的待约束原子的原子速度参与下一次计算;
其中,构成约束原子对的两个原子k1和k2各自对应的原子速度修正量为:
Figure BDA0003220248380000145
Figure BDA0003220248380000146
其中,所述约束原子对的约束力
Figure BDA0003220248380000147
Figure BDA0003220248380000148
为矩阵A逆矩阵的第k行,Bk表示约束序号为k的约束原子对间的约束力,也是所有约束原子对对应的N阶矢量约束力集合B的第k个元素,Δt为单步迭代的时间,mk1为原子k1的原子质量,mk2为原子k2的原子质量。
本次迭代后约束原子的速度为:
Figure BDA0003220248380000151
Figure BDA0003220248380000152
Figure BDA0003220248380000153
为本次计算后约束原子对k的约束精度。
基于本发明实施例公开的约束原子速度的装置,在本方案中,构建N×N的矩阵,即矩阵A,其中各元素的值在迭代过程中不会改变,因此不需要重新计算矩阵求逆,同时,在迭代计算过程中同时考虑相互有关联的约束原子对的相互影响,统一在矩阵中计算相互关联约束原子对内原子的作用,从而可以大大降低约束迭代次数,且通过快速收敛每个时刻刚性整体内原子速度,并将计算得到的速度用于下一时刻的原子运动,使得两两原子之间运动的相对速度与二者的相对位置正交,从而实现提高刚性模型约束原子速度的计算效率的目的。
对所公开的实施例的上述说明,使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的,本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下,在其它实施例中实现。因此,本发明将不会被限制于本文所示的这些实施例,而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (8)

1.一种约束原子速度的方法,其特征在于,所述方法包括:
确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成;
按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列;
利用所述每个约束原子对组基于相关约束条件构建N×N矩阵,以及构建N阶矢量,N的取值为约束原子对组内约束原子对的个数;
基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内两个原子对应的约束力;
在每个约束原子下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到约束精度;
判断所述约束精度是否小于预设约束精度,若否,则返回执行构建N阶矢量这一步骤进行迭代计算,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。
2.根据权利要求1所述的方法,其特征在于,所述利用所述约束原子对组基于相关约束条件构建N×N矩阵,以及构建N阶矢量,包括:
获取每个约束原子对对应的约束序号k,k∈[1,N];
获取每个约束原子对中两个关联耦合的原子k1和k2的相对矢量位置
Figure FDA0003220248370000015
以及原子速度
Figure FDA0003220248370000011
针对每一约束原子对、基于单步迭代时间Δt、原子k1和k2的相对矢量位置
Figure FDA0003220248370000012
和原子速度
Figure FDA0003220248370000013
构建N阶矢量c,其中,第k阶元素为
Figure FDA0003220248370000014
基于每一约束原子对组对应的N阶矢量和约束条件构建N×N矩阵,得到矩阵A
其中,所述矩阵A内的元素表示为:
Figure FDA0003220248370000021
其中,k和k′为约束原子对组内两个约束序号,k1和k2为约束序号k约束内的两个原子,k1′和k2′为k′约束内的两个原子,δ为狄拉克函数。
3.根据权利要求1所述的方法,其特征在于,所述基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内两个原子对应的约束力,包括:
确定每个约束原子对组的N阶矢量c;
针对每个约束原子对组,基于对应的N阶矢量c和所述N×N矩阵的逆矩阵的乘积确定组内所有约束原子对对应的N阶矢量约束力集合B,其中,B=A -1c,A -1为所述N×N矩阵的逆矩阵。
4.根据权利要求1所述的方法,其特征在于,所述在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,包括:
在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到对应的原子速度修正量;
基于所述原子速度修正量修正本次计算中的待约束原子的原子速度,得到本次计算后的待约束原子的原子速度和约束精度,在迭代未结束前,使本次计算后的待约束原子的原子速度参与下一次计算;
其中,构成约束原子对的两个原子k1和k2各自对应的原子速度修正量为:
Figure FDA0003220248370000022
Figure FDA0003220248370000023
其中,所述约束原子对的约束力
Figure FDA0003220248370000024
Figure FDA0003220248370000025
为矩阵A逆矩阵的第k行,Bk表示约束序号为k的约束原子对间的约束力,Δt为单步迭代的时间,mk1为原子k1的原子质量,mk2为原子k2的原子质量;
本次迭代后约束原子的速度为:
Figure FDA0003220248370000026
Figure FDA0003220248370000027
Figure FDA0003220248370000031
为本次计算后约束原子对k的约束精度。
5.一种约束原子速度的装置,其特征在于,所述装置包括:
确定模块,用于确定刚性模型中的原子构型,基于所述原子构型确定约束原子对,每个约束原子对由两个关联耦合的原子构成;
分组模块,用于按照约束原子对所包含的原子划分各个约束原子对,得到约束原子对组,每一约束原子对组包含不同的原子,每一约束原子对组内的约束原子对对应一个约束序号,序号从1递增连续排列;
构建模块,用于基于所述每个约束原子对组内约束原子对的约束条件构建N×N矩阵,以及构建N阶矢量,N的取值为约束原子对组内约束原子对的个数;
约束力计算模块,用于基于所述N×N矩阵和N阶矢量计算组内每个约束原子对内的两个原子对应的约束力;
处理模块,用于在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到约束精度;判断所述约束精度是否小于预设约束精度,若否,则返回构建模块,使所述构建模块构建N阶矢量,并进行迭代计算,直至当前得到的原子速度小于所述预设约束精度,结束迭代,确定当前得到的原子速度为原子最终速度,使所述待约束原子基于所述原子最终速度进行运动。
6.根据权利要求5所述的装置,其特征在于,所述构建模块,具体用于:
获取每个约束原子对对应的约束序号k,其中,k∈[1,N];
获取每个约束原子对中两个关联耦合的原子k1和k2的相对矢量位置
Figure FDA0003220248370000032
和矢量速度
Figure FDA0003220248370000033
针对每一约束原子对、基于单步迭代时间Δt、原子k1和k2的相对矢量位置
Figure FDA0003220248370000034
和矢量速度
Figure FDA0003220248370000035
构建N阶矢量c,基于每一约束原子对组对应的N阶矢量和约束条件构建N×N矩阵,得到矩阵A
其中,矢量c的第k阶元素
Figure FDA0003220248370000041
所述矩阵A内的元素表示为:
Figure FDA0003220248370000042
其中,k和k′为约束原子对组内的两个约束序号,k1和k2为约束序号k约束内的两个原子,k1′和k2′为k′约束内的两个原子,δ为狄拉克函数。
7.根据权利要求5所述的装置,其特征在于,所述约束力计算模块,具体用于:
确定对应每个约束原子对组的N阶矢量c;针对每个约束原子对组,基于对应的N阶矢量c和所述N×N矩阵的逆矩阵的乘积确定组内所有约束原子对对应的N阶矢量约束力集合B,其中,B=A -1c,A -1为所述N×N矩阵的逆矩阵。
8.根据权利要求5所述的装置,其特征在于,所述在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行迭代计算的处理模块,具体用于:
在每个约束原子对下,基于所述约束力、单步迭代时间和待约束原子的原子质量进行计算,得到对应的原子速度修正量;基于所述原子速度修正量修正本次计算中的待约束原子的原子速度,得到本次计算后的待约束原子的原子速度和约束精度,在迭代未结束前,使本次计算后的待约束原子的原子速度参与下一次计算;
其中,构成约束原子对的两个原子k1和k2各自对应的原子速度修正量为:
Figure FDA0003220248370000043
Figure FDA0003220248370000044
其中,所述约束原子对的约束力
Figure FDA0003220248370000045
Figure FDA0003220248370000046
为矩阵A逆矩阵的第k行,Bk表示约束序号为k的约束原子对间的约束力,Δt为单步迭代的时间,mk1为原子k1的原子质量,mk2为原子k2的原子质量;
本次迭代后约束原子的速度为:
Figure FDA0003220248370000051
Figure FDA0003220248370000052
Figure FDA0003220248370000053
为本次计算后约束原子对k的约束精度。
CN202110955870.4A 2021-08-19 2021-08-19 一种约束原子速度的方法及装置 Active CN113593652B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110955870.4A CN113593652B (zh) 2021-08-19 2021-08-19 一种约束原子速度的方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110955870.4A CN113593652B (zh) 2021-08-19 2021-08-19 一种约束原子速度的方法及装置

Publications (2)

Publication Number Publication Date
CN113593652A true CN113593652A (zh) 2021-11-02
CN113593652B CN113593652B (zh) 2024-02-06

Family

ID=78238629

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110955870.4A Active CN113593652B (zh) 2021-08-19 2021-08-19 一种约束原子速度的方法及装置

Country Status (1)

Country Link
CN (1) CN113593652B (zh)

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0544525B1 (en) * 1991-11-26 1999-01-27 Fujitsu Limited An apparatus for generating a constraint condition in a molecular dynamics method
US6014449A (en) * 1998-02-20 2000-01-11 Board Of Trustees Operating Michigan State University Computer-implemented system for analyzing rigidity of substructures within a macromolecule
US10489808B1 (en) * 2014-11-19 2019-11-26 VCE IP Holding Company LLC Policies based on combinations of atomic constraints
CN110866877A (zh) * 2019-11-12 2020-03-06 Oppo广东移动通信有限公司 基于约束粒子群算法的色彩校正方法、装置及终端设备

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0544525B1 (en) * 1991-11-26 1999-01-27 Fujitsu Limited An apparatus for generating a constraint condition in a molecular dynamics method
US6014449A (en) * 1998-02-20 2000-01-11 Board Of Trustees Operating Michigan State University Computer-implemented system for analyzing rigidity of substructures within a macromolecule
US10489808B1 (en) * 2014-11-19 2019-11-26 VCE IP Holding Company LLC Policies based on combinations of atomic constraints
CN110866877A (zh) * 2019-11-12 2020-03-06 Oppo广东移动通信有限公司 基于约束粒子群算法的色彩校正方法、装置及终端设备

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
李乾;徐华;张越;程芳;: "基于约束的刚体碰撞响应仿真研究与应用", 计算机仿真, no. 03 *

Also Published As

Publication number Publication date
CN113593652B (zh) 2024-02-06

Similar Documents

Publication Publication Date Title
Young et al. A numerical technique for solving stiff ordinary differential equations associated with the chemical kinetics of reactive-flow problems
Amir et al. Approximate reanalysis in topology optimization
Legoll et al. A micro-macro parareal algorithm: application to singularly perturbed ordinary differential equations
Vallabhajosyula et al. Conservation analysis of large biochemical networks
JP2012521591A5 (zh)
CN109145386B (zh) 一种用于工程结构的地震碰撞反应分析方法
Ozaki Efficient low-order scaling method for large-scale electronic structure calculations with localized basis functions
CN105975795A (zh) 一种高精度的多关节串联机械臂运动学反解解法
CN115455745B (zh) 自适应频点采样的扫频方法、系统及相关设备
CN110008635B (zh) 利用Newmark精细积分结合法对弹塑性结构地震反应分析的方法
Wu Convergence analysis of some second-order parareal algorithms
CN116166915A (zh) 一种基于改进gmres的电力系统潮流计算方法及系统
Li et al. Hyperelastic finite deformation analysis with the unsymmetric finite element method containing homogeneous solutions of linear elasticity
CN113593652A (zh) 一种约束原子速度的方法及装置
CN112434451B (zh) 一种基于分块并行计算的有限元分析方法
CN110096798A (zh) 一种多状态有限元模型修正的方法
CN108536954A (zh) 一种基于交点间断伽辽金的高精度格子波尔兹曼方法
CN104567937B (zh) 一种高精度的x射线脉冲星信号模拟方法
Li et al. ARCHITECT: Arbitrary-precision constant-hardware iterative compute
Berger et al. An optimized reconfigurable system for computing the phylogenetic likelihood function on dna data
CN111950123A (zh) 一种陀螺仪误差系数曲线拟合预测方法及系统
Lotfi et al. Two new three and four parametric with memory methods for solving nonlinear‎ equations
CN107273629A (zh) 模拟复杂流体运动的求解方法
CN104809300B (zh) 脉冲子结构与有限元刚‑弹混合连接的动力学模拟方法
CN101404005A (zh) 基于可变预处理迭代获得非比例阻尼力学系统模态频响的方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
CB02 Change of applicant information

Address after: Room 1010, No. 2, Lane 288, Kangning Road, Jing'an District, Shanghai, 200040

Applicant after: Shanghai Silang Technology Co.,Ltd.

Address before: 100176 room 506-1, 5th floor, building 6, courtyard 10, KEGU 1st Street, Beijing Economic and Technological Development Zone, Daxing District, Beijing

Applicant before: Beijing Si Lang science and Technology Co.,Ltd.

CB02 Change of applicant information
GR01 Patent grant
GR01 Patent grant