Next: Intrinsic Viscosity and Up: Main Previous: Intrinsic Viscosity and

Numerical Investigation of the Ratio in d=3

Further examples of the approximate invariance of [η] / [σ]_∞ for a variety of shapes are given in this section based on numerical (finite element) computations in combination with partial analytic results for [η] and [σ]_∞. All of the results obtained are consistent with (36).

The analogy of the elastostatic and hydrodynamic problems of fluid suspensions and solid composites [17,102,103], mentioned in Sec. 3, indicates that a modification of existing finite element programs for calculating the effective elastic properties of composite bodies can be made to also obtain η. This modification and the variational principle for obtaining Stokes' equation on which it is based, are described in Appendix E for particles with orthorhombic symmetry or higher (triaxial ellipsoids have orthorhombic symmetry). We note that Brenner's work [117] was essential in checking the consistency of this generalization, especially in the case of anisotropic elastic stiffness and viscosity tensors. We also utilize a similar finite element program for the calculation of [σ]_∞ . This finite element method is also described in Appendix E. All particles were represented by a cubic digital image, so that the elements were cubes arrayed on a simple cubic lattice. A standard lattice of size 104³ was used, which was the largest that would fit in the memory of the computer available to us and which would allow reasonable running times. Even so, the total CPU time used to compute the results in this paper was about 2000 hours on a CONVEX 3820 supercomputer.

In these numerical calculations, arbitrary shapes had to be represented by collections of pixels. Because of the overall computational cell size limit, a compromise had to be taken between using enough pixels to give a good representation of the particle, and keeping the particle small compared to the overall unit cell, so as to keep the volume fraction small enough to be in the linear regime in concentration. The size and complexity of the objects that could be treated in this fashion is necessarily limited, but a good approximation to a wide range of physically interesting objects could still be obtained.

Periodic boundaries were used in all simulations to reduce the importance of finite size edge effects. Since a cubic cell was always used, in reality all computations were really for simple cubic periodic arrays of the object considered. Exact calculations exist for the intrinsic conductivity and viscosity of rigid spheres arranged on a simple cubic lattice. This example can then be used to illustrate the effect of finite resolution, as described above, and finite system size on the accuracy of the computations.

Zuzovsky and Brenner carried out computations for the effective conductivity of cubic arrays of spheres embedded in a matrix [121], which are very useful for comparison with our numerical data. For the particular case where the spheres were 'superconducting' and the matrix was an ordinary conductor of unit conductivity, they developed an accurate formula for the effective conductivity σ of the composite medium. Subtracting one from the effective conductivity, and dividing by the sphere volume fraction φ, gives their prediction for the 'intrinsic conductivity' at any sphere volume fraction:

where φ = π (d/L)³ / 6 , d = sphere diameter, and L = size of cubic unit cell. Actually , this quantity is only equal to the true intrinsic conductivity in the limit where φ is small enough so that the expansion in eqn. (3) is applicable. Eqn. (45), however, provides a useful way to represent our numerical conductivity data.

Nunan and Keller [121] have computed the components of the viscosity tensor of the simple cubic array of rigid spheres in a fluid. There are two independent components for this symmetry (there would be three independent elastic components, but incompressibility reduces these to two), defined by Nunan and Keller as two functions of φ, p and q. Using their exact numerical results, they were able to show that an analytic expansion given by Zuzovsky [121] was accurate to within 0.2% up to φ = 0.13 (d/L = 0.63) for simple cubic sphere packings. This analytic expression equals,

where a = 0.2857 and b = −0.04655. In terms of p and q, the rotationally averaged intrinsic viscosity, [η] = (η / η_o − 1) / φ, at any volume fraction φ is given by,

Fig. 6 shows the finite element results for periodic arrays of spheres, along with the exact results, (45) and (48). Consider first the results for the intrinsic viscosity (circles). At small values of d/L, the numerical results are well above the exact result. This is due to not having enough pixels to represent the spherical shape. For example, a sphere with a diameter of 5 pixels (a pixel is considered to be part of the sphere if its center lies within a radius of the center) does not look much like a smooth continuum sphere. In fact, all the finite element results shown in Fig. 6 have been rotationally averaged, as they have cubic symmetry. As d/L increases, allowing each sphere to be represented by more pixels, resolution improves, and the numerical points approach the exact curve. There is a region, around d/L = 0.4, where the numerical results are essentially exact, but the value of [η] has not changed much from the d / L → 0⁺ limit. This is the region in which we have tried to run all the simulations: d/L high enough to give good particle shape resolution, but not high enough so that φ is out of the linear regime. Obviously for higher aspect ratio particles, it is harder to stay in this range of d/L.

Now consider the results for [σ]_∞ shown in Fig. 6. The comparison between exact formula and numerical (finite element) results is similar, except that the numerical results are consistently about 5-8% above the exact curve. This error is possibly larger for larger aspect ratio particles.

