一维PN结数值仿真原理(2)——有限差分法
本文最后更新于358 天前,其中的信息可能已经过时,如有错误请发送邮件到lysun26@163.com

有限差分法

差分格式

首先,我们需要对空间坐标进行离散化,也成为网格化(mesh)。网格化是器件仿真中十分重要的内容,合理的网格设置不仅可以帮助我们得到准确的结果,也能帮助我们减少模型的计算量,提高工作效率。在 TCAD 软件中,有许多的网格设置方法,但总归来说,都是通过一定的方式(比如说给定函数)对连续坐标进行离散化,把握好这个关键就很好理解 TCAD 中的一些网格设置。

将数据离散化之后,我们就可以将微分形式转换为差分形式。这个方法有很多,最简单的便是 Euler 法,还有 Runge-Kutta 法等,但这些方法的核心依旧是 Euler 法,感兴趣的可以去学一学微分方程数值解法。在这里,我们按下图所示的方式对数据离散化:

根据微分的定义,我们可以将 $x_k$ 处的值近似表示为:

$$
\frac{df}{dx_k}=\frac{f(x_{k+1})-f(x_k)}{d_k}
$$

这个便是其中的一种差分格式,叫做向前差分。在这里,我们考虑利用如下的差分形式。如下图所示,将中点记为 M。

令:

$$
\frac{df}{dx_k}=\frac{f(x_{M})-f(x_{M-1})}{h_k}
$$

由此可以得到二阶导为:

$$
\frac{d^{2}f}{dx_{k}^{2}}=\frac{f'(x_{M})-f'(x_{M-1})}{h_{k}}=\frac{\frac{f(x_{k+1})-f(x_k)}{d_k}-\frac{f(x_{k})-f(x_{k-1})}{d_{k-1}}}{h_{k}}
$$

由于 M 对应的点并不是我们网格上的点,因此,我们需要把 M 转化到 k 上。我们假设在每一段上线性变化,即:

$$
f(x_{M})=\frac{f(x_{k})+f(x_{k+1})}{2}
$$

代入到上面的式子中,可以得到:

$$
\frac{df}{dx_k}=\frac{f(x_{M})-f(x_{M-1})}{h_k}=\frac{\frac{f(x_{k})+f(x_{k+1})}{2}-\frac{f(x_{k})+f(x_{k-1})}{2}}{h_k}=\frac{f(x_{k+1})-f(x_{k-1})}{d_{k}+d_{k-1}}
$$

有了上面的关系,我们就可以把微分方程转化为差分方程了,如下:

$$
\begin{aligned}
&\frac{2}{d_{k-1}+d_{k}}(\frac{\varphi_{k+1}-\varphi_{k}}{d_{k}}-\frac{\varphi_{k}-\varphi_{k-1}}{d_{k-1}})=-(N_{k}+p_{k}-n_{k}) \\
& \frac{J_{p,k+1}-J_{p,k-1}}{d_{k-1}+d_{k}}=-R_{k} \\
& \frac{J_{n,k+1}-J_{n,k-1}}{d_{k-1}+d_{k}}=R_{k} \\
& \frac{n_{k+1}-n_{k-1}}{d_{k-1}+d_{k}}=\frac{J_{n,k}}{\mu_{n}}+n_{k}\cdot \frac{\varphi_{k+1}-\varphi_{k-1}}{d_{k-1}+d_{k}} \\
& \frac{p_{k+1}-p_{k-1}}{d_{k-1}+d_{k}}=-\frac{J_{p,k}}{\mu_{p}}-p_{k}\cdot \frac{\varphi_{k+1}-\varphi_{k-1}}{d_{k-1}+d_{k}} \\
& R_{k}=\frac{n_{k}p_{k}-1}{\tau_{p}(n_{k}+1)+\tau_{n}(p_{k}+1)}
\end{aligned}
$$

显然,在这个方程中,$J,R$ 都是中间量,可以消去,最后得到关于 $\varphi,n,p$ 的方程。

第一个方程可以写成:

$$
F_{\varphi}(n_{k},p_{k},\varphi_{k-1},\varphi_{k},\varphi_{k+1})=0
$$

将有关 p 的方程联立,可以得到:

$$
F_{p}(n_{k-1},n_{k},n_{k+1},p_{k},\varphi_{k-1},\varphi_{k},\varphi_{k+1})=0
$$

将有关 n 的方程联立,可以得到:

$$
F_{n}(p_{k-1},p_{k},p_{k+1},n_{k},\varphi_{k-1},\varphi_{k},\varphi_{k+1})=0
$$

