取代基和溶剂对ZnTPP-SubPc供体-受体对中能量转移与电荷转移竞争的影响:一项非绝热动力学研究

《Journal of Photochemistry and Photobiology A: Chemistry》:Substituent and solvent effects on the competition between energy transfer and charge transfer in ZnTPP-SubPc donor-acceptor dyads: A nonadiabatic dynamics study

【字体: 时间:2026年07月19日 来源:Journal of Photochemistry and Photobiology A: Chemistry 4.4

编辑推荐:

  •通过电子模拟与非绝热模拟揭示ZnTPP-SubPc激发态及超快电荷转移动力学。•这些模拟结果与实验现象一致,同时阐明了电荷转移机制。•氟取代显著稳定了低能量电荷转移状态,促进了溶液中高效的电荷分离。 1. 引言 在天然光活性系统中,太阳能向可用形式的转化依赖于一系列连续且相互

  •通过电子模拟与非绝热模拟揭示ZnTPP-SubPc激发态及超快电荷转移动力学。•这些模拟结果与实验现象一致,同时阐明了电荷转移机制。•氟取代显著稳定了低能量电荷转移状态,促进了溶液中高效的电荷分离。

1. 引言
在天然光活性系统中,太阳能向可用形式的转化依赖于一系列连续且相互关联的基本过程,包括光吸收、激发能传递以及光诱导的电子转移[1]、[2]、[3]、[4]。深入理解这些过程的协同调控及其内在联系,为构建高效人工光合作用系统以及开发新型光功能分子提供了重要的科学基础[5]。尤其是在有机供体-受体系统中,光激发往往会导致多种竞争性的激发态过程,包括局部激发、激发能传递、电子转移和电荷分离,这些过程的相对贡献直接决定了系统的最终状态及能量利用方式[6]、[7]、[8]、[9]。因此,系统地阐明D-A系统中激发态电子结构与超快动力学之间的关系,对于开发高性能光电及光化学功能材料具有重要意义。

在众多有机D-A分子中,卟啉及其衍生物因其出色的热稳定性、可调节的电子结构以及良好的光电性能,被广泛用于光伏电池、染料敏化太阳能电池以及人工光合作用领域[10]、[11]、[12]。尤其是锌卟啉衍生物,由于其较长的激发态寿命和较强的电子给予能力,已成为构建光功能D-A分子系统的重要构件[13]、[14]。然而,基于卟啉的分子也存在固有局限性。例如,它们的Q带吸收位于可见光谱的红光区域,通常具有较弱的振子强度,从而限制了对太阳光谱的高效利用[15]。为克服这一局限,人们将卟啉与在可见光区域具有更强吸收能力的“天线”色素结合,构建能够实现分子内能量传递或电荷分离的共轭体系,以此作为提高光捕获效率并控制激发态演化的重要策略[16]、[17]。已有大量研究致力于基于卟啉分子的D-A系统。例如,Flamigni等人报道了在自组装结构中,自由碱基卟啉客体向锌卟啉宿主进行能量传递,表明锌卟啉可作为高效的激发能接受体[18]。Lazarides和Leonardi等人的后续研究进一步表明,在BODIPY-卟啉二聚体以及1,2,3-三唑连接的BODIPY-Zn(II)介孔四苯卟啉阵列中,BODIPY色素可作为高效供体,实现快速且高效的福斯特共振能量传递至卟啉单元[16]、[19]。Kelber等人通过对供体-卟啉-受体三元体系的实验研究,证明含有卟啉核的D-A结构能够产生长寿命的光诱导电荷分离态[20]。

另一种与卟啉密切相关且日益受到关注的光功能单元是亚酞菁(SubPc)。作为酞菁的缩合类似物,SubPc具有独特的非平面锥形14π芳香框架,同时还具备强可见光吸收能力、高激发态能量以及显著的电子亲和力。通过引入强吸电子基团,如氟原子,其电子接受能力还可进一步增强[21]、[22]、[23]。此外,SubPc是作为供体还是受体取决于具体的分子环境,这体现了其作用的强烈依赖性。因此,其光诱导行为不能简单地类比于传统的平面受体分子[24]、[25]。先前的研究表明,当SubPc与不同的色素或受体单元结合时,可表现出丰富的能量传递和电荷转移行为[25]。例如,González-Rodríguez等人早在SubPc-C60二聚体中就发现,通过取代基效应和环境因素,可以调节SubPc的电子给予能力,从而在适当条件下形成分子内电荷分离态。在随后针对SubPc-二茂铁二聚体的研究中,同一团队又发现了光诱导电荷转移态,进一步表明在特定分子环境中,SubPc也可作为电子受体,从二茂铁那里接收电子[26]、[27]。El-Khouly等人进一步证实,在自组装的ZnNc-SubPc二聚体中,存在从锌萘酞菁到SubPc的分子内电子转移,这说明SubPc在大环共轭体系中也可作为高效受体,参与电荷分离。他们还在人工光合作用的模型系统中表明,SubPc的行为取决于其与耦合片段之间的能量匹配情况,既可以作为电子供体,也可以作为电子受体[24]、[28]。与此同时,KC等人的两项连续研究显示,通过分子修饰,可进一步调节SubPc的电子给予能力及其光诱导电荷分离行为:外围硫芳基取代可显著增强其供体特性;而在含有二级供体的SubPc-富勒烯共轭物中,改变连接二级供体与SubPc的桥键,可进一步调控电子传递及随后的空穴传递过程[29]、[30]。综上所述,这些研究结果表明,SubPc在D-A系统中并非只扮演单一角色,而是可以根据耦合伴侣以及能量级匹配情况,表现为电子供体、电子受体或中间传递体,从而产生高度可调控的光诱导电荷转移行为。

