A Comparison of Three Algorithms for the Inverse Conductivity Problem by Joshua Landrum Basic terminology: The network of resisitors is a square lattice of conductors. For a three by three network, it would be A B C | | | L - 1 - 2 - 3 - D | | | K - 4 - 5 - 6 - E | | | J - 7 - 8 - 9 - F | | | I H G where each number stands for a node where the individual conductors (the line segments) intersect. The boundary wires are listed here by letters; this is done to disitinguish them from the interior nodes. It is more con- venient to think of them as numbered, too, however. Replace the lettering with the obvious numbering. The forward problem: To each node in the lattice a number is assigned, corresponding to an electrical potential. Each of the wires is given a positive number as its conductivity. Using Kirchoff's law (the current going into a node is the same as the current going out of it) and Ohm's law (the current through the wire connecting nodes N and M is equal to the voltage drop from N to M times the conductivity of the wire) a 4nx4n matrix of boundary data, Lambda, is constructed, where Lambda(i,j) = the current out of boundary node i when a voltage of 1 is placed at boundary node j and every other boundary node is grounded. A current into the net- work is considered positive. The inverse problem: The inverse problem is to determine the values of the conductors in the network using the boundary data contained in Lambda. What I shall call the standard algorithm (described by profs. Curtis and Morrow in Determining the Resistors in a Network) works like this: 1 0 0 | | | V1 - . - . - . - 0 (0) | | | V2 - . - . - . - 0 (0) | | | V3 - . - . - . - 0 (0) | | | 0 0 0 where the parenthetical numbers are currents, the other numbers are volt- ages, and V1-V3 are voltages to be solved for. Solving this system of three equations in three unknowns, we can then find the current going into the network at boundary node 1. By inspection, it is clear that V2 and V3 are 0, as is the potential at every interior node. Therefore it is clear that the conductor from boundary node 1 to interior node 1 is equal to the current into the lattice there. Similarly, the conductor from boundary node 4n to interior node 1 is -(current in at boundary node 1)/V1. To get the rest of the conductors, assume all the conductors above the marked diagonal are known: 0 0 0 1 0 0 | | | /x | | | | |/ x | | V6 --.--.--.xx.--.--.-- 0 (0) | | /x | | | | |/ x | | | V5 --.--.xx.--.--.--.-- 0 (0) | /x | | | | |/ x | | | | V4 --.xx.--.--.--.--.-- 0 (0) /x | | | | | / x | | | | | V3 xx.--.--.--.--.--.-- 0 (0) | | | | | | | | | | | | V2 --.--.--.--.--.--.-- 0 (0) | | | | | | | | | | | | V1 --.--.--.--.--.--.-- 0 (0) | | | | | | | | | | | | 0 0 0 0 0 0 Once this system of equations is solved, all the voltages along the diagonal are easy to compute using the two laws, as all the conductors above the diagonal are known. It is also clear from the two laws that the voltage at each node below the diagonal is 0, independent of the value of the conductors, and that no current enters or leaves the network through boundary nodes 5 through 20. Therefore the current passing through the x'd conductors is easy to find; since the voltages at each end of these conductors are known, the conductors can thus be calculated. Due to prob- lems with round-off error, this procedure tends to break down for large networks as the center conductors are calculated, and results lose all sem- blence of accuracy for networks larger than 20x20 for the center 4 resist- ors, even for the case where all resitors are identical. Using a probe to calculate center resistors: If another boundary node is added to the center of the network, more accurate results can be obtained for the central resistors. The new boundary node in the diagram is node 4n+1. The given system is for finding the probe conductor. 0 0 0 0 0 | | | | | | | | | | | | | | | | | | | | V5 ---- . ---- . ---- . ---- . ---- . ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V4 ---- . ---- . ---- . ---- . ---- . ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V3 ---- . ---- . ---- . ---- . ---- . ---- 0 (0) | | |\ | | | | | \ | | | | | 1 | | | | | | | V2 ---- . ---- . ---- . ---- . ---- . ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V1 ---- . ---- . ---- . ---- . ---- . ---- 0 (0) | | | | | | | | | | | | | | | | | | | | 0 0 0 0 0 Applying the two laws, this system yields 0 0 0 0 0 | | | | | | | | | | | | | | | | | | | | V5 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V4 ---- . ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V3 ---- . ---- . ---- 0 ---- 0 ---- 0 ---- 0 (0) | | |\ | | | | | \ | | | | | 1 | | | | | | | V2 ---- . ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V1 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | 0 0 0 0 0 independent of the values of the conductors, where the internal zeros correspond to the voltages at the given internal nodes. After the system of equations is solved, the center conductor is easily found -- it is the current through the probe. To get the center cross of conductors, the following pattern is used (internal voltages which can be determined inde- pendent of the values of the conductors have been entered on this diagram for simplicity; they are not imposed directly. Note that the "determined" voltage V is clearly 0): V V1 V2 V3 0 | | | | | | | | | | | | | | | | | | | | (0) 0 ---- 0 ---- 0 ---- . ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | (0) 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 (0) | | x | | | | x | | | | x | | | | x | | (0) 0 ---- 0 ---- 0 ---- 1 ---- 0 ---- 0 ---- 0 (0) | | |\ | | | | | \ | | | | | 1 | | | | | (0) | | (0) 0 ---- 0 ---- . ---- . ---- . ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | 0 ---- . ---- . ---- . ---- . ---- . ---- 0 | | | | | | | | | | | | | | | | | | | | V8 V7 V6 V5 V4 Noting that the current through the x'd resistor is equal to the net current out of the top edge of the network, and that the voltage drop across it is one, its conductance is easy to find. Although this algorithm does not require knowledge of the central conductor to implement, slightly greater accuracy can be gained by adjusting the current through the central node so as to minimize the voltage swings for V4 through V8. For the four center nodes of a 21x21 network, whose conductors all have a value of one, the difference is: Current in at center Voltage swing (V43-V63) Computed conductivity a) 0 93514131.163309 1.3054195200356 6.625 2909274.1702191 .99921605719529 b) 0 104497107.70718 1.0393972272623 6.75 2905817.5314497 .98666315138441 c) 0 115725094.88017 .90018758095523 6.625 2877081.6031469 1.0033733145505 d) 0 98949969.326176 1.1059649597600 6.625 2851343.7527630 1.0073175724276 The program proceeded by calculating the center conductor using zero current in at the center, then upping the current by 1/8 and calculating the voltage swing as long as the voltage swing continued to decrease. It then used the calculations with the lowest voltage swing to calculate the center conductor again. To get the next ring of resistors, the following algorithm is applied: (0) (0) (0) (0) (0) 0 0 0 0 0 | | | | | | | | | | | | | | | | | | | | 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 | | | | | | | | | | | | | | | | | | | | V8 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- V1 | | " | | | | " | | | | " | | | | " | | V7 ---- 0 ---- 0 ==== 0 ==== 0 ---- . ---- V2 | | "\ | | | | " \ | | | | " 1 | | | | " | | V6 ---- . ---- 0 xxxx V ---- . ---- . ---- V3 | | | | | | | | | | | | | | | | | | | | V5 ---- 0 ---- 0 ---- 0 ---- . ---- . ---- V4 | | | | | | | | | | | | | | | | | | | | 0 0 0 0 0 (0) (0) (0) Again, the interior voltages were "calculated" using the two laws directly; these can be determined indepent of the values of the conductors. Note that V8 is clearly 0. Since the current through the central node is known, and all of it passes through the conductor to the node marked with voltage V, and the ==== resistors have known conductivity, the voltage V can be com- puted. Noting that the current through the x'd coductor is equal to the net current through the boundary nodes at V5-V8, the x'd conductor is quickly calculated. To extend the central cross, the following algorithm can be applied (or, rather, the system of equations associated with it): (0) (0) (0) (0) (0) 0 0 0 0 0 | | | | | | | | | | | | | | | | | | | | 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 | | | | | | | | | | | | | | | | | | | | 0 ---- 0 ---- 0 ==== 0 ==== 0 ---- 0 ---- V1 | " " " | | " " " | | " " " | | " " " | V7 ---- 0 ---- 0 ==== 0 ==== 0 ---- . ---- V2 | " "\ " | | " " \ " | | " " 1 " | | " \/ " " | V6 ---- . ---- 0 ==== V ==== . ---- . ---- V3 | | x | | | | x | | | | x | | | | x | | V5 ---- . ---- . ---- 0 ---- 0 ---- . ---- V4 | | | | | | | | | | | | | | | | | | | | 0 0 0 0 0 (0) (0) Since the current through the probe is known, and all of it passes down to the node marked with voltage V, and all the ==== conductors are known, V can be calculted. Using the voltage V, the current through the conductor marked with the \/ sign can be calculted, since the voltage at each end is known, as is its conductivity. Noting that the current through this conductor plus the current through the x'd conductor must be equal to the net current out bound- ary nodes 15 through 19, the current through the x'd conductor can be found, and thus its conductivity. I have not found any further algorithms to take advantage of the center conductor, though the standard algorithm can be used to finish finding the conductors by keeping the probe insulated. The resistor isolating algorithm: As an alternative to the standard algorithm for finding an arbitrary resistor, what I call the resistor isolat- ing algorithm can be applied: 1 V1 V2 0 0 | | | | | | | | | | | | | | | | | | | | (0) 1 ---- 1 ---- 1 ---- 0 ---- 0 ---- 0 ---- 0 (0) | x | | | | x | | | | x | | | | x | | | 1 ---- . ---- 0 ---- 0 ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | 0 ---- . ---- . ---- 0 ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | 0 ---- . ---- . ---- . ---- 0 ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | 0 ---- . ---- . ---- . ---- . ---- 0 ---- 0 (0) | | | | | | | | | | | | | | | | | | | | V7 V6 V5 V4 V3 Since this algorithm does not need to use the values of any other resistor to find the conductivity desired, it is significantly faster if only a few of the conductors need to be calculted. It is most accurate for finding the central conductors, since the voltage swings are smallest for these calcula- tions; there, it is more accurate than the standard algorithm. Its accuracy rapidly falls off as the edges are approached, however, so with a large net- work resistors in the corner may be computed with stunning innaccuracy using this algorithm. Happily the standard algorith is at its best where this algorithm is at its worst, so a combination of the two could easily be com- bined to get all the resistors in a network with maximum accuracy. A comparison of the accuracy of the three methods: The following numbers are the calculated conductivities of the central conductors of a 15x15 network: . a a a .bbb.ccc. d e f d e f d e f .ggg.hhh.iii.jjj. k l m k l m k l m .nnn.ooo. p p p . The above is the center 16 conductors in the network. In the first set of results, all the conductors in the network had conductivity 1. In the sec- ond, the resistors had sequential conductivity, with the one in the top left corner having conductivity 1, the one to the right of it 2, etc., moving from left to right, top to bottom. The probe was considered to be the last conductor. The standard algorithm results were calculted using a program supplied by professor Curtis. The other two were calculated using programs I wrote, using the linpac DGEFS subroutine. All three were written in FORTRAN77 and used double precision arithmetic. Actual Standard Probe Isolating a 1 1.0003155967396 .99998979795180 1.0000278793221 b 1 1.0001562569931 1.0000194063684 1.0000128308360 c 1 .99978282274008 1.0000195389202 1.0000103954070 d 1 .99975597920674 1.0000201665840 1.0000120226549 e 1 .99951990920973 .99999832005246 .99999559804263 f 1 1.0002059563182 1.0000213209924 1.0000126127107 g 1 1.0004500345630 .99999030346087 1.0000291569481 h 1 .99926468413114 .99999856113445 .99999552317763 i 1 .99928628581162 .99999833608272 .99999575194797 j 1 1.0004379629374 .99998817897260 1.0000297803636 k 1 1.0002118329522 1.0000150106983 1.0000089814608 l 1 .99934342594908 .99999854257274 .99999628247180 m 1 .99977875530022 1.0000152440153 1.0000126342989 n 1 .99978894648519 1.0000074700620 1.0000115076065 o 1 1.0001989009705 1.0000161043047 1.0000116095774 p 1 1.0004157454059 .99998692753540 1.0000273323348 a 194 193.92807643543 193.99741133526 193.99454389909 b 209 208.96368994766 208.99204933756 208.99644612895 c 210 209.98084068038 209.98340422235 209.99787103613 d 224 224.02952444796 223.99031109459 223.99632402739 e 225 225.13885746517 225.00000944652 225.00117766345 f 226 226.06973026993 225.97493247272 225.99989278085 g 239 238.91210420615 238.99509462568 238.98882198080 h 240 240.12395297055 240.00026644241 240.00183743265 i 241 240.71090840956 241.00006952616 241.00064162491 j 242 242.16809319103 242.00281280175 241.99708679208 k 255 254.95321093972 255.00060650720 254.99310471022 l 256 256.42372727799 256.00013724657 256.00116593616 m 257 257.05213817400 257.02422988191 256.99611610070 n 271 271.03992243810 271.00772223235 270.99592362475 o 272 271.91193484390 272.02008489330 272.00272830165 p 287 286.81513196042 286.99806985378 286.99661266536 Stability of the algorithm: The isolating algorithm gets decent results for small networks (on the order of 7x7 or smaller) even when the lambda matrix is rounded to 3 digits. Below, a-p refer to the same pattern of con- ductors as above. Only three of the calculated digits are given, as that is the extent of the accuracy of the input. Resistor Actual 5x5 7x7 9x9 a 1 1.16 2.19 3.25 b 1 1.10 1.06 .695 c 1 1.10 1.06 .695 d 1 1.10 1.06 .695 e 1 1.04 2.02 -3.63d-02 f 1 1.10 1.06 .695 g 1 1.16 2.19 3.25 h 1 1.04 2.02 -3.63d-02 i 1 1.04 2.02 -3.63d-02 j 1 1.16 2.19 3.25 k 1 1.10 1.06 .695 l 1 1.04 2.02 -3.63d-02 m 1 1.10 1.06 .695 n 1 1.10 1.06 .695 o 1 1.10 1.06 .695 p 1 1.16 2.19 3.25 I do not know why the center resistors in the 7x7 network were computed less accurately than most of the other resistors, whereas they are usually com- puted more accurately. In the 9x9 case, I think it can be attributed to coincidence.