Substitution of subspace collections with nonorthogonal subspaces to accelerate Fast Fourier Transform methods applied to conducting composites

01/22/2020
by   Graeme W. Milton, et al.
0

We show the power of the algebra of subspace collections developed in Chapter 7 of the book "Extending the Theory of Composites to Other Areas of Science (edited by Milton, 2016). Specifically we accelerate the Fast Fourier Transform schemes of Moulinec and Suquet and Eyre and Milton (1994, 1998) for computing the fields and effective tensor in a conducting periodic medium by substituting a subspace collection with nonorthogonal subspaces inside one with orthogonal subspaces. This can be done when the effective conductivity as a function of the conductivity σ_1 of the inclusion phase (with the matrix phase conductivity set to 1) has its singularities confined to an interval [-β,-α] of the negative real σ_1 axis. Numerical results of Moulinec and Suquet show accelerated convergence for the model example of a square array of squares at 25% volume fraction. For other problems we show how Q^*_C-convex functions can be used to restrict the region where singularities of the effective tensor as a function of the component tensors might be found.

READ FULL TEXT VIEW PDF
POST COMMENT

Comments

There are no comments yet.

Authors

page 11

06/05/2013

Quaternion Fourier Transform on Quaternion Fields and Generalizations

We treat the quaternionic Fourier transform (QFT) applied to quaternion ...
07/17/2019

Interesting Open Problem Related to Complexity of Computing the Fourier Transform and Group Theory

The Fourier Transform is one of the most important linear transformation...
07/02/2019

Robust Tensor Completion Using Transformed Tensor SVD

In this paper, we study robust tensor completion by using transformed te...
09/11/2018

A fast Fourier transform based direct solver for the Helmholtz problem

This paper is devoted to the efficient numerical solution of the Helmhol...
07/20/2018

Fast transforms over finite fields of characteristic two

An additive fast Fourier transform over a finite field of characteristic...
09/25/2020

On the Compression of Translation Operator Tensors in FMM-FFT-Accelerated SIE Simulators via Tensor Decompositions

Tensor decomposition methodologies are proposed to reduce the memory req...
08/19/2018

Fourier analysis perspective for sufficient dimension reduction problem

A theory of sufficient dimension reduction (SDR) is developed from an op...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

Subspace collections with orthogonal subspaces obviously have a lot of relevance to physical problems. Here we give an example to show that subspace collections with nonorthogonal subspaces are an important tool in analysis, reaching beyond their connection with rational functions of several complex variables. This chapter is mostly self-contained. While it uses nonorthogonal subspace collections, it is not necessary for the reader to have read Chapter 7 of [24]. However, it is recommended for the reader to look at Chapters 1 and 2 in that book.

The specific problem we consider is a two-component conducting medium with conductivities and , phase being the inclusion phase. The dimension is or . Without loss of generality we can take . Assume for simplicity that the composite has square or cubic symmetry so that the effective conductivity tensor is isotropic, taking the form . [2] made the pioneering observation that the function is an analytic function of and has all its singularities on the negative real axis. He made some assumptions which were not valid [17]. These assumptions could be circumvented by approximating the composite by a large resistor network ([18]). Subsequently Golden and Papanicolaou ([9], [10]) gave a rigorous proof of the analytic properties. Here we make the additional assumption (only true for some geometries) that the function has all its singularities confined to an interval on the negative real axis, and (see 1). This information results in tighter bounds on the conductivity and complex conductivity, or equivalently the complex dielectric constant, ([4]; [35]; [8]), and these bounds may be inverted to yield information about and from experimental measurements ([31]). By contrast, here our goal is to utilize the information to improve the speed of convergence of Fast Fourier Transform methods for computing the fields in composites and their associated effective tensors. These methods were first introduced by Moulinec and Suquet ([27], [28]): see [26] for a recent review. One important application has been to viscoplastic polycrystals ([14], [15]). For materials with a linear response, previous significant advances in the acceleration of these schemes were made by [7] (see also the generalization in section 14.9 of [21], and in particular equation (14.38)), and by [39] and [38].