由于这些方程都是对于内部的点,而非边界点,因此方程的个数共有 $3(N-2)$,未知量个数也是 $3(N-2)$ 个,因此理论上是可以解的。

方程个数这么多,得到理论解几乎是不可能的,因此在实际中常常采用迭代法。迭代的方法有很多,常用的是牛顿迭代法

牛顿迭代法

只有一个方程的情形

$$
f(x)=0
$$

牛顿迭代法的原理如图所示:

首先,我们需要给定一个初值 $x_0$,然后过这个点作切线。我们可以得到切线方程为:

$$
y=f'(x_{0})(x-x_{0})+f(x_{0})
$$

令 $y=0$,可以求得该切线与 x 轴的交点坐标。可以得到:

$$
x_{1}=x_{0}-\frac{f(x_{0})}{f'(x_{0})}
$$

从图中可以发现,该点距离最终的结果更近了一步。就这样迭代下去,便可以得到接近目标值的解。其迭代方程为:

$$
x_{k+1}=x_{k}-\frac{f(x_{k})}{f'(x_{k})}
$$

方程组情形

$$
\begin{aligned}
f_{1}(x_{1},x_{2},\cdots ,x_{n})&=0 \\
f_{2}(x_{1},x_{2},\cdots ,x_{n})&=0 \\
\vdots \\
f_{n}(x_{1},x_{2},\cdots ,x_{n})&=0
\end{aligned}
$$

首先,我们把自变量和因变量都写成向量的形式,令:

$$
\mathbf{X}=(x_{1}, x_{2},\cdots, x_{n})^T
$$

$$
\mathbf{F}=(f_{1}, f_{2}, \cdots, f_{n})^T
$$

则上述方程组可以写成:

$$
\mathbf{F(X)}=0
$$

当只有一个方程时,我们得到:

$$
x_{k+1}=x_{k}-\frac{f(x_{k})}{f'(x_{k})}
$$

我们可以对其进行推广,这是数学研究中十分常用的方法。首先,在矩阵中,除法表现为逆,故将 $f'(x_{k})$ 推广为一个矩阵的逆,记为 $\mathbf{F'(X^{k})}^{-1}$,然后我们就能够得到:

$$
\mathbf{X^{k+1}}=\mathbf{X^{k}}-\mathbf{F'(X^{k})}^{-1}\mathbf{F(X^{k})}
$$

而这个便是方程组时的牛顿迭代法的递推公式。接下来的问题就是这个 $\mathbf{F'(X^{k})}$ 是什么,要怎么求。

首先,我们知道 $\mathbf{X,F}$ 均为 $n\times 1$ 的矩阵,故 $\mathbf{F'(X^{k})}^{-1}$ 为 $n\times n$ 的矩阵,所以 $\mathbf{F'(X^{k})}$ 也为 $n\times n$ 的矩阵,这个矩阵便是雅可比矩阵(Jacobian matrix)。雅可比矩阵的形式如下:

$$
\mathbf{J} = \begin{bmatrix}
\frac{\partial f_1}{\partial x_1} & \cdots & \frac{\partial f_1}{\partial x_n} \\
\vdots & \ddots & \vdots \\
\frac{\partial f_n}{\partial x_1} & \cdots & \frac{\partial f_n}{\partial x_n}
\end{bmatrix}
$$

从这个具体形式中可以看出,雅可比矩阵描述了多元函数的局部线性近似的信息。如何理解雅可比矩阵是怎么导出的呢?可以这么简单理解:

考虑一个二元方程组,即:

$$
\begin{aligned}
f_1(x_1,x_{2})=0 \\
f_{2}(x_{1},x_{2})=0
\end{aligned}
$$

仿照着一元方程的做法,考虑 $f_1(x_1,x_2)$,根据全微分的知识,我们知道这个二元函数的变化率为:

$$
df_{1}^k=\frac{\partial f_1^k}{\partial x_1}dx_{1}^k+\frac{\partial f_1^k}{\partial x_2}dx_2^k
$$

$$
f_{1}^{k+1}=f_{1}^{k}+df_{1}^{k}=f_{1}^{k}+\frac{\partial f_1^k}{\partial x_1}dx_{1}^k+\frac{\partial f_1^k}{\partial x_2}dx_2^k
$$

这可以看作是 $x^k$ 点处的一小段直线的表达式,我们可以将其扩展成一个直线方程:

$$
f=f_{1}^{k}+\frac{\partial f_1^k}{\partial x_1}(x_1-x_{1}^k)+\frac{\partial f_1^k}{\partial x_2}(x_2-x_{2}^k)=f_{1}^{k}+\begin{bmatrix}
\frac{\partial f_1^k}{\partial x_1}&\frac{\partial f_1^k}{\partial x_2}
\end{bmatrix}
\begin{bmatrix}
x_1-x_{1}^k\\x_2-x_{2}^k
\end{bmatrix}
$$

同样的,我们考虑 $f_2$,就可以得到另一个方程,然后将两个方程联立,写成矩阵形式,便可以得到解了,就可以发现导数矩阵就是雅可比矩阵。

牛顿迭代法的流程

综上所述,我们可以得到求解方程组的流程,如下:

  1. 选择初始猜测向量 $\mathbf{X^{0}}=[x_1^0,\cdots,x_n^0]$.
  2. 计算函数值 $\mathbf{F(X^0)}$ 和雅可比矩阵 $\mathbf{J(X^0)}$.
  3. 计算雅可比矩阵的逆 $\mathbf{J(X^0)}^{-1}$.
  4. 计算下一个近似解 $\mathbf{X^1}=\mathbf{X^{0}}-\mathbf{J(X^{0})}^{-1}\mathbf{F(X^{0})}$.
  5. 重复步骤 2~4,直到满足收敛条件为止。

总结

通过牛顿迭代法,我们就可以设计程序,利用计算机对上面的 $3(N-2)$ 个方程求解。

有问题可以留言哦~ 觉得有帮助也可以投喂一下博主,感谢~
文章链接:https://www.corrain.top/pndiode-simulation2/
版权声明:本博客所有文章除特别声明外,均采用 CC BY-NC-SA 4.0 许可协议。转载请注明文章地址及作者
暂无评论

发送评论 编辑评论


				
|´・ω・)ノ
ヾ(≧∇≦*)ゝ
(☆ω☆)
(╯‵□′)╯︵┴─┴
 ̄﹃ ̄
(/ω\)
∠( ᐛ 」∠)_
(๑•̀ㅁ•́ฅ)
→_→
୧(๑•̀⌄•́๑)૭
٩(ˊᗜˋ*)و
(ノ°ο°)ノ
(´இ皿இ`)
⌇●﹏●⌇
(ฅ´ω`ฅ)
(╯°A°)╯︵○○○
φ( ̄∇ ̄o)
ヾ(´・ ・`。)ノ"
( ง ᵒ̌皿ᵒ̌)ง⁼³₌₃
(ó﹏ò。)
Σ(っ °Д °;)っ
( ,,´・ω・)ノ"(´っω・`。)
╮(╯▽╰)╭
o(*////▽////*)q
>﹏<
( ๑´•ω•) "(ㆆᴗㆆ)
( ゜- ゜)つロ
_(:з」∠)_
(⌒▽⌒)
( ̄▽ ̄)
(=・ω・=)
(*°▽°*)八(*°▽°*)♪
✿ヽ(°▽°)ノ✿
(¦3【▓▓】
눈_눈
(ಡωಡ)
_(≧∇≦」∠)_
━━━∑(゚□゚*川━
(`・ω・´)
( ̄3 ̄)
✧(≖ ◡ ≖✿)
(・∀・)
(〜 ̄△ ̄)〜
→_→
(°∀°)ノ
╮( ̄▽ ̄)╭
( ´_ゝ`)
←_←
(;¬_¬)
(゚Д゚≡゚д゚)!?
( ´・・)ノ(._.`)
Σ(゚д゚;)
Σ(  ̄□ ̄||)<
(´;ω;`)
(/TДT)/
(^・ω・^)
(。・ω・。)
(● ̄(エ) ̄●)
ε=ε=(ノ≧∇≦)ノ
(´・_・`)
(-_-#)
( ̄へ ̄)
( ̄ε(# ̄) Σ
(╯°口°)╯(┴—┴
ヽ(`Д´)ノ
("▔□▔)/
(º﹃º )
(๑>؂<๑)
。゚(゚´Д`)゚。
(∂ω∂)
(┯_┯)
(・ω< )★
( ๑ˊ•̥▵•)੭₎₎
¥ㄟ(´・ᴗ・`)ノ¥
Σ_(꒪ཀ꒪」∠)_
٩(๛ ˘ ³˘)۶❤
(๑‾᷅^‾᷅๑)
😂
😀
😅
😊
🙂
🙃
😌
😍
😘
😜
😝
😏
😒
🙄
😳
😡
😔
😫
😱
😭
💩
👻
🙌
🖕
👍
👫
👬
👭
🌚
🌝
🙈
💊
😶
🙏
🍦
🍉
😣
Source: github.com/k4yt3x/flowerhd
颜文字
Emoji
小恐龙
花!
小黄脸
热词系列一
tv_小电视
上一篇
下一篇