高通量计算筛选的材料设计接口应用软件①

冯岳璐,王 浩,杨小渝,靳 凯,万 萌

1(中国科学院 计算机网络信息中心,北京 100190)

2(中国科学院大学 计算机科学与技术学院,北京 100049)

3(南京航空航天大学 材料科学与技术学院,南京 210016)

复合材料是由两种或两种以上不同性质的材料,通过物理或化学的方法,在宏观上组成具有新性能的材料.高分子基复合材料由基体、增强体和界面3 部分组成,相比于传统材料,由于其具有独特优异的性能,例如:轻质、高强度、高模量,被广泛用于航空航天等领域.研制高分子基复合材料时,基体扮演着重要角色,组成基体的不同性质材料的配比如何确定,及不同配方的基体对复合材料的宏观性能又有哪些影响,一直是一个十分具有挑战性的难题[1-5].基于传统实验“试错法”,成本高,周期长,而计算模拟、大数据和人工智能的方法,可以帮助更快地进行材料的配方设计并预测材料的性质.

材料计算的主流软件如Material Studio[6,7]虽然提供了不同尺度的计算模块(CASTEP、Forcite Plus、Mesocite 等),并且该软件也支持同一尺度下多种计算任务(如Forcite Plus 模块支持结构优化、模拟退火、动力学平衡等).然而使用起来非常不方便,不同尺度及不同任务间的计算需要每步单独开展.尤其是,计算结果数据得不到有效保存和管理,极易丢失,非常不便于人工智能时代的材料设计和材料数据挖掘.

高通量计算和筛选方法可以生成大量参考的配方,自动计算,自动保存计算结果,根据筛选条件选择出符合要求的配方.可以快速排除许多不合理的配方,减少了材料研发周期和材料开发中的代价.通过计算模拟进行复合材料的配方设计和性质预测,往往涉及多个尺度的计算(如微观-介观-宏观)[8].即便在同一尺度下,还涉及不同的计算任务(如结构优化-模拟退火-动力学平衡等).因此,使用高通量计算和筛选方法进行复合材料的配方设计,应该是这样的一种思路:顺着“结构-分子运动-物理性质” 研究路线,如何高通量的产生大量不同的复合材料配方,并在这些配方的结构中通过一系列动力学模拟后,计算每一种配方下材料的结构性质和力学性质,并按照一定的筛选条件选择出合理的配方.通过筛选出的结构,分析其结构性质与物理性质的联系,得出哪些结构性质是决定其物理性质的关键因素.一般考虑运动单元(链节运动、链段运动、侧基运动、支链运动、晶区运动以及整个分子链运动等) 的变化以及运动方式(键长、键角的变化,侧基、直链、链节的旋转和摇摆运动) 的变化等因素.采用这种方法和模式进行材料配方设计和性质预测,涉及大规模计算,大量的数据处理,以及同一尺度下的不同计算任务之间的计算结果如何传递,和不同尺度之间的计算结果转换和传递方法.

因此,我们需要开发这样一种基于材料计算软件Materials Studio 的支持高通量自动流程计算和筛选的材料设计接口应用软件,它能自动生成大量的复合材料配方,根据用户需要添加计算模块、设置计算参数,计算作业以及参数设置完成之后提交作业;计算环境包含计算作业的调度和材料计算软件,并能根据用户此次的筛选需求自动存储相应的计算结果,利用统计学方法以及机器学习的方法对数据进行二次处理,从而得到相应的用于结构筛选所需的筛选因子数据,并根据筛选因子筛选出符合条件的配方.

1 需求分析和软件架构设计

1.1 软件需求分析

本着“一键式、自动流程、高通量筛选” 等特点,该软件使用了客户端/服务器模式的架构.客户端可以设置计算作业、设置计算参数、监控作业状态等;服务器端负责在用户提交计算作业并开始计算之后,负责作业间的调度与计算结果的保存于分析,并将其存入数据库中.提供计算资源的Materials Studio 及其计算模块可以是任意计算节点,小到PC 机,大到计算中心.该软件的主要特点如下:

(1)自行设计每次计算机作业流程、添加计算模块.

(2)高通量的生成材料配方.

(3)按需设置各个计算模块中的参数.

