电流辅助烧结过程的多物理场耦合体系模拟

电流辅助烧结（electriccurrent-assistedsintering，ECAS）是在粉体颗粒间通入电流辅助烧结（electriccurrent-assistedsintering，ECAS）是在粉体颗粒间通入交流或直流脉冲大电流（上千甚至上万安培电流）实现快速烧结和致密化。相比于粉末激光烧结技术，电流辅助烧结产品（例如陶瓷部件）具有更好的力学性能，在降低成本、减小能耗的同时能够将烧结时间缩短到十几分钟，展示了诱人的工业应用前景。到目前为止，电流辅助烧结在材料制备方面已经取得了诸多进展，许多特殊的高性能材料已经能够利用电流辅助烧结技术制造。为了进一步提高电流辅助烧结技术的功用和效率，需要对烧交流或直流脉冲大电流（上千甚至上万安培电流）实现快速烧结和致密化。相比于粉末激光烧结技术，电流辅助烧结产品（例如陶瓷部件）具有更好的力学性能，在降低成本、减小能耗的同时能够将烧结时间缩短到十几分钟，展示了诱人的工业应用前景。到目前为止，电流辅助烧结在材料制备方面已经取得了诸多进展，许多特殊的高性能材料已经能够利用电流辅助烧结技术制造。为了进一步提高电流辅助烧结技术的功用和效率，需要对烧结过程中的作用机制进行更深入地理解和定量化描述。[1][2][3-5]，涉及应力（应变）场、温度场、电场和浓度场之间的相互作用，整个烧结过程被看作是一个多物理场相互作用的耦合过程。在此指导思想下，许多科研人员提出了多物理场耦合方法[6-9][10][11][12-14]以及电流对激活能的影响三个方面。在传统烧[15-17]，但已有不少研究表明，在电流活[18-22]Kondo[21]研究发现，对于碳化锰复合陶瓷，当电流密度小于600A·cm-2时，活化能保持不变，当电流密度大于600A·cm-2Garay等[22]关于镍钛复合材料的研究发现，当电流密度为1527A·cm-2292kJ·mol-12036A·cm-2时，活化能降低到86kJ·mol-1。因此在焦耳热效应和电迁移作用的基础上，还需要考虑电流对激活能的影响。本文结合相场法，推导了一组耦合了热、力、电和扩散场变量且适用于非等温扩散问题的相场模型控制方程。根据前人的实验研究基础[本文结合相场法，推导了一组耦合了热、力、电和扩散场变量且适用于非等温扩散问题的相场模型控制方程。根据前人的实验研究基础[21-22]，假设了碳化硅（SiC）材料电流与激活能之间的函数关系，通过将扩散系数引入耦合方程，在模拟中实现了活化作用可实时影响扩散行为。1.控制方程在非平衡体系（V）中，在单相多晶粒情况下不考虑化学反应和塑性变形，可将热-力-电-扩散耦合系统总自由能扩展为式（1）[10,23,24]。FB=∫VfB(T,ε(r),D,c,ξi,∇c,∇ξi)dVFB=∫VfB(T,ε(r),D,c,ξi,∇c,∇ξi)dV的梯度项。这里比自由能（f）是单位体积的自由能，可以由式（2）得出。BfB(T,ε(r),D,c,ξi,∇c,∇ξi)=f(0)B(T,ε(r),D,c)+f(1)B(c,ξi)+f(2)B(∇c,∇ξi)fB(T,ε(r),D,c,ξi,∇c,∇ξi)=fB(0)(T,ε(r),D作[10]，f(0),ε(r),f(1),ξ)和f(2),▽ξ)可以由式式表示。BBiBif(0)B(T,ε(r),D,c)=f00+a(1)ijεij+a(2)iDi+a(3)c+a(4)T+b(1)jiklεjiεkl+b(2)jiDjDi+b(3)c2+b(4)T2+b(12)jikεjiDε(r),D,c)=f00+aij(1)εij+ai(2)Di+a(3)c+a(4)T+bjikl(1)εjiεkl+bji(2)DjDi+b(3)c2+b(4)T2+bjik(12)εjiD式中：式中：F为系统总自由能；f为比自由能；T、ε(r)和D分别为绝对温度、应变张量的BB弹性部分和电位移（电通量密度弹性部分和电位移（电通量密度ξ分别为物质单位质量浓度和表征晶粒取向的序参i数，ξ0~1连续变化，其中i=和▽ξ为c和ξi i i式中：f式中：f(0)(T,ε(r),D,c)是体系在（拟）平衡状态下的比自由能；f(1)(c,ξ)表示包含BBi(2),▽ξ)是c和ξ的梯度场。参考之前的工Biif(1)B(c,ξf(1)B(c,ξi)=−A12(c−12)2+A24(c−12)4+A34(c−1)4+A44c4+⎡⎣−A52c2(∑i=1Nξi)2+A64(∑i=1Niξi)4⎤ξi)=−A12(c−12)2+A24(c−12)4+A34(c−1)4+A44c4+[−A52c2(∑i=1Nξi)2+A64(∑i=1Niξi)4]+∑f(2)B=κc2|∇c|2+∑i=1Nκξi2|∇ξi|2fB(2)=κc2|∇c|2+∑i=1Nκξi2|∇ξi|2数；κ和κ为梯度能量系数。按照广义的Gibbs-Duhem公式可得式。cξi−Bσ=B)=Bμc=Bμξi=Biδδξi=δfBδξi应的广义力。将式（3）~式（5）带入式（6）可得式（7）~式（11）。σij=a(1)ij+2b(1)jiklεkl+b(12)jikDk+b(13)jic+b(14)jiTσij=aij(1)+2bjikl(1)εkl+bjik(12)Dk+bji(Ei=a(2)i+2b(2)ijDj+b(12)kεkl+b(23)ic+b(24)jTEi=ai(2)+2bij(2)Dj+bk(12)εkl+bi(23)c+−s=a(4)+2b(4)T+b(14)jiεji+b(24)jDi+b(34)c−s=a(4)+2b(4)T+bji(14)εji+bj(24)Di+bμc=a(3)+2b(3)c+b(13)jiεij+b(23)iDi+b(34)T−A1(c−12)+A2(c−12)3+A3(c−1)3+A4c3+(−c=a(3)+2b(3)c+bji(13)εij+bi(23)Di+b(34)T−A1(c−12)+A2(c−12)3+A3(c−1)3+A4c3+(−式中：式中：f(0)为自由能关于变量、和c的泰勒展开式并保留到二阶；参数Ba(1)a(1),a(2),a(3),a(4)，b(12)，b(13)，b(14)，b(23)，b(24)，b(34)为待定常数，具体物理意义与取ijijikjijiij值将在后文讨论为室温下稳态时的自由能密度（i=1,2,3,···,7）为现象学参00 i、、μ和μξ对cξii在随后的数值计算中，参考Sherief的工作在随后的数值计算中，参考Sherief的工作[24]，式（7）~式（11）中的待定参数分μξiμξi=−A5c2ξi+A6(ξi)3+A72∑i=1,j≠iN(ξi)2ξj−κξi2ξiμξi=−A5c2ξi+A6(ξi)3+A72∑i=1,j≠iN(ξi在参考状态下σ、E、s和μ的取值对以后的讨论没有影响，所以f =B ij i c 00ij0，a(1)=0，a(2)=0，a(4)=0，a(3)=0。ij在弹性范围内，有

)2C=2δjδ+δδ+δδkδδkl+μδikδjl+μδilδjk)。

ijkl与D、μsT共轭量之间线性关系的影响可忽略，因i i c此有b(2)=0，b(3)=0，b(4)=0。ij ijσ与D、ET之间的线性关系可忽略，所以取b(12)=0，b(24)=0。一般来ij j i jik j说，b(13)=-βδ，b(23)=εE，b(14)=-。ij 2 ij i ei ij ij式中：a为热扩散效应度量，β=(3λ+2μ)α，β=(3λ+2μ)α，δ

为克罗内克符1 T 2 c ij号，λ、μ为拉梅常数，αT

为线性热膨胀系数，αc

为线性扩散膨胀系数。引入式（12）~式（15）几何方程，将几何方程式（12）带入式（7）中，结合动量方程式（13），可得位移表示的动量方程式（16）。εij=(ui,j+uj,i)/2εij=(ui,j+uj,i)/2动量方程:σij,j+Fi=ρu¨i动量方程:σij,j+Fi=ρu¨i商方程:Ts˙=−qi,i+ρμcc˙商方程:Ts˙=−qi,i+ρμcc˙:qi=−kt:qi=−ktT,i热传导定律:qi=−ktT,iμui,jj+(λ+μ)uj,ij−β1T,i−β2ρc,i+Fi=ρu¨iμui,jj+(λ+μ)uj,ij−β1T,i−β2ρc,i+Fi=ρu¨i式中：体力F可由三部分组成，如式（17）所示。iFi=ρeEi+PjEj,i+F′iFi=ρeEi+PjEj,i+Fi′可得热传导微分方程式（18）。ρCεT˙+β1Tε˙ji+aTρc˙1=(ktT,i),i+r˙+ρμcc˙ρCεT˙+β1Tε˙ji+aTρc˙1=(ktT,i),i+r˙+ρμ场的影响，则电场控制方程，如式（19）所示。∇⋅(γ∇V)=0∇⋅(γ∇V)=0式中：V为电势。相场法是通过引入一系列保守场变量和非保守场变量来表征微结构演化。保守场变量浓度（c）的演化可用Cahn-Hilliard（C-H）非线性扩散方程来描述[25-26]，并考虑电迁移的作用，可得浓度场控制方程，如式（20）所示。c˙=∇⋅M(∇μc−Z∗|e|∇V/Ω)c˙=∇⋅M(∇μc−Z∗|e|∇V/Ω)式中：式中：ρ表示体电荷密度，ρE表示与电场相关的体力，P为极化矢量，F'为需要考eeiji虑的其他形式的体力项。将式（9）带入熵方程式（14）中，结合热传导定律式（15），虑的其他形式的体力项。将式（9）带入熵方程式（14）中，结合热传导定律式（15），式中：式中：r˙=j2i/γr˙=ji2/γ为焦耳热源，γ为电导率，j为电流密度分量。此外，忽略磁i因此，总的扩散驱动力可以定义如式（21）所示。F=∇因此，总的扩散驱动力可以定义如式（21）所示。F=∇μc−Z∗|e|∇V/ΩF=∇μc−Z|e|∇V/Ω如式（22）所示。ξ˙i=−Lξiμξiξ˙i=−Lξiμξi（23）所示。[1−c3(10−15c+6c2)]+Msurfc(1−c)+MGB∑ξiξjM−15c+6c2)+Mvap[1−c3(10−15c+6c2)]+Msurfc(1−c)+MGB∑ξiξjMk=DkVm/此外，非保守场变量（此外，非保守场变量（ξ）的演化可用Ginzburg-Landau（G-L）方程来描述[25,27]，i[27]，如式i式中式中为有效电荷数为电子电荷为电势Mvolvap、Msurf和MGB分别为体积扩散、蒸发凝结、表面扩散和晶界扩散迁移率，可简化取为10Mvap=Mvol=0.1MGB=M [27]；D为各扩散形式对应扩散系数为摩尔体积为气体常数surfkmk0为各扩散形式对应扩散激活能。bk综上所述，式（16）、式（18）、式（19）、式（20）和式（22）构成了一整套热-力-电-扩散强耦合控制方程。式（16）为考虑温度场、电场、浓度场影响的动量方程，热传导方程式（18）包含应变场、电场和浓度场的影响，扩散方程式（20）考虑了温度场、电场和应力场的影响。该强耦合控制方程忽略了温度场、应变场和浓度场对电场分布的影响。模型与参数选取宏观烧结模拟分布，为后续微观模拟提供合理的边界条件。模拟采用COMSOL多物理场耦合的商业软Manièred[28]电流辅助烧结技术的烧结装置由石墨模具、硬质合金压头和SiC试样整体装配组成，图1所示，材料属性选取自COMSOL材料库。为了方便后期与实验对比，参考Bothara的实验工作[29]，设置模具内径为20。在考虑焦耳热和表面对流换热的基[30]。为了与实验条件相符，取模具上的参考点A作为温度反馈点，调节输入电流，并以内部中心点B作为后续微观结构演化模型的边界条件。图 1烧结模型示意图Figure 1.Schematicdiagramofthesinteringmodel下载:全尺寸图片幻灯片从20℃开始，分别在、601001700℃，随后保温10000。模拟结果表明，在不同升温速率下，稳定状态下A点和B点之间的温度相差约13℃，如图2（a）所示。由于SiC的电阻率随温度升高而逐渐降低，所以在温升率一定时，粉体中电流密度在升温期间逐渐增大，如图2（b）所示。图 2不同加热速率下温度和电流密度的变化B点温度变化点电密度变化Figure 2.Changesinthetemperatureandcurrentdensityatthedifferentheatingrates:(a)thetemperaturechangesatpointAandpointB;(b)thecurrentdensitychangesatpointB下载:全尺寸图片幻灯片模型与边条为了模拟烧结过程中的微结构演化，我们在计算域（宽度w=239nm，高度h=300nm）中建立平面四颗粒模型，球体半径R=60nm，初始接触颈部半径为r/2=15nm，对应初始相对密度约为0.60，见图3。温度和电流的边界条件取自宏观热-力-电模拟结果，如图2（a）和图2（b）所示。相应的边界条件为：四周温度（T）随时间（t）变化（T(t)）如图2（a）所示，左右两侧轴向压力为30MPa；电流密度（j）随时间t（t）变化（j(t)）如图2（b）所示，从左侧流入，右侧流出，上下两端固定且电绝缘；t初始温度为室温20℃，初始电流密度为0。图 3二维模型及边界条件示意图Figure 3.Schematicdiagramoftwo-dimensionalmodelandboundaryconditions下载:全尺寸图片幻灯片计算参数0SiC陶瓷粉末的材料参数如下所示：密度ρ=3.14×103kg·m-3，弹性模量E=0430GPa，泊松比ν=0.17，热容Cε

=0.68×103J·(kg·K)-1k

=170ttW·(m·K)-1，热膨胀系数αtt

=4.1×10-6K-1，扩散膨胀系数αc=2.68×10-6kg·m3，摩尔体积Vm

=1.246×10-5m3·mol-1Z=20D

=5m·s-e2，热扩散系数α=4.1448×109m2·(s2·K)-1，相对介电常数εe

=7.51，参考温度T0

=293K。现象学系数取如下：A1

=7.78×104J·m-3，A

=3.89×104J·m-3，A

=A=42353.4038×104J·m-3，A235

=1.167×104J·m-3，A

=2.334×103J·m-3，A

=7.7867×104J·m-3，晶界迁移系数L67

=2×10-6m3·(J·s)-1κc

=7.78×10-2J·m-1，κ =6.2224×10-19J·m-1。电导率的选取与激活能对扩散系数的影响电导率与温度变化的关系eSiC作为一种半导体材料，其电导率受温度的影响较大。当温度低于500℃时，SiC几乎绝缘；当温度高于500℃时，SiC电导率会随温度增加而增加[31]。为便于计算，采用以下拟合函数表示电导率（σ，S·m-1）与温度（T）的关系，如式（24）所示。eσe={0.9617573.23−0.34904T+4.09007×10−4T2T⩽500℃T>500℃σe={0.96175T 500℃73.234T2T>500℃2.4.2 激活能对扩散系数的影响2.4.2 激活能对扩散系数的影响许多研究表明，电流活化烧结中，电流会大幅降低激活能，极大地促进烧结的进行[17-20]。Kondo等[21]和Garay等[22]分别研究了电流对Mo2C和NiTi激活能的影响，结果如图4（a）所示。由图可知，当电流密度较低时，激活能基本保持不变；当电流密度高于阈值时，激活能将迅速降低。根据Ray等[17]的研究，取无电流情况下激活能Q=10200kJ·mol-1。在此基础上，假设最小活化能为420kJ·mol-1，约为Q的2/5，这与上述两种0材料的变化范围一致。据此假设拟合SiC活化能与电流密度关系为Q=-300tanh((J-4500)/1150)+720，如图4（b）所示。图 4不同材料激活能随电流密度的变化规律

C和NiTi；（b）SiC2Figure 4.Changesoftheactivationenergywiththecurrentdensity:(a)MoC2andNiTi;(b)SiC下载:全尺寸图片幻灯片结果与讨论电流辅助烧结致密化曲线以热-力-电-扩散强耦合理为基础，建立二维四颗粒模型作为代表单元，运用COMSOL的偏微分方程模块模拟烧结中微结构演化过程。定义颈长率（X），相对密度（ρ）和致密化参数（ψ）如式（25）~式（27）所示。r RX=r/2RX=r/2Rρr=1/[1+(1−X)3]ρr=1/[1+(1−X)3]ψR=(ρr−ρ0)/(1−ρ0)ψR=(ρr−ρ0)/(1−ρ0)式中：r式中：r为图3中两个相邻颗粒接触颈部的直径，由此可以得到致密化参数（ψ）随R时间变化情况。30MPa不变的情况下，取不同升温速率60100K·min-1），图5（a）所示，其中实验部分在900℃前无有效数据。由图可知，温升率越高，致密化所需时间越短，模拟与实验致密化曲线变化趋势基本[30]基本一致。由于颗粒颈部的初始曲率半径极大，在模拟烧结依然有一定致密化，随后致密化参数增长缓慢；在约1400℃处SiC致密化速率开始迅速增加，在约1600图5（b）所示；当致密化参数接近所不同的是，实验结果显示在60K·min-1和100K·min-1温升率下烧结后期（即致密化参数达到0.9之后）几乎没有差异，而模拟结果没有显示这种特征，对于这种差异需要作进一步研究。总的来说，模拟的致密化曲线与实验结果基本吻合，且均在1400℃时开始加速致密化，这在低温升率20K·min-1情况下最为明显。图 5不同温升率下模拟与实验致密化曲线与致密化速率随温度变化曲线Figure 5.Simulationandexperimentalresultsofthedensification(a)andthedensificationratecurveswithdifferenttemperatures(b)atdifferentheatingrates下载:全尺寸图片幻灯片以温升率60为例，观察电流辅助烧结中的微结构演化。烧结过程中电流密图6图6（a）所示，烧结颈部存在明显的电流集中效应，在快速致密化阶段（对应时间为1450~1680尤为明显。如图6（b）所示，在电流激活作用下，颈部扩散迁移率远大于其他区域，最多可高出约8个数量级，使得颈部加速生长，极大的促进烧结致密化。由浓度分布图6（c））及等值线形状演化图7）可看出，随着烧结的进行，接触颈部半径逐渐增大且颈部曲率减小，整体物质向内收缩，孔隙缩小逐渐接近圆形，最终消失，达到完全致密化。从烧结致密化曲线的演化可看出，在1450s1400℃后，由于此时颗粒在颈部位置的曲率较大，开始迅速致密化。随着烧结过程的进行，当孔隙逐步收缩并接近圆形s）时，颈部生长速度降低，致密化曲线逐渐平缓图7所示。图 6烧结过程中电流密度、扩散迁移率、物质浓度在不同时刻的分情况（60K·min-1）Figure 6.Distributionofthecurrentdensity(a),diffusionmobility(b),andconcentration(c)atdifferenttimeperiodsat60下载:全尺寸图片幻灯片图 7温升率为60K·min-1的电流辅助烧结致密化曲线及微结构演化Figure 7.DensificationcurveofECASandthemicrostructureevolutionat60K·min-1下载: