引入
我们早就学过了计算串并联电阻的等效电阻的公式:
针对混合的电阻网络,有些时候我们依然可以反复运用这两个公式来计算等效电阻。在有些电阻网络中,这样是无法得到答案的:
如果你熟悉星三角变换 ,你应该能很快算出答案。(或者你能看出这就是一个平衡桥,那就能口算了)
星三角变换固然是一种手工计算等效电阻的很好方法,但是对于计算机而言,利用这个方法计算电阻会很麻烦。下面,我将使用高斯消元法来计算两点之间的等效电阻。
建模
我们可以把两点之间的纯电阻网络看作是一张有向图,节点表示导线,边表示电阻,边权表示电阻的大小。要求源点和汇点之间的等效电阻大小。
为什么要建成有向边呢,因为虽然电阻是无向的,但是电流是有向的,我们一开始并不知道电流的方向,但是可以任意地假定一个方向,然后给每一条边赋上一个电流I i I_i I i ,当I i I_i I i 为正,表示电流方向与箭头同向,当I i I_i I i 为负,表示电流方向与箭头反向。然后我们再假定汇点的电势为0 0 0 ,其余每个节点赋上一个电势φ i \varphi_i φ i 。
若要求源点到汇点的等效电阻,只需要知道源点的电势φ s \varphi_s φ s 和汇点的电势φ t \varphi_t φ t ,以及汇入汇点的总电流I I I ,就可以计算出等效电阻R = φ s − φ t I R = \frac{\varphi_s - \varphi_t}{I} R = I φ s − φ t 。为了计算这些量,我们要把整个网络的电流和电势都算出来。
设这张图共有n n n 节点m m m 条边,则我们一共有m m m 个I i I_i I i 未知,n − 2 n-2 n − 2 个φ i \varphi_i φ i 未知。一共m + n − 2 m+n-2 m + n − 2 个未知量,需要m + n − 2 m+n-2 m + n − 2 个方程才能解出来。我们可以利用所学的物理知识来找到这些方程。
欧姆定律
最简单的就是安培定律了,对于编号为i i i 的边( u , v ) (u, v) ( u , v ) ,有
φ u − φ v = I i R i \begin{align*}
\varphi_u - \varphi_v = I_i R_i\\
\end{align*}
φ u − φ v = I i R i
这一共有m m m 个方程。
基尔霍夫电流定律(KCL)
除了源点和汇点外,每个结点流入的电流都等于流出的电流。
对于编号为u u u 的节点,有
∑ i ∈ A d j ( u ) I i = 0 \sum_{i \in Adj(u)} I_i = 0
i ∈ A d j ( u ) ∑ I i = 0
这 n − 2 n-2 n − 2 个方程,再加上前面的 m m m 个,正好有 m + n − 2 m+n-2 m + n − 2 个方程,我们可以用高斯消元法来解出所有的电流。
高斯消元法
现在不妨令 φ s = 1 \varphi_s = 1 φ s = 1 ,φ t = 0 \varphi_t = 0 φ t = 0 ,则有
{ φ u − φ v = I i R i ∀ ( u , v ) ∈ E ∑ i ∈ A d j ( u ) ± I i = 0 ∀ u ∈ V − { s , t } φ s = 1 φ t = 0 \begin{cases}
\varphi_u - \varphi_v = I_i R_i && \forall (u, v) \in E \\
\sum_{i \in Adj(u)} \pm I_i = 0 && \forall u \in V - \{s, t\} \\
\varphi_s = 1 \\
\varphi_t = 0 \\
\end{cases}
⎩ ⎨ ⎧ φ u − φ v = I i R i ∑ i ∈ A d j ( u ) ± I i = 0 φ s = 1 φ t = 0 ∀ ( u , v ) ∈ E ∀ u ∈ V − { s , t }
有了这 m + n m+n m + n 个方程,我们就可以用高斯消元法解出所有的电势和电流,进而求出等效电阻了。
高斯消元法是一种求解线性方程组的算法,在计算机中的实现也很简单。此处复杂度O ( ( m + n ) 3 ) O((m+n)^3) O (( m + n ) 3 )
此处的优化空间:从上面可以看出,这些方程一般是很稀疏的,如果改用计算稀疏线性方程组的算法可以达到更高的效率。如LU分解、bfgs、lbfgs
代码实现
使用纯C实现
Github地址:https://github.com/cyrus28214/resistance
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 #include <stdio.h> #include <stdlib.h> #include <math.h> double **array2d (int rows, int cols) ;void gauss (double **a, int n) ;int main () { int n, m, s, t; scanf ("%d%d%d%d" , &n, &m, &s, &t); int k = n + m + 1 ; double **a = array2d(k, k+1 ); for (int i = 0 ; i < m; i++) { int u, v; double r; scanf ("%d%d%lf" , &u, &v, &r); a[i][m+u] = 1 ; a[i][m+v] = -1 ; a[i][i] = -r; if (u != s && u != t) a[m+u][i] = -1 ; else if (u == t) a[m+n][i] = -1 ; if (v != s && v != t) a[m+v][i] = 1 ; else if (v == t) a[m+n][i] = 1 ; } a[m+s][m+s] = 1 , a[m+s][k] = 1 ; a[m+t][m+t] = 1 , a[m+t][k] = 0 ; a[m+n][m+n] = -1 ; gauss(a, k); double I_total = a[k-1 ][k]; printf ("%6e\n" , 1 / I_total); } double **array2d (int rows, int cols) { double **arr = malloc (rows * sizeof (double *)); for (int i = 0 ; i < rows; i++) { arr[i] = calloc (cols, sizeof (double )); } return arr; } void gauss (double **a, int n) { for (int i = 0 ; i < n; i++) { for (int j = i + 1 ; j < n; j++) { if (fabs (a[j][i]) > fabs (a[i][i])) { double *temp = a[i]; a[i] = a[j]; a[j] = temp; } } double c = a[i][i]; for (int k = 0 ; k <= n; k++) a[i][k] /= c; for (int j = 0 ; j < n; j++) { if (i == j) continue ; c = a[j][i]; for (int k = 0 ; k <= n; k++) a[j][k] -= c * a[i][k]; } } }
验证
对于开头的这个例子,我们可以验证一下:
输入
1 2 3 4 5 6 7 4 5 0 1 0 2 20 0 3 20 1 2 60 1 3 60 3 2 50 4.000000e+01
输出4.000000e+01
,和手算结果一致
再来一个例子,经典高中物理题,每个电阻都是1 Ω 1\Omega 1Ω ,求对角线电阻:
用我们的程序算出:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 8 12 0 6 0 1 1 1 2 1 2 3 1 3 0 1 4 5 1 5 6 1 6 7 1 7 4 1 0 4 1 1 5 1 2 6 1 3 7 1 8.333333e-01
而这一题的答案是5 / 6 Ω 5 / 6 \,\Omega 5/6 Ω ,答案正确。
一些私货
qzqz机房的hxd们可以在cofun上看到我出的这一题,快去做(cofun服务器好像又爆炸了)
如果你不知道cofun是什么,请忽略这段私货。