(4)自动进行作业的调度.

(5)自动提取计算结果,使用统计学和机器学习方法分析数据.

(6)可下载所有的计算结果文件或下载经过筛选的符合条件的配方文件.

1.2 软件架构设计

该软件主要分为流程设计、流程监控、作业调度、计算结果分析、结构筛选5 个部分.流程设计、作业调度及结构筛选3 部分是整个系统的核心,主要包含了高通量生成计算任务,高通量自动流程的实现,以及高通量结构筛选,如图1所示.

图1 软件架构图

流程设计及流程监控部分中,可自行添加和删除计算作业、在每个作业中可根据用户需要添加计算模块、在每个计算模块中也可根据计算需求进行参数设置、用户添加完计算模块及设置完计算参数后可以提交作业.之后,用户可对计算作业的运行状态进行监控.

作业调度部分包含作业的自动调度及计算节点.作业的自动调度中,软件服务器根据用户设置的计算参数和计算模块,将每个计算模块中的计算任务发送至计算节点,计算节点计算完成之后,将计算结果发送至提交计算任务的服务器.

计算结果分析部分从计算节点计算完成并发送至提交计算任务服务器的计算结果中,提取本次计算所需的数据,直接存储或者进行处理后存储,生成筛选因子.

结构筛选部分中,用户可以根据添加的筛选因子进行结构筛选,对满足筛选条件的结构,用户可以根据需要选择其中的某些配方,并下载这个配方下的所有计算结果文件,或者直接下载所有配方的计算结果文件.

1.3 流程设计

本节主要介绍使用本软件进行高通量微观建模以及微观尺度的分子动力学模拟两部分的计算流程,具体的流程介绍如图2所示.

在微观尺度的模拟中,首先需要产生大量的候选配方并建立模型,系统提供两种产生候选模型的方法,分别是用户上传和高通量生成大量候选配方.模型上传或高通量建模完成之后需要对每个配方进行分子动力学模拟、结构交联、结构分析、力学性质计算和结构筛选,用户可以自行选择本次作业中需要添加的计算模块.

在进行流程设计时,用户首先需要选择使用高通量建模或上传已经建立好的模型,有了微观尺度的计算模型,就需要对结构进行初步优化.优化完成之后可以选择是否进行交联反应建立交联模型,如果选择了交联反应模块,则在交联模型建立完成之后对结构再次进行动力学模拟,最后进行结构分析、力学性质计算和筛选.不选择交联反应模块则直接进行结构分析、力学性质计算和筛选.

2 接口应用软件设计

软件基于客户端-服务器架构,主要基于Materials Studio 的核心模块进行开发.不同尺度的材料计算模块(如Amorphous cell、Forcite Plus、Mesocite、CASTEP等) 运行于计算节点.服务器可以高通量生成作业,进行参数设置,以及作业提交等.将每次计算任务发送至计算节点,一旦计算结束,计算节点会将计算结果返回至服务器端进行分析、筛选以及可视化呈现.

图2 微观尺度模拟计算流程图

目前该软件主要支持Materials Studio 的Amorphous Cell 模块和Forcite Plus 模块,通过Amorphous Cell 模块生成微观尺度的复合材料分子模型,使用Forcite Plus 模块中的不同计算任务对模型进行分子动力学模拟.自动流程是整个软件的核心,自动流程的核心算法实现参照了MatCloud[9-13].

2.1 Materials Studio 核心模块简介

进行性质预测首先需要建立高分子模型,Amorphous Cell 就是一个高分子建模模块,能够进行塑料、玻璃、树脂等无定型聚合物的建模.Forcite Plus 是Materials Studio 提供的一个分子动力学计算模块,它包含各种通用力场和支持5 种任务,每个任务可以是一个单独的模块.

实际应用中,用户可以根据自己的需求选择不同的任务组合.可以选择上传建立好的模型或上传单分子模型,根据环境参数进行无定形建模.在分子动力学模拟中可以选择多种任务,比如几何优化-模拟退火-动力学平衡、几何优化-动力学平衡等多种流程.

在进行上述操作时,每个计算任务均需要单独地完成,并将计算结果存储下来.当前一个任务(如几何优化) 完成后,需要将其计算结果读取,再开始下一个任务的计算(如模拟退火),非常不便捷.