尽管已有大量实验和理论研究探讨了卟啉分子及SubPc与其他色素结合时的激发态过程,但对于它们直接结合形成的共轭体系,尤其是关于低能级激发态的排序、局部激发与电荷转移状态之间的竞争,以及外部环境如何调控超快动力学等方面,仍缺乏系统而连贯的理论理解[29]。先前的理论研究表明,D-A系统中的电荷转移过程不仅受供体和受体片段内在电子性质的影响,还与分子构型、取代基修饰以及溶剂环境有关[31]、[32]。例如,Ghosh等人使用密度泛函理论及受限密度泛函理论,研究了锌四苯卟啉/C70复合物中的电荷转移态[33]。Basurto等人则对内嵌富勒烯-锌四苯卟啉二聚体的电子结构及电荷转移态进行了密度泛函理论分析,发现受体种类的变化可显著影响电荷转移的激发能量[34]。对于基于SubPc的系统,Farhat等人通过密度泛函理论及时间依赖性自洽场计算表明,通过修饰受体端基,可有效调节SubPc衍生物的电子和光学特性以及电荷转移能力[35]。此外,Tang等人将受限全线性响应时间依赖性密度泛函理论与非绝热分子动力学模拟相结合,发现SubPcH在PTZ-SubPcH和PTZ-SubPcF体系中表现为电子受体,而在PDI-SubPcH体系中,SubPcH则表现为电子供体[36]。总体而言,这些研究表明,分子结构修饰和环境因素可显著改变D-A系统中局部激发与电荷转移状态的能量排序及其后续的动力学行为。

在此背景下,由卟啉与SubPc单元通过桥键连接而成的直接耦合D-A共轭体系,因其兼具优良的光物理性质和可调节的电子结构,为研究能量传递与电荷传递过程之间的相互作用提供了理想平台。Bressan等人使用飞秒瞬态吸收光谱和二维电子光谱技术,研究了在甲苯溶液中的锌四苯卟啉-亚酞菁二聚体(ZnTPP-SubPcH)的激发态动力学,其中SubPcH为供体,ZnTPP为受体。实验表明,在560纳米波长激发下,由于SubPcH的吸收贡献更为显著,该体系主要表现为从SubPcH向ZnTPP的能量传递[37]。然而,相关激发态的排序、不同的松弛路径及其背后的微观机制仍有待阐明。为此,我们基于Bressan等人的研究,通过整合静态电子结构计算与基于受限全线性响应时间依赖性密度泛函理论框架下的非绝热分子动力学模拟,对ZnTPP-SubPcH二聚体的光诱导动力学进行了全面的理论研究。为进一步探究取代基的影响,我们还在此基础上设计了氟代衍生物ZnTPP-SubPcF。通过这项研究,我们旨在阐明取代基修饰和溶剂环境如何共同影响局部激发与电荷转移状态的能量排序,以及能量传递与电荷转移在ZnTPP-SubPc体系中的竞争关系,从而深入理解卟啉-亚酞菁D-A体系中的激发态调控机制。

2. 方法
采用Perdew–Burke–Ernzerhof(PBE)函数结合Grimme的D3色散修正,对ZnTPP-SubPcH和ZnTPP-SubPcF的最小能量结构进行了几何优化。所有原子的计算均采用6-31G*基组[38]、[39]、[40]、[41]、[42]。由于两种体系都涉及电荷转移激发态,因此随后利用经优化的基态几何结构,通过结合经优化调整的筛选范围分离杂化函数(OT-SRSH函数)的线性响应时间依赖性密度泛函理论,计算了它们的激发态性质(详细内容见补充材料)。通过以下两步流程优化了两个关键参数:范围分离参数ω和长程极限参数β[43]。首先,根据电离势定理优化ω参数,得到ZnTPP-SubPcH和ZnTPP-SubPcF的ω值分别为0.1765和0.1049。在优化后的ω值基础上,接着在极化连续模型中,以甲苯和苯腈作为溶剂,优化β参数。由此得到的ZnTPP-SubPcH在甲苯中的β值为0.1241,在苯腈中的β值为-0.1921;而ZnTPP-SubPcF在甲苯中的β值为0.1871,在苯腈中的β值为-0.2442。在本研究的所有计算中,长程修正版PBE0函数中的α参数均采用默认值0.25。对于溶液中的体系,长程精确交换分数被设定为1/ε,其中ε表示溶剂的介电常数。据此,可通过公式β=1/ε?α来确定β值。对于介电常数为2.3741的甲苯,计算得到的β值为0.1712,与ZnTPP-SubPcH和ZnTPP-SubPcF的优化值0.1241和0.1871较为接近。同样,对于介电常数为25.592的苯腈,计算得到的β值为-0.2109,也与相应优化值较为接近。在优化了ω和β参数之后,结合PCM溶剂模型,使用OT-SRSH函数进行了线性响应时间依赖性密度泛函理论计算,以研究溶剂化体系中的激发态性质。包含这些优化参数的范围分离函数被称为LC-PBE0*,以区别于原有的LC-PBE0函数。这种方法在多项先前研究中已被证明具有可靠性[43]、[44]、[45]、[46]、[47]。

随后,利用GTSH程序中的NAMD模拟技术,基于线性响应时间依赖性密度泛函理论水平上的最少开关表面跳跃加经典路径近似方法,研究了ZnTPP-SubPcH和ZnTPP-SubPcF的激发态动力学[48]、[49]。所需的激发态能量、时间导数非绝热耦合系数等其他量,均采用OT-SRSH函数,即LC-PBE0*函数来计算。这两种体系中电子态之间的时间导数非绝热耦合系数的数值计算,则采用了我们团队最新开发的算法[49]、[50]、[51]、[52]。在最少开关表面跳跃加经典路径近似框架内,采用了Granucci等人提出的经验性能量依赖量子退相干修正,修正参数设定为0.1原子单位[53]、[54]、[55]。此外,NAMD计算是在经典路径近似下进行的,这种方法已被证明能够在没有显著的结构重排或键的形成与断裂情况下,有效描述激发态的松弛过程[49]、[50]、[51]、[56]、[57]、[58]。在经典路径近似中,仅需一个基态轨迹即可。因此,用于最少开关表面跳跃方法的态间能量和时间导数非绝热耦合系数只需计算一次。之后,可以沿着这个已计算出的基态轨迹,采样多个初始条件,这些初始条件代表表面跳跃轨迹的不同起始点。由于无需实时计算能量、梯度和时间导数非绝热耦合系数,这种近似方法大大提高了表面跳跃模拟的计算效率[56]、[59]、[60]、[61]、[62]。NAMD模拟的具体步骤如下:首先将ZnTPP-SubPcH和ZnTPP-SubPcF的最小能量结构加热到300开尔文,然后在该温度下让体系平衡约1皮秒,此过程中采用了链长为五的Nosé-Hoover链式恒温器[39]、[42]、[63]、[64]。之后进行2皮秒的非绝热动力学模拟,根据垂直跃迁能量及其对应的振子强度,随机生成1000个初始条件。所有从头算分子动力学模拟均采用PBE + D3水平进行,时间步长为1飞秒,所得到的NVE基态轨迹随后被用于FSSH-CPA模拟。对于每个初始条件,沿着相同的预定义核轨迹生成300条独立的FSSH轨迹,积分步长同样为1飞秒,总积分时间根据所研究过程的特征时间尺度来确定。最终统计结果是通过对所有1000×300条轨迹的平均值得出的。为在经典路径近似中实现详细平衡,采用了玻尔兹曼因子缩放的跃迁率来替代传统的跃迁拒绝机制与速度重标定方法:pkjt=pkjtbkjtEj>Ek#,其中bkj=exp?Ej?EkkBT。这一改进方法已在多种材料的从头算NAMD模拟中得到广泛应用[49]、[50]、[56]、[57]、[64]。几何结构优化及LR-TDDFT计算均使用GAUSSIAN 09软件完成[65]。ω和β参数则是通过Lu等人开发的经局部修改的optDFTω程序进行优化的[66]。AIMD模拟在CP2K软件包中开展,使用了QUICKSTEP模块,以及DZVP-MOLOPT-SR-GTH基组与Goedecker-Teter-Hutter(GTH)伪势,依然处于PBE + D3水平[67]、[68]、[69]、[70]、[71]、[72]。本研究中所有的NAMD模拟都是通过GTSH软件包完成的[48]。相应的电子密度和空穴密度则是用Multiwfn 3.8软件计算得出的[73]、[74]。更多理论细节可参见补充材料以及我们之前的研究。