Figure 1: We assume that in the complex -plane, the function has all its singularities on the negative real -axis, between and . Here the zig-zag line denotes a possible branch cut.

The singularities in the interval could be poles or could be a branch cut (if the inclusion has sharp corners, or if there is some randomness in the geometry). Note it is not only the effective conductivity which has this analytic form but also the electric field and current field at fixed , and fixed applied field, with varying. If we didn’t know anything about and there could potentially be singularities anywhere along the negative real axis. Let’s first look at the case , .

The original Fast Fourier Transform scheme of Moulinec and Suquet ([27], [28]) is based on series expansions, such as (2.35) in [24], which for conductivity takes the form

(1.1)

for the effective conductivity , and electric field , with being the applied (average) electric field. [To obtain (1.1) from (2.35) in [24] we use the fact that .] Their key and beautiful idea was that since the action of is readily evaluated in Fourier space, while the action of is readily evaluated in real space, both and can be readily calculated from these series by going back and forth between real and Fourier space, using Fast Fourier Transforms to do so.

In a two-phase medium with isotropic conductivities and , and with the choice of used by Moulinec and Suquet, takes the value in phase 1 and the value in phase 2. So it is clear that this scheme gives an expansion of the form

(1.2)

The transformation maps the right half of the complex plane to the unit disk (see 2).

Figure 2: If we want a series expansion which converges in the entire right half of the -plane, when a singularity is located at the point

, then we make a fractional linear transformation which takes the right half of the

-plane to the unit disk in the -plane, and find an series expansion in powers of . The scheme of Moulinec and Suquet ([27], [28]) provides such an expansion.

So with the series will converge if , i.e., in the entire right half of the complex plane, . (Recall that it is the distance from the origin to the nearest singularity in the plane which determines the radius of convergence in the plane, and hence the region of convergence in the plane.)

With the accelerated scheme of [7] (see also the generalization in section 14.9 of [21], and in particular equation (14.38)), one has the expansion

(1.3)

The transformation where maps the complex plane minus the slit along the negative real axis to the unit disk (see 3).

Figure 3: If we want a series expansion which converges in the entire the -plane minus the negative real -axis, then we first make a square root transformation which maps the cut complex -plane to the right half of the -plane, followed by a fractional linear transformation which takes it to the unit disk in the -plane, and find an expansion in powers of . The scheme of [7] provides such an expansion.

So with the series will converge if , i.e., in the entire complex plane minus the slit along the negative real axis. These arguments show that when the accelerated scheme should have a larger region of convergence, and by the same line of reasoning a faster rate of convergence when both schemes converge. (Note, however, that if and the scheme of Moulinec and Suquet could outperform that of Eyre and Milton for small and very large values of , since the scheme of Moulinec and Suquet should then converge for sufficiently small or sufficiently large negative real values of .)

If we know and (or bounds for them) then it makes sense to use a transformation which maps the complex plane minus the slit to the unit circle . Since the transformation

(1.4)

maps the interval to and takes to , it is clear that the transformation

(1.5)

maps the complex plane minus the slit to the unit disk (see 4).

Figure 4: If we want a series expansion which converges in the entire the -plane minus the cut on -axis between and then we first use a fractional linear transformation which maps this cut to the entire negative real axis in the -plane, followed by a square root transformation mapping the cut complex -plane to the right half of the -plane, followed by a fractional linear transformation which takes it to the unit disk in the -plane, and find an expansion in powers of . The scheme developed in this chapter provides such an expansion.

So we want to get a Fast Fourier Transform (FFT) method associated with an expansion

(1.6)

It looks rather formidable but really it is just a matter of substituting

(1.7)

in the [7] scheme associated with

(1.8)

But of course one wants to do this substitution at the level of the underlying Hilbert space, not just in the conductivity function.

