Two basic assumptions implicit in the present algorithm are: 1) the microstructure of interest can be accurately represented as a 2-D or 3-D digital image, where each pixel of the image is treated as a homogeneous single phase, and 2) the impedance of each individual phase in the composite is known at arbitrary frequencies, either experimentally as a list of numbers, or by using a fitted circuit model to represent the single phase properties. In this paper, we show examples of use of the algorithm on models where a simple RC circuit, consisting of a perfect resistor and perfect capacitor placed in parallel, is used to represent single-phase properties. Because of the random geometry in composite materials, complicated composite IS behavior is commonly encountered, even when the single phase properties are given by simple RC circuits.
Figure 1 shows a portion of a two dimensional image, in which electrodes are simulated by placing an extra layer of pixels on opposite sides of the image. The resistances (R) and capacitances (C) for each phase are uniquely determined by that phase's resistivity and dielectric constant, respectively. Different electrical properties may be assigned to the sample-electrode interface, if necessary, to simulate electrode effects.
The R and C values used in the circuits for each pixel are:
| R = (σA/D) −1, | (1) |
| C = krεoA/D | (2) |
where σ is the conductivity of the specified phase, kr = ε/εo is the dielectric constant, εo is the permittivity of free space, D is the length between the center of the pixel and the edge of the pixel (1/2 unit), and A is the area of the pixel face (1 unit2). It is possible to assign a scale to the digital images, e.g. pixel length = 1 mm, but for the purpose of this description each pixel has a length of one unit.
Figure 1: Part of a 2-D image showing how the resistor-capacitor network is mapped onto the digital image of the microstructure.
The admittance 
, which is the reciprocal of impedance, is the AC analog of
conductance. For each RC circuit at each applied frequency, the admittance is then given by:
| ψ = 1/Rj + iωCj, | (3) |
where the subscript (j) refers to the jth phase of the microstructure, R and C are as defined above, ω is the applied angular frequency, and i = (−1)1/2. The sample size is much smaller than the wavelength associated with the applied frequency, so that the same value of ω is used at each pixel. Each pixel has four (2D) or six (3D) impedances extending from its center to the pixel boundaries (see Figure 1). Neighboring pixels are connected together by joining the appropriate impedances into a single bond. Perfect bonds are usually assumed between different phases. However, special internal-interface elements could be easily inserted if necessary.
The end result is the creation of a two or three dimensional electrical network with a node at
the center of each pixel. The solution of Kirchoff's law on this impedance
network is
mathematically equivalent to a finite-difference solution of
Laplace's equation, 
· j
= 2V = 0, where j = E, E = −
V, and
is the admittance [18]. This
equation comes from the static (time-independent) limit of the continuity equation [19]. The
well-known "correspondence principle" relating static (
= 0) elastic problems to single-
frequency viscoelastic problems [20] has its analogy for electrical problems as well, which
we have exploited in developing this algorithm. The mathematics is the same for the single-frequency complex case as for the DC case, but all quantities are allowed to be complex and
therefore dependent on frequency.
The impedance calculation begins by applying a voltage across the microstructure, i.e. the field of pixels. This is accomplished, using a unit voltage step, by setting the voltage of the electrode pixel nodes on one side of the image to 1, and those on the opposing side of the image to 0. Since the applied voltage step is taken to have a zero phase angle, the voltage at the electrodes is real and is not a function of frequency. However, the voltages at the pixel-nodes within the sample will be complex, in general, and are iterated until the net current at each node is zero, satisfying Kirchoff's law [19]. Knowing the voltage at each pixel, it is possible to calculate all the local currents, as well as the total current across the sample, and hence, the composite impedance. A similar technique was employed [8] for DC studies of a three dimensional cement paste model. Since the nodes in this model are connected by impedances instead of only resistors, the network now has a frequency-dependent response.
Periodic boundaries are maintained in directions perpendicular to the applied field by connecting pixels located on a surface of the image to the pixels on the opposite surface. This reduces statistical errors associated with using finite-size digital images. If sample-electrode phenomena are not of interest, then periodicity can also be used in the field direction, while still maintaining the applied field [8]. Otherwise, special interface resistors and capacitors can be inserted between the electrode and the sample pixels to mimic the sample-electrode interfacial response.
Solving Kirchoff's law at each of N nodes leads to an independent set of N linear equations that can be arranged into a matrix of the type Ax = b, where Aij is the admittance between nodes i and j, x is the voltage vector of length N, and b is a constant vector, also of length N, that comes from the constraints in the system. The constraints in this system are periodic boundary conditions perpendicular to the applied field, and the constant electrode voltages. A matrix equation of this type for the nodal voltages lends itself well to a conjugate gradient iteration technique [21]. The iteration process stops when the magnitude of the real and imaginary parts of the residual current going in to a node, averaged over all nodes, are below a prescribed cutoff, chosen low enough to ensure an accurate solution, but not so low that computer time is wasted on achieving more than the necessary accuracy. Since there is one equation and one unknown voltage for each node, a three dimensional, 1003 image, for example, will contain one million equations and unknowns. In general, the matrix A has elements connecting each node to every other node. Fortunately, in a 3-D system, only 7 of the A coefficients for each node are non-zero. These 7 elements come from the six connections from each node to its nearest neighbors to a central site, plus a self term.
Memory and computational speed are the limiting constraints on the size of models that can be used. The memory requirements are on the order of 100 bytes per pixel, when all real numbers are double precision (eight bytes per real number). For each applied frequency, a typical 1003 system is solved in approximately 1600 seconds on a Cray-YMP supercomputer. Usually, a single spectrum requires at least 20 different frequencies to be sampled. Other versions of the conjugate gradient algorithm may give better performance [22].