3. 结果与讨论
ZnTPP-SubPcH和ZnTPP-SubPcF的最小能量基态结构是采用PBE + D3水平优化得到的。图1展示了这两种物质的骨架结构以及优化后的三维几何形态。在这两种结构中,SubPc中的中心硼原子都通过一个氧原子与ZnTPP直接相连。为进一步探究色散相互作用对分子几何结构的影响,这两种结构在相同的函数和基组设置下进行了再次优化,但此次未考虑D3色散修正,优化后的结构显示在图S1中。对比有无色散修正情况下的优化结构可以发现,两种体系在几何形态上存在明显差异,这说明色散相互作用对其平衡几何结构有着显著影响。此外,ZnTPP-SubPcH与ZnTPP-SubPcF在结构对比中的偏差程度并不相同,这表明取代基的不同会进一步影响分子内部片段间的非共价相互作用以及由此产生的构象变化。基于此,后续关于激发态性质与动态特性的分析均采用了经过PBE + D3优化后的结构作为计算基础。

下载:下载高分辨率图像(240KB)
下载:下载全尺寸图像
图1. 采用PBE + D3水平优化后的ZnTPP-SubPcH(上)与ZnTPP-SubPcF(下)的化学结构及三维几何形态。

基于优化后的结构,我们采用了经最优调谐的LC-PBE0*功能项的LR-TDDFT方法,计算了ZnTPP-SubPcH和ZnTPP-SubPcF在甲苯和苯腈溶剂中10个能量最低的单重激发态。为简化表述,两种结构中的氧原子都归入了ZnTPP片段。这样处理的依据是,作为连接桥的氧原子仅起到连接ZnTPP与SubPc的作用,对整体π共轭体系的贡献可以忽略不计。电子跃迁密度矩阵分析(见图S2)显示,跃迁密度主要分布在ZnTPP和SubPc的发色团上。对于ZnTPP-SubPcH,所有激发态下桥接氧原子上都没有明显的跃迁密度。而在ZnTPP-SubPcF中,桥接氧原子上仅有少量跃迁密度,占比约为0.1。其余的跃迁密度则主要集中在ZnTPP片段上。由于桥接氧原子的贡献相对较小,且不会改变跃迁密度的整体空间分布,因此在后续分析中也将该氧原子归入ZnTPP片段。如图1所示,蓝色区域代表ZnTPP片段,红色区域则代表SubPcH或SubPcF片段。基于激发态计算结果,我们再用Multiwfn 3.8软件分析了电子密度和空穴密度分布。ZnTPP-SubPcH与ZnTPP-SubPcF在甲苯和苯腈中的相应结果分别展示在图2和图3中。空穴密度以蓝色表示,电子密度以红色表示。我们首先研究了ZnTPP-SubPcH在甲苯和苯腈中的激发态特性(见图2)。从空间上的电子密度和空穴密度分布来看,甲苯中的10个激发态大致可分为四类。第一类包括S1和S2态,在这些状态下,电子和空穴都位于ZnTPP片段上,属于局限于ZnTPP上的局域激发态。第二类是S3和S4态,此时电子和空穴都集中在SubPcH片段上,属于局限于SubPcH上的局域激发态。第三类包含的状态数量最多,为S5至S8态,这类状态下电子位于SubPcH上,空穴位于ZnTPP上,属于从ZnTPP向SubPcH的电荷转移态。第四类是S9和S10态,此时电子主要在ZnTPP上,空穴主要在SubPcH上,属于从SubPcH向ZnTPP的电荷转移态。此外,这些激发态的振子强度与吸收峰强度密切相关。像S1(0.0543)、S2(0.0439)、S3(0.4250)和S4(0.5788)这样的局域激发态的振子强度,远远高于S5(0.0009)和S6(0.0002)这样的电荷转移态的振子强度。尤其值得注意的是,SubPcH的局域激发态的振子强度甚至高于ZnTPP的对应值。在苯腈中,10个激发态同样可以大致分为四类。第一类也是数量最多的一类,属于从ZnTPP向SubPcH的电荷转移态,包括S1、S2、S5和S6。第二类是S3和S4,属于局限于ZnTPP上的局域激发态。第三类是S7和S8,属于位于SubPcH上的局域激发态。第四类是S9和S10,属于从SubPcH向ZnTPP的电荷转移态。与之前类似,包括S3(0.0729)、S4(0.0665)、S7(0.4067)和S8(0.5498)在内的局域激发态的振子强度,仍远高于S1(0.0005)、S2(0.0000)、S5(0.0074)、S6(0.0037)、S9(0.0015)和S10(0.0025)等电荷转移态的振子强度。这些结果表明,局域激发态对应于吸收光谱中的亮特征,而电荷转移态则多为暗特征。

下载:下载高分辨率图像(517KB)
下载:下载全尺寸图像
图2. 在甲苯和苯腈的弗兰克-康登点上,采用TD-LC-PBE0*水平计算的ZnTPP-SubPcH的10个能量最低的单重激发态的电子(红色)和空穴(蓝色)密度分布,以及对应的垂直激发能。(关于图中颜色标注的含义,可参阅本文的网络版本。)