2.2 软件核心功能

为了解决Material Studio 多个计算任务不够便捷的问题,我们开发了一个基于Material Studio 多个计算模块的高通量材料计算接口应用软件[14-25],它能够高通量生成作业,能够进行不同任务间的协调和自动流程,并将作业发送至计算节点开展计算,基于计算结果进行筛选.我们以环氧树脂基复合材料的配方设计,进行该软件的应用验证.

目前该软件的1.0 版,主要支持微观尺度的计算,支持Material Studio 的Amorphous Cell 和Forcite Plus 模块.软件主要有以下的核心模块:(1) 无定形建模模块(Amorphous Cell);(2) 上传模型模块(Upload-Model);(3) 动力学平衡模块(Forcite Plus Dynamic);(4) 几何优化模块(Forcite Plus Geometry Optimization);(5) 模拟退火模块(Forcite Plus Anneal);(6) 结构分析模块(Forcite Plus Analysis);(7) 力学性质计算模块(Forcite Plus Mechanical Properties);(8) 结构交联模块(XLink);(9) 结构筛选模块(Pickout).

2.2.1 按不同配比,高通量方式生成大量的分子结构

无定形建模模块:支持用户上传单分子模型,可以根据用户选择的配比,通过接口调用服务器端的Amorphous Cell 模块,生成大量的复合材料分子模型.这些输出的高分子模型,就是Forcite Plus 各个计算模块的输入文件.

高分子基复合材料的配方设计,涉及组成基体的组分种类及其配比.例如,以环氧树脂基复合材料的基体为例,包含3 种组分,分别是:环氧树脂、固化剂、增韧剂.而环氧树脂组分可以有两种分子结构,分别是E54 环氧树脂分子,AG80 环氧树脂分子,固化剂又可以有双氰氨,二氨基二苯酚等.究竟在哪个配方下,才能使这种复合材料有着更优的性能,就是配方设计问题.

确定基体中不同组分的分子结构数之比,如表1.

表1 基体中不同组分的分子结构数比例

确定环氧树脂组分中不同的分子的分子结构数,如表2所示(此例中选择两种环氧树脂分子结构).

表2 环氧树脂组分中不同分子的分子结构数(单位:个)

从表2中任意选取一个组合作为环氧树脂组分的配方.比如,本例中选择50 个E54 分子、50 个AG80分子,组成100 个环氧树脂组分配方.选择双氰胺、二氨基二苯酚、甲基硫氢苯酐3 种固化剂分子,选择PES、PEI、PSF 这3 种增韧剂分子.则按表1的基体中不同组分的分子结构比例计算,则会产生1326 种固化剂配方,1326 种增韧剂配方,则本例中共有1758 276(1326×1326) 个基体配方.如表3和表4所示为固化剂和增韧剂组分的配方.

表3 生成的固化剂组分配方(单位:个)

表4 生成的增韧剂组分配方(单位:个)

无定形建模模块的主要功能,就是完成搭建具有多种组分,如上面的例子中包含两种环氧树脂分子,3 种固化剂分子,3 种增韧剂分子,在不同配方下的高分子共混模型.每一个分子模型的构建都是Materials Studio 的Amorphous Cell 模块采用的蒙特卡洛方法进行无定形模块搭建,每个分子模型都是随机产生的,所以在本项目中的无定形建模模块对每一个配方只生成一个分子模型,每个配方的后续计算都是在这个分子模型的基础上进行的.

2.2.2 上传搭建好的模型

上传模型模块:在微观尺度的模拟中本软件支持用户上传两种类型的文件,一种是无定形建模模块搭建好的模型,一种是单分子模型,两种结构类型支持的是两种流程,本模块支持用户上传无定形建模模块搭建好的模型,可以直接作为Forcite Plus 的各个计算模块的输入文件,直接开展分子动力学模拟.

2.2.3 微观尺度的分子动力学模拟

分子动力学模拟的目标是对高分子结构中的每一个小分子进行优化,优化之后找出体系在变化过程中总能量最低的构型,再对其进行特定条件下的模拟运动,并分析其各种性质.包含以下几种计算任务.