Judging by the results for spheres shown in Fig. 6, we can expect that both the numerically computed intrinsic conductivity and viscosity will be systematically high, with the ratio [η] / [σ]_∞ probably somewhat low, as the intrinsic conductivity seems to overshoot the true result slightly more than the intrinsic viscosity. Slightly increasing all the ratios involving a numerical computation of [σ]_∞ in Table 9, Table 10, and Table 11 would improve the agreement with (29).

We have run other tests on the intrinsic viscosity using non-spherical shapes. Having exact results for non-trivial shapes is crucial in order to check the accuracy of numerical simulations. We have computed [η] for a spherical dumbbell with r_p = 1.526, giving [η] ≈ 4.9 , as compared to the exact value [η] = 4.89 (r_p = 1.5431) . An ellipsoid of revolution with an aspect ratio of 3 gave [η] = 3.91 which is 6.1% higher than the exact value of 3.685. The numerical result for two touching spheres was 3.62, about 5% larger than the exact result of [η] = 3.45 given in eqn. (25).

The important physical example of a right circular cylinder is next considered. Very precise analytical calculations of the polarizability α_e have been made for the cylinder [122,122]. Values of [σ]_∞ calculated from these results are given in Table 9. Finite element calculations of [η] and [σ]_∞ for several aspect ratios, which have comparable accuracy to the sphere and touching sphere test cases, are also given in Table 9. The ratio [η]/ [σ]_∞ obtained from this combination of numerical and analytical calculations accords rather well with eqn. (29) and the general correlation (36). Note the difference between the columns marked "Numerical/Exact" and "Numerical/Numerical" giving the results for the ratio [η]/ [σ]_∞. ( The term 'numerical' refers to an estimate obtained by finite element calculation while 'exact' refers to analytic results. The 'exact' results often involve a non-trivial numerical evaluation of the integral expressions which define the analytic results, however.) Using the finite element estimates for the intrinsic conductivity instead of the exact results gave a value somewhat closer to the prediction of (29), since similar systematic computational errors for [η and [σ]_∞ probably compensate.

We next examine the case of rectangular parallelepipeds which is summarized in Table 10. Simulation results closely parallel the exact analytical calculations for the ellipsoid of revolution case discussed in Sec. 3. Again the ratio [η] / [σ]_∞ is shown to be nearly invariant with respect to shape. The oblate result for a very large aspect ratio (marked "Infinity" in Table 10) corresponds to the case of a square plate and such objects are discussed more fully in Sect. 6. The value for [σ]_∞ is approximately 7% higher than the best known experimental value (see Sec. 6) and we expect that our estimate of [η] is too large by about the same amount. The [η]/ [σ]_∞ ratio tends to decrease as aspect ratio increases. These results parallel the analytic results for ellipsoids of revolution in the prolate and oblate limits.

Next, we illustrate a simple means to increase [η] and [σ]_∞ to large values without making a very extended or flat object. We consider a cube of unit edge length, in which a square 'channel' is cut through the center of each face, which passes completely through the cube. The parameter m is taken to be the edge length of the cutout face in units of the cube edge length. We obtain a rigid cubic wire frame when m approaches 1. Notice that cutting out the center, which makes the particle more 'sponge-like', has a very large effect on [η] and [σ]_∞, as can be seen in Table 11. It would be interesting to push the effect to the extreme in a different way by generating a Menger sponge [123] 'fractal' by a repeated decimation of the cube at different scales so that [η] and [σ]_∞ would diverge in a characteristic fashion related to the fractal dimension of the 'sponge'. The memory capacity of the computer was not large enough to allow us to consider more than one or two generations of such an iteratively constructed 'diffuse' object, so we presently confine ourselves to the first generation wire frame structure as described above. The rapid increase of [η] and [σ]_∞ when large holes are cut out is noted. In the limit where m goes to one the intrinsic conductivity and viscosity appear to scale roughly quadratically in (1−m).

In a similar vein we also consider a flat square 'ring' where each side has length 25 and the ring has a unit cross-section. The effect on the virials is large, as in the 'sponge' case (see Table 11). Results for a square cross-section tube with a hole of width 1/3 of the side length are given in Table 11 and a less pronounced effect is found.

We consider other strategies of particle modification in Table 11. For example, instead of decimating the structure we introduce protuberances onto our object. Specifically, we poke three rectangular parallelipipeds orthogonally through a sphere to create a 'jack-like' structure. Taking the width of these parallelipipeds as a unit of length, we take their length as 15 units and the sphere diameter as 9. The increase of [η] and [σ]_∞ is not as dramatic as the 'sponge' case, but the effect is still appreciable. We next consider the case of separated and aligned cubes ('dice'), connected by a rigid, conducting wire of vanishing thickness to maintain object connectivity. The results in Table 11 show that such a tethering of 'bulky' groups gives a very large effect on the values of [η] and [σ]_∞. The values of [η] for 'dice' are similar to [η] values for the spherical particle dumbbell (see Table 7) at comparable separations r_p.