2 Substitution at the level of the Hilbert space

Substitution at the level of the Hilbert space requires the substitution of subspace collections where the subspaces do not satisfy the orthogonality property. At a discrete level it is easy to get an idea of how this can be done. If we consider the composite as a resistor network, in phase with resistors having resistance , while in phase with resistors having resistance , so e.g., the composite is replaced by a network (see 5(a)).

Figure 5: A composite containing two phases, as for example a periodic array of disks with a unit cell as shown in (a), can be approximated by a discrete network of resistors having resistors and . At this discrete level we are free to replace each resistor having resistance by the combination of resistors in figure (b). This transformation at the Hilbert-space level, corresponds to a fractional linear transformation in the -plane.

Then at this level we could replace each resistor by the compound resistor (see 5(b)) with real constants. That is the resistance gets replaced by

(2.1)

when , and so gets replaced by

(2.2)

which is a fractional linear transformation of . In other words, in phase the number of fields is multiplied by (replacing one resistor by ). Note that with real positive values of , and the transformation (2.2) maps the nonnegative real -axis onto an interval on the positive real axis, and the negative real -axis onto its complement, which cannot be the desired interval . Essentially, instead of making a transformation which has the effect of lengthening the branch cut as desired in the map at the top of 4, we have a transformation which presumably (in some suitably defined metric) shortens the branch cut. This is why it is necessary to substitute nonorthogonal subspace collections, rather than orthogonal ones [see also the discussion below (4.17)], since, as we will see, they only act to lengthen the branch cut.

How to do this in general in the continuum case is described in section 29.1 of [21], bottom of page 621 and page 622, for orthogonal subspace collections, and in the second half of Section 7.8 of [24], without assuming orthogonality. That is a bit abstract so let us go through it for the case in question.

3 The original subspace collection

Our starting point is the Hilbert space of fields (real or complex

dimensional vector fields), that are cell-periodic and square-integrable in the unit cell. We denote

as the projection onto all fields which are zero in phase , as the space of constant vector fields, as the space of gradients which derive from periodic potentials (i.e., if and ) and as the space of divergence free fields with zero average value (i.e., if and ). In our Hilbert space the inner product is taken to be

(3.3)

where the overline denotes complex conjugation. With respect to this inner product the spaces and are orthogonal. The projection onto we denote by . In Fourier space

(3.4)

The field equations are solved once we have found electric fields and current fields such that

(3.5)

The effective tensor is then defined via

(3.6)

where is the projection onto . We assume the composite is isotropic so that is a scalar, although the method is easily generalized to anisotropic materials where the effective conductivity is a second order tensor.

4 The vector subspace collection we substitute into the original subspace collection

Now consider a dimensional subspace collection consisting of component vectors with inner product

(4.7)

where the overline denotes complex conjugation. The projection projects onto the one dimensional space of fields proportional to the unit vector where and are given constants such that . The ’s could be complex but we do not mean . Thus is a projection but not an orthogonal projection when the ’s are complex, as then is not Hermitian. We take the following:

(4.8)

The field equations become

(4.9)

where the constant will be chosen so the associated “effective modulus” is . That is

(4.10)

where is the projection onto , so that

(4.11)

Without loss of generality we can choose so the field equations become

(4.12)

From the middle equation we get

(4.13)

which with gives

(4.14)

So we have

(4.15)

which with is satisfied with

(4.16)

where

(4.17)

Note that (respectively ) is obtained by substituting (respectively ) in (4.16). Given real we need to choose and so that these equations are satisfied. This will necessitate complex solutions for and since otherwise will be negative. Note that with and being complex, is no longer Hermitian, even though it is a (non self-adjoint) projection, and so we have a subspace collection which is not orthogonal: is not orthogonal to . Also from the field equations with we get

(4.18)

i.e.,

(4.19)

5 The subspace collection after the substitution

Now consider the Hilbert space consisting of all periodic fields of the form