几何优化模块:优化原子坐标和晶胞参数,支持原子笛卡尔坐标和晶胞参数的限定,可添加外应力.输入文件是原子或分子结构,输出是几何优化后的结构.

模拟退火模块:基于不同温度点的动力学模拟,实现体系的反复升降温过程,辅助扫描势能面,寻找最优的分子构象,吸附构象等.该模块的输入文件是分子结构,输出是从不同温度条件下的分子结构中找出的5 个能量最低的构象.

动力学平衡模块:基于牛顿运动方程,研究原子核在特定系综(NVE、NVT、NPT、NPH) 条件下的运动,并结合分析模块计算各种性质.本模块可以对模型在单一温度条件下进行动力学平衡,也可以对模型在一个温度范围内按一定步长进行多次的动力学平衡,用以计算其在温度变化过程中模型的各种性质变化.

2.2.4 结构交联反应

交联反应是指2 个或者更多的分子(一般为线型分子) 相互键合交联成网络结构的较稳定分子(体型分子) 的反应.反应分为物理交联和化学交联,分别是由氢键和共价键结合而成的.在本系统中我们使用的实验理论基于化学交联,如橡胶的硫化、不饱和聚酯通过链式聚合反应的固化、环氧树脂与固化剂的反应等.线形聚合物经过适度交联,在力学强度、弹性、尺寸稳定性、耐溶剂性或化学稳定性等方面均有改善,所以交联反应常被用于聚合物的改性.

结构交联模块:本软件中的交联模块是在对模型进行初步的分子动力学模拟之后,找出每个模型的总能量较低的构型,需要对其进行交联,产生组合键,该模块中也可以设置交联度以及其他参数.初步的分子动力学模拟之后的优化结构文件作为输入,输出为多个分子相互键合交联之后的较稳定分子.

2.2.5 结构分析和力学性质计算

结构分析/力学性质计算模块:这两个模块主要对结构进行力学性质计算和结构分析,并将所需的计算结果存入数据库.

在结构分析模块中,用户可以自行上传计算和分析脚本,也可以从提供的16 种性质中选择一些作为本次计算的分析对象.默认的16 种性质包含与结构相关的性质(径向分布函数等),与动力学相关的性质(温度分布曲线等),与统计学相关的性质(密度分布曲线等)等.在基本性质计算完成之后可以对计算结果进行分析从而得出其他性质,以玻璃化转变温度为例,如果一个配方在不同的温度条件下的动力学平衡计算次数不满足条件(小于4 次),则不会计算玻璃化转变温度.如果不同温度条件下动力学平衡次数满足条件(大于等于4),则会计算如玻璃化转变温度等性质.其中玻璃化转变温度的计算以及计算结果存储如下:一个配方在不同温度条件下进行动力学平衡之后,会继续计算其密度,并在温度密度两个维度下分段拟合出两条最符合温度密度变化趋势的直线,两直线的交点则为该配方的玻璃化转变温度.将计算出的两条直线的斜率、截距还有交点处所对应的玻璃化转变温度预测值存入数据库.

力学性质模块和结构分析模块的计算对象相同,需要对每个动力学平衡计算完成之后的结构进行力学性质计算,力学性质计算完成之后需要提取相关的筛选因子,在本文中此模块选择的是杨氏模量作为筛选因子,提取完成之后将其存入数据库.

这两种模块的输入为分子动力学模拟计算完成之后的文件,输出为各种分析结果,包含表格、文本文档等.在设计计算流程时,可以选择这两个模块的任意一个或两个都选择.

2.2.6 结构筛选

结构筛选模块:该模块主要是展示所有的计算结果,也可以进行筛选、下载等功能.在本模块中,会根据用户在结构分析和力学性质计算两个模块中选择的需要计算、提取和存储的性质,将所有配方的计算结果以图表的形式展示出来.页面也有筛选条件设置框和筛选按钮,用户可以根据需要设置筛选参数,筛选出符合条件配方的计算结果进行展示,也可以进行下载.比如,在一个流程中,我们添加了两个计算模块,分别进行结构分析以及力学性质计算,结构筛选模块中读取以上两个模块的计算结果,先将所有的计算结果展示出来,之后用户可以自己设置筛选参数,筛选出满足条件的配方.