下载:下载高分辨率图像(800KB)
下载:下载全尺寸图像
图3. 在甲苯和苯腈的弗兰克-康登点上,采用TD-LC-PBE0*水平计算的ZnTPP-SubPcF的10个能量最低的单重激发态的电子(红色)和空穴(蓝色)密度分布,以及对应的垂直激发能。(关于图中颜色标注的含义,可参阅本文的网络版本。)

ZnTPP-SubPcF的激发态特性如图3所示。在甲苯和苯腈中,这10个激发态都可以大致分为三类。首先在甲苯中,第一类包括S1(0.0003)、S2(0.0006)、S3(0.0007)、S4(0.0006)、S9(0.0019)和S10(0.0015)态,这类状态下电子主要位于SubPcF片段上,空穴位于ZnTPP片段上,属于从ZnTPP向SubPcF的电荷转移态。第二类是S5(0.0725)和S6(0.0623)态,此时电子和空穴都主要位于ZnTPP部分上,属于ZnTPP单元的局域激发态。第三类是S7(0.5123)和S8(0.5066)态,属于位于SubPcF单元上的局域激发态。在苯腈中,第一类同样属于从ZnTPP向SubPcF的电荷转移态,主要包括S1(0.0001)、S2(0.0005)、S3(0.0005)、S4(0.0001)、S9(0.0016)和S10(0.0013)态。第二类仍是S5(0.0868)和S6(0.0935)态,属于局限于ZnTPP上的局域激发态。第三类是S7(0.5141)和S8(0.4804)态,属于位于SubPcF上的局域激发态。由于振子强度决定了吸收峰的强度,这些结果表明,那些位于SubPcF上的局域激发态是造成强吸收现象的主要原因。

对比ZnTPP-SubPcH与ZnTPP-SubPcF在甲苯和苯腈中的电子-空穴密度分布可以发现,两种体系之间最显著的差异在于低能量电荷转移态的稳定性,以及这些态对环境极性的反应方式。对于ZnTPP-SubPcH来说,在极性较弱的甲苯中,能量最低的S1至S4态主要是局域激发态,而更高能量的激发态则以电荷转移态为主,这说明含有氢取代基的体系依然具有明显的激发态特征,低能量区域以局域激发态为主,高能量区域则以电荷转移态为主。当溶剂极性升高到苯腈时,电荷转移态的稳定性显著增强,能量也相应降低,而局域激发态则处于中间能量区域。这表明ZnTPP-SubPcH的激发态排序对周围环境极为敏感,具有明显的溶剂依赖性。相比之下,ZnTPP-SubPcF的电荷转移态更倾向于向更低能量移动。在甲苯和苯腈中,低能量激发态的分布基本保持不变:S1至S4以及S9至S10态主要以电荷转移态为主,而S5至S8态则为局域激发态。这一结果进一步说明,即便在极性较弱的甲苯中,氟取代也能有效稳定最低能量区域的电荷转移态。这意味着氟取代显著增强了受体基团的吸电子能力,从而使电荷转移态在不同环境中的稳定性得到提升。与这一趋势一致的是,ZnTPP-SubPcF在甲苯中的S1态激发能分别为1.77电子伏特和1.56电子伏特,均低于ZnTPP-SubPcH对应的2.26电子伏特和1.94电子伏特。这些结果表明,氟取代会导致最低激发态及低能量激发态谱带的整体红移,这一效应在溶液环境中尤为明显。