(5.20)

Fields in take the form

(5.21)

Fields in take the form

(5.22)

where . Fields in take the form

(5.23)

where . The space consists of all vectors of the form

(5.24)

and consists of all fields of the form

(5.25)

Also consists of all vectors of the form

(5.26)

and consists of all fields of the form

(5.27)

The inner product on is defined to be

(5.28)

where the overline denotes taking a complex conjugate. With this inner product, the subspaces and are mutually orthogonal. We define , i.e.,

(5.29)

where is the identity matrix which defines even if and are complex. Now the field equations are

(5.30)

These are easy to solve given periodic solutions and to the equations in the Hilbert space , i.e.,

(5.31)

We take (with )

(5.32)

Note that we have and . Also, with , we have

(5.33)

Finally if is the projection onto we have

(5.34)

and since we deduce that

(5.35)

That is, is still the effective tensor. Now the idea is to apply, either the basic Fast Fourier Transform scheme of Moulinec and Suquet ([27], [28]) to the Hilbert space that is associated with the expansion

(5.36)

or the accelerated Fast Fourier Transform method of [7] as generalized in section 14.9 of [21] to the Hilbert space that is associated with the expansion

(5.37)

The operator is easily evaluated in real space. The operator which projects onto is easily evaluated in Fourier space since

(5.38)

where in Fourier space

(5.39)

Hence the accelerated Fast Fourier Transform method of Eyre and Milton can be directly applied in the Hilbert space .

6 Proof of acceleration

Let us suppose is fixed such that . The rate of convergence will be determined by the magnitude of

(6.1)

Since , when , the convergence will be quicker the smaller is. Now as increases,

(6.2)

decreases monotonically from the value at , to the value at . So the larger the value of , the faster the rate of convergence. Similarly as decreases, decreases monotonically from the value at to the value of at . So the smaller the value of , the faster the rate of convergence. More generally when is complex, to see the rate of convergence one should plot the contours in the complex plane. It might be instructive to do this for particular values of and , and compare the contours with .

7 The numerical example of Moulinec and Suquet

This numerical example is due to Hervé Moulinec and Pierre Suquet (to be published) and compares in a model example the speeds of convergence of results in the Hilbert spaces and for the Fast Fourier Transform scheme first proposed by Moulinec and Suquet ([27], [28]). Then the speeds of convergence of results are compared in the Hilbert spaces and for the accelerated scheme of [7] [see also the generalization in section 14.9 of [21]]. In both cases the Fast Fourier Transform schemes converge substantially faster in the Hilbert space .

It is to be emphasized that while these numerical results are only for the effective tensor, the real interest in the new algorithm is for obtaining results for the fields inside the body: accelerated rates of convergence for the effective conductivity alone could easily be obtained by using Padé approximants (

[1]) or (almost equivalently) the associated bounds which use series expansion coefficients up to a given order ([19]; [16]; [23]: see also Chapter 27 in [21] and Chapter 10 in [24], and references therein). Using the new method I expect there will be a similar acceleration of the convergence rates for the fields. Indeed, bounds on the norm of the difference between the actual field and the field obtained by truncating the series expansion show these improved convergence rates. However this does not guarantee pointwise convergence, and in any case it needs to be numerically explored. Padé approximants methods could also be used for the fields, but this approach has the disadvantage that one needs to simultaneously store a lot of information: not just the fields, but also the fields that appear up to a given order in the perturbation expansion for a nearly homogeneous medium.

The model example is a regular array of squares of conductivity occupying a volume fraction of 25% in a matrix of conductivity , as illustrated in 6.

Figure 6: A periodic array of squares at a volume fraction of 25%, as illustrated here, provides a benchmark for testing the theory as, thanks to [30], there is an exact formula for its effective conductivity, given in (7.1).

An exact formula for the effective conductivity of this array (and for the interior fields) was discovered by [30],

(7.1)