3 应用案例

下面给出一个例子,来说明本软件如何产生大量配方并进行计算和筛选的.案例描述:这里假设一种双酚A 型环氧树脂分子(DGEBA),一种固化剂分子(DDS)作为基础组分(如图3所示),按不同的比例产生配方.因此我们需要确定在哪种配方下该种树脂基复合材料的玻璃化转变温度在指定区间内(如1100 K到1200 K).单个模型进行初步优化之后的计算基本步骤为,首先建立交联模型,对模型进行动力学平衡,之后进行结构分析和力学性质计算,最后对配方进行筛选.

图3 环氧树脂基体不同组分的分子模型

3.1 流程设计

生成不同配方的环氧树脂基复合材料模型,并计算不同配方下结构的杨氏模量和玻璃化转变温度.然后基于杨氏模量和玻璃化转变温度进行筛选,从而筛选出符合条件的配方.

本次计算在首次动力学平衡之后,对每个模型中的大量分子进行交联.因为要计算不同配方的玻璃化转变温度,所以需要对每个配方下的模型在不同的温度条件下进行动力学平衡.可以通过设置起始温度和终止温度,并设置温度变化的步长,实现每个配方的模型在不同温度条件下进行动力学平衡计算.添加力学性质计算模块主要是为了计算筛选结构所需的杨氏模量.具体的计算流程如图4所示.

图4 微观尺度的计算和筛选流程

3.2 参数设置

本例中,采用两种分子混合,分别是一种环氧树脂分子和一种固化剂分子(如图3),因为要进行交联反应,基于这两种分子的特性,对其配比为大于1:1 即可,本例中固化剂分子采用5 至20 个结构分别进行混合,共16 种配方.如图5所示.

图5 配方设计

在交联任务前,需要对每种配方的结构在单个温度条件下进行动力学平衡,故将本次计算的起始、终止温度设置为一样,温度步长大于0 即可.如图6所示.

图6 单一温度下进行动力学平衡

在交联计算中,设置本次计算的交联度为50.如图7中参数Conversion 所示.

图7 交联反应参数设置

在第二次动力学平衡中,将起始温度设置为298 K,终止温度设置为1900 K,温度的步长为100 K,系综为NPT.温度设置如图8所示.

图8 同一模型进行多次动力学平衡温度设置

在分析模块中,本次的计算中采用温度密度曲线来预测玻璃化转变温度,所以在本任务中选择密度作为本次分析的性质.如图9所示,选中的性质颜色会变化.

图9 选择需要进行分析的性质

3.3 计算结果

计算流程设计和参数配置完成之后,点击“Start”键,启动该作业.整个作业会自动化运行,不需人工的干预.计算结束后,可进入结构筛选模块的View 界面,展示所有配方的计算结果,分别是配方的分子结构数,摩尔质量比,该配方的玻璃化转变温度以及最接近这个温度的两个温度下的杨氏模量.

可见,共有16 种候选样本空间;无定形建模,共生成16 种复合材料配方分子的模型;第一次几何优化,输出16 种复合材料配方分子模型优化后模型;第一次动力学平衡,对16 个配方的模型进行动力学平衡;交联反应,输出为16 个配方交联后的模型;第二次几何优化,对16 个交联后的模型进行优化;第二次动力学平衡,16 个优化后的模型中每个结构都需要在不同的温度条件下进行动力学平衡,共16×17 次动力学平衡计算任务,同样也有16×17 次的结构分析和力学性质计算任务.从建模到获取杨氏模量,总共运行了2560 分钟 (42.7 h),平均每个配方计算时间是160 分钟(2.7 h).

如图10所示,由于计算结果太多(16 个配方的计算结果),此处只展示一部分(前两个配方).

图10 计算结果

3.4 筛选结果

本次的筛选条件为玻璃化转变温度在1100 K 到1200 K 之间的所有配方.筛选结果如图11所示,筛选出的结果5_20 的配方,该配方的玻璃化转变温度为1152 K.