基于激发态的垂直激发能及对应的振子强度,我们采用了Barbatti等人提出的半经典方法[75]、[76],并结合从基态动力学模拟中获取的1000个结构,计算了ZnTPP-SubPcH和ZnTPP-SubPcF在甲苯和苯腈中的紫外-可见吸收光谱。该方法能够考虑核热运动对吸收光谱的影响,图4展示了计算得到的光谱结果。此外,通过对激发态性质的片段级分析,还可以详细识别出在特定波长激发时所产生的不同类型激子的贡献。具体而言,局域激子被标记为SubPcHF?(蓝色)和ZnTPP?(绿色),而电荷转移激子则被标记为SubPcHF+ZnTPP?(黄色)和SubPcHF?ZnTPP+(橙色)。在甲苯中的ZnTPP-SubPcH体系中,可以观测到一条吸收带,其最大吸收波长为523纳米(2.37电子伏特)(见图4(a))。在该波长附近进行激发时,主要会产生局域激子,其中SubPcH?的贡献最为显著,其次是ZnTPP?,而两种电荷转移激子的贡献则相对较小。在苯腈中(见图4(b)),吸收带的最大吸收波长为532纳米(2.33电子伏特)。在该波长附近激发时,仍然主要以SubPcH?和ZnTPP?这两种局域激子为主,电荷转移激子的贡献依旧很小。与ZnTPP-SubPcH相比,ZnTPP-SubPcF在两种溶剂中的吸收峰都出现了更为明显的红移现象。在甲苯中(见图4(c)),ZnTPP-SubPcH的523纳米吸收带在ZnTPP-SubPcF中红移到了545纳米(2.27电子伏特)。在该波长附近激发时,主要产生的依然是局域激子,其中SubPcF?的贡献比ZnTPP?更强。在苯腈中(见图4(d)),532纳米的吸收带红移到了552纳米(2.25电子伏特)。在该波长附近激发时,产生的仍然是SubPcF?和ZnTPP?这两种局域激子。为了进一步验证我们理论方法的准确性,我们还在图4(a)中绘制了ZnTPP-SubPcH在甲苯中的实验吸收光谱(灰色虚线)。在470至650纳米的范围内,实验光谱也在560纳米处出现了一个峰值(2.21电子伏特),该值与计算得到的523纳米处的峰值(2.37电子伏特)仅相差约0.16电子伏特。图4. 使用LC-PBE0*函数,在OT-SRSH能级上,基于从基态分子动力学轨迹中抽取的1000个构型,模拟了ZnTPP-SubPcH在甲苯(a)和苯腈(b)中的吸收光谱,以及ZnTPP-SubPcF在甲苯(c)和苯腈(d)中的吸收光谱。图中还展示了相应的激子型贡献。虽然通过静态电子结构计算可以获取这两种体系激发态性质的重要信息,但其背后的光诱导动态过程仍不十分清楚,因此需要通过NAMD模拟进一步研究。在实验中,560纳米的激发光能够作用于SubPcH的强吸收带,而ZnTPP在此波长下的吸收较弱,因此以SubPcH的激发为主。在随后的弛豫过程中,系统主要发生从SubPcH到ZnTPP的正向能量转移。因此,为了与实验结果进行比较,我们采用了LR-TDDFT框架下的FSSH-CPA方法,对ZnTPP-SubPcH和ZnTPP-SubPcF在甲苯和苯腈中的情况进行了NAMD模拟,由于它们在可见光区域的吸收较强,初始激发主要集中在SubPcH或SubPcF单元上。首先,在FSSH-CPA框架中,模拟所包含的激发态是根据它们的垂直激发能量和振子强度来选择的。以位于SubPcH上、与ZnTPP-SubPcH在甲苯中523纳米(2.37电子伏特)吸收带相关的LE态为例,用于模拟的初始态必须满足两个条件:(1)激发能量要在2.37±0.16电子伏特范围内;(2)该态的振子强度要高于从0到1之间的随机值。相应的分布如图5(a)所示。采用相同的程序,确定了532纳米(ZnTPP-SubPcH在苯腈中,图5(b))、545纳米(ZnTPP-SubPcF在甲苯中,图5(c))以及552纳米(ZnTPP-SubPcF在苯腈中,图5(d))激发时的初始态。对于在甲苯中的ZnTPP-SubPcH(图5(a)),523纳米附近的激发主要使S3-S5态被占据,其中S3的贡献最大,其次是S4和S5,而S1、S2和S6的贡献较弱。在苯腈中(图5(b)),占据态主要分布在S3-S10之间,其中S6-S9的贡献最强,其余态的贡献较弱。对于在甲苯中的ZnTPP-SubPcF(图5(c)),545纳米的激发主要使S4-S8态被占据,其中S8的贡献最大。在苯腈中(图5(d)),552纳米的激发以S8为主,其次是S5,S6和S7的贡献相对较小且相近。图5显示了在FSSH-CPA模拟中用于ZnTPP-SubPcH在甲苯(a)和苯腈(b)中,以及ZnTPP-SubPcF在甲苯(c)和苯腈(d)中的选定激发态的初始分布情况。需要注意的是,图5中所显示的态分布代表选定的初始条件,而图6中所展示的种群数量则对应于动力学开始时的平均绝热态种群数量。从选定的1000个初始条件中,我们为每个条件模拟了300条FSSH-CPA轨迹。对于ZnTPP-SubPcH,在甲苯和苯腈中的总模拟时间均为2000飞秒。对于ZnTPP-SubPcF,在甲苯中的总模拟时间同样为2000飞秒,而在苯腈中则延长至10皮秒。所有模拟均采用1飞秒的时间步长。报告的结果是基于1000×300条轨迹的集合平均值得出的。图6展示了两种体系相关绝热电子态的随时间变化的种群数量。我们首先讨论ZnTPP-SubPcH的激发态种群分布情况。在甲苯中(图6(a)),在523纳米处激发后,初始种群主要分布在S3-S10之间,只有少量来自S1-S2。在1000飞秒内,S3-S10的种群数量从0.91降至0.14,而同一时期S1-S2的种群数量从0.09上升至0.86。在苯腈中(图6(b)),在1500飞秒内,S3-S10的初始种群数量从1.00降至0.02,而S1-S2的种群数量从0.00上升至0.98。对于在甲苯中的ZnTPP-SubPcF(图6(c)),初始种群主要集中于S4-S10态。在1000飞秒内,S4-S10的种群数量从1.00降至0.00,而S1-S3的种群数量从0.00上升至1.00。在苯腈中(图6(d)),在10皮秒内,S5-S10态的种群数量从1.00降至0.03,而S1-S4的种群数量从0.00上升至0.97。上述结果表明,经过光激发后,两种体系的种群演化都经历了一个从较高能级激发态向较低能级激发态的顺序弛豫过程。为了定量描述这种逐步弛豫行为,根据上述描述,所涉及的10个激发态被分为两组:最初被归入M流形的那些种群数量较多的态,以及被归入Q流形的最低能级激发态。随后,使用逐步的M→Q动力学模型来分析时间依赖的种群动态,其中种群数量表示为:Mt=M0e?k1t,Qt=M0+Q0?Mtin,其中M0和Q0分别对应于t=0时M流形和Q流形的初始种群数量。k1则对应于M→Q过程的速率常数。在这个顺序动力学模型中,我们得到了两种体系在甲苯和苯腈中的相应时间常数。对于ZnTPP-SubPcH,其在甲苯中的τ1值为162.1飞秒,在苯腈中为350.9飞秒。相比之下,对于ZnTPP-SubPcF,其在甲苯中的τ1值为162.9飞秒,而在苯腈中则显著增加至2.9皮秒。实验测得的ZnTPP-SubPcH在甲苯中的弛豫时间为6.81皮秒,而本NAMD模拟得到的相应弛豫时间为162.1飞秒。尽管这两个时间尺度之间存在定量差异,但这种差异可能源于理论模型的固有近似以及模拟条件与实验条件之间的差异。不过,无论是实验研究还是理论研究都一致认为,能量从SubPcH向ZnTPP的转移是该体系中主要的 photoinduced过程。因此,尽管绝对时间常数存在偏差,但实验观察结果与理论模拟在本质机制上是一致的。图6展示了ZnTPP-SubPcH在甲苯(a)和苯腈(b)中,以及ZnTPP-SubPcF在甲苯(c)和苯腈(d)中,相关绝热激发态的随时间变化的种群数量。总体而言,两种体系在光激发后都表现出从较高能级激发态向较低能级激发态的级联式弛豫行为。对于ZnTPP-SubPcH,苯腈溶剂在一定程度上减缓了这一弛豫过程。而对于ZnTPP-SubPcF,苯腈溶剂则引发了明显的动力学延迟效应。这些结果表明,溶剂环境显著影响了两种体系的激发态弛豫动力学,而且ZnTPP-SubPcF在苯腈中的弛豫受到更强烈的抑制。为了探究这两种体系行为背后的机制,我们还分析了ZnTPP-SubPcH和ZnTPP-SubPcF中相邻绝热态之间的能隙分布,以及不同绝热态之间的平均绝对非绝热耦合强度,如图7所示。众所周知,非绝热跃迁既受能隙的影响,也受非绝热耦合的影响,非绝热耦合向量与能隙之间的关系可表示为:dij=Ψi?H?RΨjEj?Ei。与状态i和j相关的非绝热耦合向量记为dij,其中Ei和Ej分别是它们的能量值。随后,定义Λij=v×dij为连接绝热态的TDNAC。因此,能隙越小,平均非绝热耦合强度往往越大,非绝热跃迁也就越快。我们首先讨论ZnTPP-SubPcH的情况。在甲苯中(图7(a)),相邻态之间的能隙主要分布在0–6千卡/摩尔范围内,平均值为1.75千卡/摩尔,而所有相邻电子态的平均TDNAC值为15.18皮秒?1。在苯腈中(图7(b)),相邻态之间的能隙主要分布在0–4千卡/摩尔范围内,平均值为1.15千卡/摩尔,所有相邻电子态的平均TDNAC值为15.66皮秒?1。这些发现表明,苯腈中较小的能隙和稍大的TDNAC有利于非绝热跃迁,不过实际的弛豫时间仍然由激发态的排序以及不同弛豫途径之间的竞争共同决定。对于在甲苯中的ZnTPP-SubPcF(图7(c)),相邻态之间的能隙同样主要分布在0–6千卡/摩尔范围内,平均值为1.59千卡/摩尔。同时,所有相邻电子态的平均TDNAC值为17.41皮秒?1。在苯腈中(图7(d)),相邻态之间的能隙主要分布在0–10千卡/摩尔范围内,平均值为2.38千卡/摩尔,所有相邻电子态的平均TDNAC值为15.79皮秒?1。与ZnTPP-SubPcH相比,ZnTPP-SubPcF在甲苯中的弛豫动力学更快。对于ZnTPP-SubPcF,高能区域存在的某种分布可能是其弛豫过程相对较慢的主要原因。这些结果表明,弛豫动力学并非仅由能隙或TDNAC的数值大小决定,还受到激发态的顺序以及是否存在高效弛豫途径的影响。通过NAMD模拟,我们还获得了两种体系在两种环境中的不同片段上的电子和空穴权重随时间的变化情况,不同激子组分的时间演化情况,以及电子-空穴分离指数(D指数)。我们首先讨论ZnTPP-SubPcH的情况。如图8(a)和(b)所示,在甲苯中1000飞秒内,SubPcH上的空穴权重从0.88降至0.13,而ZnTPP上的相应权重则从0.12升至0.87,表明空穴从SubPcH转移到了ZnTPP。在同一时期,SubPcH上的电子权重从0.93降至0.17,而ZnTPP上的相应权重则从0.07升至0.83,说明电子也同样从SubPcH转移到了ZnTPP。由于空穴和电子共同向同一个片段移动,这种行为被归因于从SubPcH到ZnTPP的光诱导能量转移,这与实验解释在定性上是一致的。如图10(a)中各种激子类型的随时间变化权重所示,对于在甲苯中的ZnTPP-SubPcH,最初占主导地位的LE激子SubPcH?,再加上CT激子SubPcH?ZnTPP+的次要作用,会在1000飞秒内逐渐转变为LE激子ZnTPP?。与此同时,D指数(图S3(a))在初始阶段迅速下降,最终趋于0,表明电子和空穴并没有持续的空间分离,而是通过同方向的协同转移在两个片段之间重新分配激发能量。因此,在甲苯中的ZnTPP-SubPcH中,主要的photoinduced过程是从SubPcH到ZnTPP的能量转移,而非真正的长寿命电荷分离。这一计算结果在定性上支持了实验提出的从SubPcH到ZnTPP的单重态能量转移机制。然而,在苯腈中,ZnTPP-SubPcH的激发态动力学则明显不同。如图8(c)和(d)所示,500飞秒内空穴就从SubPcH转移到了ZnTPP。而电子的轨迹则有所不同:其在SubPcH上的权重在200飞秒内从0.91降至0.53,然后在1500飞秒内恢复到0.98;而在ZnTPP上的权重则在最初的200飞秒内从0.09升至0.47,随后在1500飞秒内降至0.02。这些结果表明,从SubPcH到ZnTPP的能量转移发生在最初的200飞秒内,之后电子转移的方向与空穴转移的方向相反。空穴权重在超快时间尺度上迅速从SubPcH转移到ZnTPP,之后一直保留在ZnTPP片段上。相比之下,电子权重则先发生局部重新分布,然后逐渐重新转移到SubPcH上,最终主要分布在SubPcH片段上。也就是说,在苯腈中,电子和空穴最终分别定域在SubPcH和ZnTPP两个片段上,从而形成了明显的空间分离。图10(b)的分析表明,当BN中的SubPcH片段被激发后,最初形成的物种以LE激子SubPcH?为主,同时还有少量CT激子SubPcH+ZnTPP?的贡献。在200飞秒内,这些物种逐渐转变为LE激子ZnTPP?和CT激子SubPcH?ZnTPP+。在接下来的1500飞秒里,ZnTPP?激子以及其他最初存在的物种持续衰减,最终所有的激子种群都转化为CT激子SubPcH?ZnTPP+。与此同时,D指数(见图S3(b))从较小的数值开始持续上升,最终稳定在相对较高的水平,进一步证实了BN中的系统从LE态向CT态演变。因此,对于ZnTPP-SubPcH而言,溶剂极性显著影响了激发后的弛豫路径:在极性较弱的甲苯中,系统的动态过程以能量转移为主;而在极性较强的BN中,系统则更容易转变为稳定的CT态。我们还研究了ZnTPP-SubPcF的激发态动态,如图9所示,其动态行为与ZnTPP-SubPcH有本质不同。在甲苯中(见图9(a)和(b)),SubPcF上的空穴权重从0.85降至0.00,而ZnTPP上的相应值则从0.15升至1.00,表明在500飞秒内空穴就迅速从SubPcF转移到了ZnTPP上。同时,电子也在最初的100飞秒内从SubPcF转移到ZnTPP,但在之后的500飞秒内又发生了反向转移,从ZnTPP回到SubPcF。这些结果表明,早期阶段的动态特征是激发能从SubPcF向ZnTPP重新分配,之后系统进入CT通道,最终形成以CT为主的最终状态。对不同类型激子的时变权重分析(见图10(c))显示,在甲苯中,激发最初主要形成LE激子SubPcF?。在最初的100飞秒内,LE激子ZnTPP?的权重有所增加,但大约在500飞秒时,ZnTPP?和SubPcF?都转变为CT激子SubPcF?ZnTPP+,后者最终成为主导成分。更重要的是,D指数(见图S3(c))在几百飞秒内迅速上升并稳定在较高数值,表明在甲苯中系统已能够实现明显的电子-空穴空间分离。这些结果表明,氟取代使得系统即使在极性较弱的溶剂中也能超越单纯的激发能转移,有效进入电荷转移通道,最终形成以CT为主的最终状态。下载:下载高分辨率图像(379KB)下载:下载全尺寸图像图9. 甲苯中ZnTPP-SubPcF的空穴权重时变曲线(a)和电子权重时变曲线(b),以及苯腈中ZnTPP-SubPcF的空穴权重时变曲线(c)和电子权重时变曲线(d)。下载:下载高分辨率图像(402KB)下载:下载全尺寸图像图10. 甲苯中ZnTPP-SubPcH的激子权重时变曲线(a)和苯腈中的曲线(b),以及甲苯中ZnTPP-SubPcF的曲线(c)和苯腈中的曲线(d)。图9(c)和(d)显示,在BN中,系统同样表现出明显的电荷转移行为,尽管其动态过程比在甲苯中更为缓慢。SubPcF上的空穴权重首先从0.86迅速降至0.40,然后再降至0.00,而ZnTPP上的权重则从0.14急剧上升到0.60,最终达到1.00,表明系统在皮秒时间尺度上进入了超快电荷转移状态。随后,空穴继续从SubPcF片段转移到ZnTPP片段,而电子则从ZnTPP反向转移到SubPcF片段。对图10(d)的分析显示,激发后,主要的贡献来自SubPcF?激子,其次是ZnTPP?。此外,ZnTPP?的权重在皮秒时间尺度上从0.14先上升到0.57,然后与SubPcF?一起转变为CT激子SubPcF?ZnTPP+。D指数(见图S3(d))呈单调上升趋势,并在较长时间内逐渐接近饱和值,表明BN有助于形成更完整、更稳定的CT态。值得注意的是,尽管ZnTPP-SubPcF在甲苯和BN中最终都会趋向于以CT为主的最终状态,但在BN中的演化更为平滑,且持续时间更长,这表明强极性环境不仅有利于CT态的稳定,还改变了激发态弛豫的时间尺度和路径。综合来看,两种系统在两种环境中的动态行为表明,取代基效应和溶剂环境共同决定了激发后的演化路径。对于ZnTPP-SubPcH,在甲苯中以SubPcH向ZnTPP的能量转移为主,而在BN中则进一步向电荷转移方向发展。相比之下,对于ZnTPP-SubPcF,氟取代显著增强了受体片段的吸电子性,使得系统即使在甲苯中也能快速进入电荷转移通道,并在BN中形成更稳定的CT态。换句话说,H取代系统的动态主要受溶剂极性控制,而F取代系统则表现出更强的内在电荷转移倾向。为进一步阐明非绝热动力学过程中的激发态电荷转移和能量转移机制,我们基于电子和空穴种群引入了两个基于片段层面的描述符。第一个描述符定义为qCTt=hAt?eAt,其中hAt和eAt分别表示片段A上的空穴和电子种群数量。这一数值反映了激发态演化过程中电子-空穴分离的程度。当qCTt接近0时,电子和空穴主要位于同一片段上,对应于没有明显电荷分离的LE态。相反,qCTt的增加表明电子和空穴逐渐定位在不同的片段上,体现了CT态的形成。qCTt的符号还可以进一步指示电荷转移的方向。如图S4所示,甲苯中的ZnTPP-SubPcH在整个模拟过程中qCTt都接近0,表明在激发态演化过程中没有发生显著的电荷转移。而BN中的ZnTPP-SubPcH以及两种溶剂中的ZnTPP-SubPcF的qCTt则持续上升,最终接近1,表明形成了空间分离的CT态。此外,为了进一步识别整个激子的迁移情况,我们还定义了LAt=hAt+eAt2,它表示片段A上电子-空穴对的平均位置。与qCTt不同,这一数值反映的是激子的整体迁移,而非电子和空穴之间的空间分离。因此,这两个描述符为激发态动力学提供了互补的信息。如图S5所示,甲苯中的ZnTPP-SubPcH的LAt存在显著变化,激子种群逐渐从SubPcH片段转移到ZnTPP片段,而qCTt在整个模拟过程中始终接近零。这表明电子和空穴是协同迁移的,没有明显的空间分离,这是能量转移的特征。相比之下,BN中的ZnTPP-SubPcH以及两种溶剂中的ZnTPP-SubPcF,LAt逐渐向0.5靠拢的同时,qCTt也持续上升。这种行为是由于电子和空穴逐渐分离并形成CT态所导致的,而非单纯的激子迁移。因此,这些系统中的激发态动力学主要是由电荷转移而非能量转移所主导。综合qCTt和LAt的分析,可以明确区分能量转移和电荷转移过程。具体而言,甲苯中的ZnTPP-SubPcH表现出典型的能量转移过程,而其余系统在光激发后则以电荷转移为主。基于上述计算结果,图11总结了甲苯和BN中ZnTPP-SubPcH和ZnTPP-SubPcF的光诱导动态演化机制。对于ZnTPP-SubPcH,在甲苯中约523纳米(2.37电子伏特)的光激发下,初始激发主要集中在SubPcH片段上,如图4(a)中的LE态所示,对应于图2中的S3和S4态,而更高能级态(S5-S10)的贡献相对较小。随后,NAMD模拟显示,激发态迅速弛豫,同时电子密度在各个片段之间重新分布,逐渐形成集中在ZnTPP部分的低能级、类LE的电子分布。对不同片段上电子和空穴种群的时变分析进一步表明,这一过程本质上是从SubPcH到ZnTPP的光诱导能量转移(PEnT)。在BN中约532纳米(2.33电子伏特)的光激发下,初始激发同样主要对应于集中在SubPcH片段的LE态,如图4(b)所示,主要涉及图2中的S7和S8态,S9和S10态的贡献较小。在早期动态过程中,系统表现出能量重新分配和电荷转移加剧的耦合演化,导致电子密度逐渐在ZnTPP片段上积累,形成以ZnTPP为主的类LE的瞬态电子分布。这种中间构型代表了沿连续弛豫路径的瞬态电子结构,而非孤立的稳定状态。随着动态过程的推进,系统逐渐演变为CT态SubPcH?ZnTPP+,对应于图2中的S1和S2态。因此,这一过程可视为在PEnT辅助下的光诱导电子转移(PET)途径,随着CT态的建立,PET机制在后期逐渐占据主导地位。对于ZnTPP-SubPcF,在甲苯中约545纳米(2.27电子伏特)和BN中约552纳米(2.25电子伏特)的光激发下,其动态行为与BN中的ZnTPP-SubPcH极为相似。初始激发同样主要集中于SubPcF片段的LE态,如图4(c)和(d)所示,对应于图3中的S7和S8态,更高能级态的贡献相对较小。在早期动态过程中,系统表现出能量重新分配和电荷迁移趋势的同步激活,逐渐形成以ZnTPP为中心的类LE的瞬态电子分布。随后,系统进一步演变为CT态SubPcF?ZnTPP+,最终达到以CT为主的弛豫终点,其中BN环境显著提升了CT态的稳定性和种群数量。总体而言,该系统的激发态动力学可统一描述为一种协同演化过程,即由PEnT驱动的早期弛豫之后是主导电荷转移阶段的PET过程,最终共同形成CT态。下载:下载高分辨率图像(245KB)下载:下载全尺寸图像图11. 甲苯中ZnTPP-SubPcH(上)和苯腈中ZnTPP-SubPcF(下)的光诱导动力学过程。为进一步评估CPA的潜在影响,我们使用了经过优化的ωB97XD*函数,对甲苯和BN中的ZnTPP-SubPcH和ZnTPP-SubPcF进行了额外的S?几何结构优化。选择ωB97XD*函数是因为LC-PBE0*缺乏 properly参数化的D3色散校正,而我们之前的研究已经表明,经过优化的ωB97XD*函数能够得到与优化后的LC-PBE0*函数类似的结果[44]。优化后的S?极小值与相应S0结构之间的均方根偏差(RMSDs)在甲苯中分别为ZnTPP-SubPcH为1.19埃,ZnTPP-SubPcF为1.80埃;在BN中则为ZnTPP-SubPcH为1.61埃,ZnTPP-SubPcF为1.81埃。这些结果表明,激发态势能面上存在结构弛豫现象。虽然这说明CPA所依赖的假设并未完全满足,但结构变化的幅度仍然较为适中。不过,考虑到对大型分子系统进行完全耦合的非绝热动力学模拟需要巨大的计算成本,FSSH-CPA在当前的计算条件下实现了计算效率与模拟精度之间的合理平衡。因此,本研究中呈现的模拟结果主要是为了揭示系统中激发态弛豫、能量转移和电荷转移的主要动态特征及背后的物理机制,而动态时间常数的定量数值则应在CPA适用范围内进行解读。通过分析S?和S?优化几何结构下的激发态特性(见图S7),进一步评估了结构弛豫对激发态电子结构的影响。对于ZnTPP-SubPcH,弛豫后相关激发态的电子特性和能量排序基本保持不变。相比之下,ZnTPP-SubPcF的较高能级激发态谱系出现了变化,包括几种LE态和CT态之间的部分重排。不过,负责后续电荷转移过程的最低能级CT态仍保持了其电子特性,从最初以SubPc为中心的LE态到低能级CT谱系的整个光诱导弛豫序列依然未变。这些观察结果表明,CPA近似主要会影响高能级激发态之间种群重分布及弛豫时间尺度的定量描述,而不太可能改变本研究中主要的电荷转移路径或定性机制结论。需要指出的是,这类给体-受体体系中的激发态能量/电子转移过程不能简单地归类为传统的F?rster共振能量转移(FRET)或Dexter能量转移机制,因为后者是基于费米黄金法则建立的,最适合初始态和终态较为明确的弱耦合给体-受体体系。而在当前研究的体系中,多个低能级激发态几乎处于简并状态,且局部激发与电荷转移特性之间存在强烈混合。在这种条件下,微扰弱耦合模型可能不再适用。因此,研究采用了NAMD模拟来精确追踪激发态种群及电子-空穴分布的超快演化过程。结果表明,激发态的弛豫涉及强耦合的LE/CT态混合以及随时间变化的电荷重分布,而非纯粹的FRET或Dexter型能量转移过程。将这一过程简单归为FRET或Dexter能量转移是一种过度简化,对于当前这种强耦合的激发态体系而言甚至可能具有误导性。

