[发明专利]一种激光与物质非线性相互作用的并行化仿真方法

权利要求书

1.一种激光与物质非线性相互作用的并行化仿真方法，其特征在于，包括如下步骤：

S01：将含有激光与物质相互作用中非线性光学效应的等离子体速率方程引入麦克斯韦方程中，作为额外的微分方程，以获得激发的自由电子等离子体密度；

S02：将材料对光产生的克尔效应和多光子吸收非线性光学效应表达为克尔极化强度P_kerr(t)和多光子吸收极化强度P_MPI(t)；并将德鲁德模型反映的自由电子等离子体的产生对物质光学性质的影响表达为德鲁德极化强度P_Drude(t)，根据克尔极化强度P_kerr(t)、多光子吸收极化强度P_MPI(t)、德鲁德极化强度P_Drude(t)和电通量密度D(t)构建本构方程D′(t)，本构方程可表示为：

D′(t)＝D(t)+P_kerr(t)+P_MPI(t)+P_Drude(t)

D(t)＝ε₀ε_∞E(t)

其中：

P_Drude(t)由以下微分方程得到：

P_Kerr(t)由瞬时克尔效应引起，表达为：

P_MPI(t)由多光子吸收引起，表达式为：

式中：

D′(t)为电位移矢量；P_Drude(t)为德鲁德极化强度；P_kerr(t)为克尔极化强度；P_MPI(t)为多光子吸收极化强度；E(t)为电场强度；ε₀为真空介电常数；ε_∞为未激发材料的相对介电常数；ω_p为等离子体频率；为三阶极化率；Γ为电子碰撞频率；

S03：利用极化电流和极化强度的关系对麦克斯韦方程变形：

极化电流J：

变形后的麦克斯韦方程为：

式中：

E(t)为电场强度；μ₀为真空磁导系数；H(t)为磁场强度；J(t)为极化电流；

S04：由激光参数和物质参数根据变形后的麦克斯韦方程、本构方程和微分方程确立耦合方程组，对仿真区进行网格化，设在x方向上有i_max网格，在y方向上有j_max网格，z方向上有k_max网格；

S05：设置各向异性完全匹配层作为吸收边界，并采用时域有限差分法数值求解耦合方程组，对计算区和边界区进行统一数值仿真，由t时刻的电场E(t)求解t时刻的磁场H(t)，由t时刻的磁场H(t)求解t时刻的电位移矢量D′(t)；

S06：通过牛顿迭代法并行化数值求解本构方程，由t时刻的电位移矢量D′(t)更新t时刻的电场E(t)；

S07：采用龙格库塔方法求解自由电子等离子体密度速率方程，由电场E⁽ⁿ⁺¹⁾(t)更新激光激发的等离子体密度n_e、等离子体频率ω_p、极化强度P_drude(t)和P_MPI(t)；

以Δt为时间间隔，按照时间步进往前推进，重复步骤S05～S07，直到t＞最大时间t_max，结束仿真。

2.根据权利要求1所述的激光与物质非线性相互作用的并行化仿真方法，其特征在于，获得激发的自由电子等离子体密度的微分方程为：

式中：

n₀为未激发材料的折射率；E_g为材料的带隙；n_a为饱和电子密度；n_e为激光激发的自由电子等离子体密度；σ为逆韧致辐射截面；σ^(k)为多电子电离系数；I为激光光强；k是材料同时吸收的激光光子的个数；τ_R为电子-空穴复合时间。

3.根据权利要求1所述的激光与物质非线性相互作用的并行化仿真方法，其特征在于，由t时刻的电位移矢量D′(t)更新t时刻的电场E(t)，具体为：

S06.1：确定所有空间网格上依赖于三维空间坐标的光场矢量E的三个分量E_x、E_y、E_z矩阵：

E_x(i，j，k)→i＝1～i_max

E_y(i，j，k)→j＝1～j_max

E_z(i，j，k)→k＝1～k_max

式中：i为x方向上的网格序号；j为y方向上的网格序号；k为z方向上的网格序号；i_max为x方向上最大网格数；j_max为y方向上最大网格数；k_max为z方向上最大网格数；

构建一维列矩阵A，共3N行，N为空间网格的数目N＝i_max×j_max×k_max，一维列矩阵A中每三行对应于一个空间网格的光场三分量：

S06.2将每个空间网格上由本构方程构造的三个函数[X Y Z]的函数值写成一维行矩阵，把这个行矩阵重复3次，由此得到3行3列的矩阵，然后把所有N个空间网格上的所有3行3列的矩阵在列的方向上堆叠为一个3N行3列的矩阵B：

其中

式中：Δt为时间网格尺寸，n为时间网格数；D′⁽ⁿ⁺¹⁾(t)为t＝(n+1)Δt时刻的电位移矢量，E⁽ⁿ⁾(t)为t＝nΔt时刻的电场强度矢量，E⁽ⁿ⁺¹⁾(t)为t＝(n+1)Δt时刻的电场强度矢量，为t＝nΔt时刻的德鲁德极化强度矢量，为t＝(n-1)Δt时刻的德鲁德极化强度矢量，为t＝nΔt时刻的多光子吸收极化强度矢量，E_x为x方向的电场分量，E_y为y方向的电场分量，E_z为z方向的电场分量；

S06.3：确定所有空间网格各自对应的3*3的雅各比矩阵的逆矩阵：

其中，为雅各比矩阵的逆矩阵的9个分量；

构建一个3N行3列的矩阵C，每三行对应一个空间网格的雅各比矩阵的逆矩阵：

以.g为迭代次数，令A^g代表t＝(n+1)Δt时刻第g次迭代的电场值E^g，根据A^g+1＝A^g-C^g.*B^g，确定A^g+1，令A^g+1代表t＝(n+1)Δt时刻第g+1次迭代的电场值E^g+1，如果|E^g+1-E^g|＞eps，则采用E^g＝E^g+1迭代，并采用该电场值更新整合后的雅各比逆矩阵C和整合后的矩阵B，增加迭代次数，迭代过程直至|E^g+1-E^g|＜eps，获得新时刻下的电场E⁽ⁿ⁺¹⁾(t)＝E^g+1；其中eqs为误差范围。