在所有计算结果和根据条件筛选出的结果展示界面中,用户都可以选择下载某一配方或所有配方下的计算结果文件.比如本次筛选出的符合条件的配方只有5_20,那么只能下载这一种配方的计算结果文件,如果根据筛选条件筛选出的配方有两个,那么用户可以选择下载这两个配方的其中一种或者两种都下载.

4 讨论

4.1 创新点

通过试验方法开展复合材料的设计,往往周期较长,效率不高.不使用本软件,仅使用Material Studio,用户首先需要建立庞大的文件系统,每次计算任务都需要单独设置配方和计算参数,在计算完成之后,只能手动提取需要分析的计算结果并保存,还需要有机器学习领域知识并会使用相关方法进行数据拟合,拟合结果记录、筛选等都需要耗费大量的人力和时间.

使用本软件,通过高通量筛选的办法,开展复合材料设计,在如下的关键点上进行了创新.

(1) 供筛选样本空间的自动生成.本软件使用递归算法生成了大量的配方(以本文第2.2.1 节为例,生成了超过106种配方),用户可以在预生成的大量配方中选择保留哪些配方,删除哪些配方(以本文第2.2.1 节为例,用户可以在候选配方中选择性删除一部分,如果环氧树脂组分对某一种固化剂敏感度较低时,可以将复合材料的各个配方种该固化剂分子的数量间隔增大一些).高通量方式快速生成复合材料配方的候选空间,极大地提高了计算效率,降低了实验的盲目性,缩短复合材料的研发周期.

(2) 多个步骤的自动化、流水式计算(如无定形建模-几何优化-模拟退火-动力学平衡-结构交联-力学性质计算-结构分析-筛选) 通过工作流一键式的完成.如果没有该软件,操作的流程是:首先需要进行建模,建模完毕,用户需要先打开几何优化模块,对结构进行优化得到优化后的输出文件,然后重新打开模拟退火模块,加载上一步几何优化得到的输出文件,再进行模拟退火.以此类推,模拟退火,动力学平衡的计算,也是同样的操作模式,非常不便捷.而计算模块的流程化设计,实现了数据的自动流动,避免了上述逐步的操作流程,提高复合材料计算仿真的灵活性和可操作性.

(3) 材料关键性数据的自动获取并存储,实现材料数据的重复利用和可回收.一旦计算结束,数据会抽取和自动归档.如果没有该软件,这些物性数据(计算参数如温度,动力学计算中的密度、径向分布函数,力学性质如杨氏模量等) 是隐藏在各个模块的输入输出文件中.而本软件,通过与结构算法的结合,为多尺度计算到组织性能-工艺和设计的相互连接打下良好的基础.

(4) 集成人工智能方法,自动进行筛选.根据关键性质使用机器学习方法生成配方的筛选因子.以玻璃化转变温度为例,本软件采用二分单变量线性回归的方法对数据进行了分段拟合,找出了各个配方的玻璃化转变温度.

(5) 配方的计算环节中文件归档和管理,各个配方下计算结果文件的下载.根据筛选结果,用户可以选择性下载部分满足条件配方的所有计算结果文件.用户可以在任何地方访问服务器,下载计算结果文件并开展材料研究.

图11 筛选结果

4.2 下一步的工作及难点

下一步我们会继续开发介观尺度的耗散粒子动力学模块Mesocite 的接口应用,在筛选出微观尺度的配方结构之后,需要对其进行介观建模,而如何对微观复合材料分子结构中进行粗粒化,以及建立好粗粒化分子之后填充介观模型,有很多种方式,而从中选择最能体现材料本质的模型也是一个难题.

5 结束语

航空复合材料作为21世纪主流材料,采用高效设计方法研发新型复合材料,减少复合材料研发成本依然是未来急需解决的问题.高通量材料计算和多尺度模拟方法和理念,尤其适合于复合材料配方的理论设计.为此,我们研发了一个基于Materials Studio 的支持高通量计算和筛选的材料设计接口应用软件.该软件具有“一键式、自动流程、高通量筛选” 等特点,目前该软件仅支持微观尺度的高通量建模和分子动力学模块的自动流程筛选,下一步我们将实现基于微观尺度的分子动力学和介观尺度的耗散粒子动力学的跨尺度计算模拟和高通量自动流程筛选,开展环氧树脂基复合材料的配方设计.