and clearly has a branch cut between and . [Interestingly, for the four-phase checkerboard [25] conjectured a formula, that was later independently proved by [5] and [20]

.] Taking a wider estimate for the branch cut (with

and ), Moulinec and Suquet used the algorithm described in this chapter and found that the acceleration provided by the new method was generally substantially improved as shown in their 7 and 8.

Figure 7: The results of Moulinec and Suquet (to be published) comparing the convergence of numerical results to the formula (7.1) for the effective conductivity when and , for their basic scheme (Moulinec and Suquet [27], [28]), and for the basic scheme applied in the Hilbert space . Printed with permission of Hervé Moulinec and Pierre Suquet.
Figure 8: The results of Moulinec and Suquet (to be published) comparing the convergence of numerical results to the formula (7.1) for the effective conductivity , when and for the accelerated scheme of [7] (see also the generalization in section 14.9 of [21]), and for the accelerated scheme applied in the Hilbert space . Printed with permission of Hervé Moulinec and Pierre Suquet

8 Estimating the parameters and

[4] has derived rigorous lower bounds on and upper bounds on , for suspensions of separated spheres in a medium. Actually we need only estimates on and . For instance if an estimate on is slightly too large, then the transformations will look like those depicted in 9. The series will still converge but only for , where will be close to (but less than ). However, the support of the measure will be highly dependent on small details of the microstructure: a tiny sharp cusp on the surface of an otherwise smooth inclusion will in general dramatically change the support of the measure: see [13], page 378 of [21], and equation (15) of [12]. As remarked on page 378 of [21], this is related to the fact that these sharp corners behave as sinks of energy in the mathematically equivalent dielectric problem with the conductivities and being replaced by electrical permittivities and , where is real and positive, while is almost real and negative with a tiny positive imaginary part. Related observations and analysis include those of [34], [33], [6], and [3]. Also there is a close connection with the essential spectrum of the Neumann–Poincaré operator on planar domains with corners ([11]; [32]).

Figure 9: The new accelerated scheme should still work although not as effectively, even when we make a small error in the estimation of (and/or ). With this error a singularity at , gets mapped to a point slightly greater than zero in the -plane, and then gets mapped to a point slightly inside the unit disk in the -plane. There is a consequent reduction in the radius of convergence, and an associated decrease in the rate of convergence, for the series in powers of .

9 Bounds on the support of the measure using -convex functions

For many other equations in periodic composite materials robust things can be said about the measure which restrict its support. We assume we have a Hilbert space of square integrable periodic fields with an inner product

(9.1)

where the angular brackets denote an average over the unit cell of periodicity. We also assume has a decomposition into three orthogonal subspaces

(9.2)

We are interested in the equations

(9.3)

where

(9.4)

and where the local tensor field takes the form

(9.5)

in which

is the orthogonal matrix (satisfying

) associated with a rotation acting on elements in the tensor space, while is a field of rotation matrices giving the local orientation of each phase, and is the indicator function that is in phase and zero elsewhere.

First, following section 14.8 of [21] let us show that there are series expansions for the effective tensor that converge even when the local tensor is not everywhere positive definite. To obtain the series expansion we introduce a constant reference tensor and an associated operator : we say

(9.6)

Explicitly, is given by the operator

(9.7)

where the inverse is to be taken on the subspace . Introducing the polarization field

(9.8)

we see it satisfies , and hence we have the identity

(9.9)

which gives

(9.10)

Averaging both sides yields a formula for the effective tensor, and associated series expansion:

(9.11)

where the angular brackets denote a volume average.

A sufficient condition for convergence is that the operator has norm less than . We take a reference tensor of the form

(9.12)

where the tensor is assumed to be bounded and self-adjoint, but need not be positive definite. We are interested in seeing whether the expansion for the effective tensor converges when is very large. Expanding the operator in powers of gives

(9.13)

with a remainder term

(9.14)

When this remainder term has norm satisfying the bound