4. 结论

在本研究中,通过结合LR-TDDFT分析与NAMD模拟,系统研究了ZnTPP-SubPcH和ZnTPP-SubPcF在Tol和BN溶剂中的激发态电子结构、吸收特性以及光诱导动力学行为。计算结果显示,两种体系中的主要亮态均源自LE态,而CT态的振子强度较低,因此多为暗态。与含氢取代的体系相比,氟取代显著降低了低能级CT态的能量,增强了电荷分离的趋势,使得ZnTPP-SubPcF在不同环境中更易进入电荷转移通道。动力学模拟进一步表明,在Tol溶剂中,ZnTPP-SubPcH的主要过程为从SubPcH到ZnTPP的PEnT传递,这与实验提出的SubPc→ZnTPP能量转移路径定性一致。在BN溶剂中,该体系会进一步演化为一种独特的CT态。而对于ZnTPP-SubPcF而言,其在BN溶剂中更容易形成稳定的CT态,同时也能在Tol溶剂中进入电荷转移通道,这说明氟取代显著提升了其电荷转移能力。总体而言,取代基的修饰与溶剂环境对ZnTPP-SubPc体系的激发态性质及超快动力学行为具有显著的协同调控作用。这项工作不仅为相关实验观察提供了微观理论解释,还为卟啉-SubPc D-A体系及先进有机光电材料的发展奠定了重要的理论基础。

作者贡献说明:
梅杰·江:数据整理,初稿撰写。
熊思宇:形式分析。
龙晓雪:形式分析。
刘丽霞:软件应用。
张舒:软件应用。
杨佳佳:监督指导、项目管理、资金获取,审稿与编辑工作。
相关新闻
生物通微信公众号
微信
新浪微博

热点排行

    今日动态 | 人才市场 | 新技术专栏 | 中国科学人 | 云展台 | BioHot | 云讲堂直播 | 会展中心 | 特价专栏 | 技术快讯 | 免费试用

    版权所有 生物通

    Copyright© eBiotrade.com, All Rights Reserved

    联系信箱:

    粤ICP备09063491号