## Abstract

We study the problem of optimal design of a two-phase composite material to maximize the sum of thermal and electrical conductivity. We obtain a characterization of an optimal configuration in terms of a free boundary type H̅ problem. In the process of obtaining this result, we provide a new proof of the relevant case of Bergman's cross-property bound. We use our characterization to argue that the optimal interface is not, as has been suggested, a periodic minimal surface.

## 1. Introduction

In this paper, we consider periodic two-phase composites. If a composite is made of two homogeneous materials with thermal (or electrical) conductivities *σ*_{1} and *σ*_{2}, then the effective conductivity *σ*_{e} of the composite depends on the microstructure and can be computed by solving the so-called *cell* problem (see §2).

Let *ϕ*_{1} be the proportion of the *n*-dimensional unit cube covered by the first material and *ϕ*_{2} the proportion covered by the second. Assume that the configuration is cubically symmetric, so that the resulting composite is isotropic. When the thermal conductivities of the two materials are *σ*_{1} and *σ*_{2} and the corresponding electrical conductivities are *λ*_{1} and *λ*_{2}, Bergman (1978) derived the following cross-property bound:where *λ*_{e} is the effective electrical conductivity and 〈*λ*〉=*ϕ*_{1}*λ*_{1}+*ϕ*_{2}*λ*_{2}. From the above equation, the value of *σ*^{*} is obtained and gives an upper bound for *σ*_{e} as long as (*λ*_{2}*σ*_{1}−*λ*_{1}*σ*_{2})/(*σ*_{1}−*σ*_{2})>0. In particular, if *σ*_{1}=1, *σ*_{2}=*ϵ*, *λ*_{1}=*ϵ*, *λ*_{2}=1, and *ϕ*_{1}=*ϕ*_{2}=1/2, Bergman's upper bound says that(1.1)

A question left open by Bergman's work is whether or not this upper bound is achievable. In recent works by Torquato *et al*. (2002, 2003) two microstructures where the phases are separated by a minimal surface (either the Schwartz *P* surface (figure 1) or the diamond *D* surface) were tested numerically. Both examples were shown to achieve the upper bound (1.1) within three digits of accuracy when *ϵ*=0.1. Based on this computation, it was natural to hypothesize that these structures achieve equality in Bergman's bound (1.1). However, it is not clear whether this microstructure is truly optimal or just produces an effective conductivity that is so close to optimal that it cannot be distinguished by the numerical approximation. The same structure was shown by Torquato & Donev (2004) to give at least a very close to optimal value for the sum of the electrical (or thermal) conductivity and the bulk modulus.

In this paper, we find a new proof of the relevant case of Bergman's upper bound (1.1). Our proof is very simple and elementary, which allows us to derive explicit optimality conditions. If the optimal structure has a minimal surface as its interface, our conditions have implications that can be easily tested numerically. We report some tests for the Schwartz *P* surface which strongly suggest that it is *not* an optimal structure.

Whether or not the Schwartz *P* surface gives an optimal effective conductivity may not be interesting for most practical engineering purposes if the difference is less than 0.001. However, it is interesting to answer this question in order to understand the theory better and, in particular, to learn whether there is a connection (as suggested by Torquato *et al*. (2002, 2003) and Torquato & Donev 2004) between *multifunctionality* and the mean curvature of the interface.

The structure of the paper is as follows. In §2, we set up the problem and notation in a precise way. In §3, we present a new proof of (1.1) and provide explicit conditions for equality. In §4, we consider the possibility that an interface would satisfy our optimality condition and also have mean curvature zero. We derive some consequences (in any space dimension) and test them numerically for the Schwartz *P* surface (in three space dimensions). The numerical results indicate that the Schwartz *P* surface cannot be exactly optimal. Finally, in appendix A we show a new proof of the well-known Hashin–Shtrikman bounds.

The new proof of (1.1) that we present in §3 is based on the observation that in the optimal configurations the solutions of the cell problem agree with the directional derivatives of a potential function. This idea may have some interest in itself. It was motivated by looking at the known proofs of the Hashin–Shtrikman bounds and especially at the constructions that achieve the bounds. After this observation is made, the proof comes quickly after very straightforward elementary computations. The same idea can be used to prove the classical Hashin–Shtrikman bounds. We include the proof in appendix A; it is different from (and more elementary than) the previously known arguments.

## 2. Preliminaries

We start by recalling the definition of the effective conductivity of a periodic composite.

Let *Q* be the unit cube in ^{n}: *Q*=[0,1]^{n}. For a set *A*⊂*Q*, we defineand extend *a* periodically to ^{n}. This corresponds to a periodic two-phase composite where one phase is a good conductor (of heat or electricity) and covers a subset *A* of the cube, and the other phase is a poor conductor and covers the rest of the cube.

The effective conductivity *A*_{eff} is a tensor obtained by solving the cell problem. The value of *A*_{eff} is characterized by the following expression:(2.1)A standard source for the general theory of composites is Milton (2002).

If the set *A* is left invariant by any rigid motion that preserves the unit cube *Q* (i.e. *A* is cubically symmetric), then *A*_{eff} must be a scalar matrix *A*_{eff}=*a*_{eff}*I*.

For two-phase composites, it is also interesting to define the problem switching the phasesand(2.2)with *B*_{eff}=*b*_{eff}*I* for symmetric domains. Note that if *A*+*e*=*Q*\*A* (modulo the unit cube), for *e*=(1/2, 1/2, …, 1/2), then automatically *A*_{eff}=*B*_{eff}. This follows by a simple change of variables in (2.1).

We study the problem considered by Torquato *et al*. (2002, 2003), namely to find the sets *A* that are cubically symmetric and maximize *a*_{eff}+*b*_{eff}.

In order to avoid the complications of restricting our analysis to cubically symmetric sets *A*, we consider the quantity (tr *A*_{eff}+tr *B*_{eff}) instead of (*a*_{eff}+*b*_{eff}). In case *A* is cubically symmetric, both quantities coincide (by a factor *n*).

Note that the two-phase configuration separated by the Schwartz *P* surface is cubically symmetric and also satisfies *A*+*e*=*Q*\*A* for *e*=(1/2, 1/2, …, 1/2). The value of the sum (tr *A*_{eff}+tr B_{eff}) is invariant by all those symmetries.

## 3. Upper bound for tr *A*_{eff}+tr *B*_{eff}

For any dimension *n*, we will prove an upper bound for tr *A*_{eff}+tr *B*_{eff} independent of the set *A*. Then, we will check what conditions *A* must satisfy in order to achieve this bound. Our upper bound coincides with equation (1.1) for cubically symmetric composites, since in that case tr *A*_{eff}=*na*_{eff} and tr *B*_{eff}=*nb*_{eff}.

*Given A*_{eff} *and B*_{eff} *defined by* *(2.1)* *and* *(2.2)**, we have the upper bound*(3.1)

In the case when the conductivities are *σ*_{1} and *σ*_{2} instead of 1 and *ϵ*, the same proof gives the estimate

We start writing tr *A*_{eff}+tr *B*_{eff} in a long form

If we write *U*=(*u*_{1}, …, *u*_{n}) and *V*=(*v*_{1}, …, *v*_{n}), the above expression takes the cleaner form(3.2)

In order to find an upper bound, we restrict the minimum to a smaller set. In this case, it will be those *U* and *V* such that there exists a periodic potential function *p* so that(3.3)(3.4)

The key idea in this proof is that once we restrict our set of test functions *U* and *V* to this smaller set of gradients of a potential, the value of the minimum in (3.2) will *magically* become independent of the set *A*. We do the computations below.

We rewrite the expression (3.2) in terms of *p*,

The value of Δ*p* can be *any* function we choose as long as its average in the unit cube *Q* is zero. We can take the best choice of Δ*p*(*x*) pointwise, which is the minimum of the respective quadratic function inside the integral in the expression above. If *x* is in *A*, the minimum of would be achieved for . On the other hand, for *x* in *Q*\*A*, the best choice would be Δ*p*(*x*)=*θ*. Since we have |*A*|=(1/2)|*Q*|, this choice for Δ*p* has average zero in *Q*. Note that the value of *θ* is independent of *A* and we obtain exactly(3.5) ▪

Now, we want to check the conditions for the upper bound in proposition 3.1 to be achieved. We follow the proof and see that the upper bound is achieved if the optimal *U* and *V* for (3.2) happen to have the form given by (3.3) and (3.4).

For a given set *A*, the function *p* is easily computed by the equation

Equivalently, we can first solve the equation(3.6)then set *p*=−*θq*. We also have and .

In theorem 3.3, we provide necessary and sufficient conditions for a set *A* to be optimal.

*A set A*⊂*U with a smooth boundary* ∂*A realizes the upper bound of* *proposition 3.1* *if the corresponding function q defined in* *(3.6)* *has the following behaviour near each point x*∈∂*A*:(3.7)*for some matrix M*(*x*) *depending on each point x*∈∂*A such that M*(*x*).*ν*=0.

Note that *D*^{2}*q*(*x*) has a jump on ∂*A*; that is why we give two different values on each side of the surface.

To check if the computed estimate in proposition 3.1 is achieved for a given set *A*, we must check whether the computed *U* and *V* satisfy the Euler–Lagrange equation of (3.2). In other words,(3.8)(3.9)

Since Δ*q* is constant in *A* and *Q*\*A*, the above equalities are automatically satisfied within those sets. The only place to check is their boundary ∂*A*, where we must have the following *jump* in the normal derivative (at least if ∂*A* is smooth):(3.10)(3.11)where denotes the normal derivative to ∂*A* on the *A* side and denotes the normal derivative to ∂*A* on the *Q*\*A* side.

Note that from constructions *DU*=−*θD*^{2}*q* and *DV*=*θD*^{2}*q*, a jump in *DU* corresponds to a jump in *D*^{2}*q*. From (3.6), Δ*q* is a discontinuous function across ∂*A*. The Hessian matrix *D*^{2}*q* will be a continuous function on both sides and up to ∂*A*, but on ∂*A* it has two different values from each side. Note also that *q*∈*C*^{1,α} for every α<1 and *q*∈*C*^{1,1} if ∂*A* is smooth. Assuming that ∂*A* is smooth, we rewrite the equations (3.10) and (3.11) as(3.12)(3.13)

The vector field is well defined since *q*∈*C*^{1,α}. Taking derivatives tangentially on the surface ∂*A*, we see that *D*^{2}*q*.*τ* must have the same value on both sides *A* and *Q*\*A*. Since *D*^{2}*q* is self-adjoint, the tangential component of *D*^{2}*q*.*ν* must be equal on both sides of ∂*A*. But from the identities above, this implies that *D*^{2}*q*.*ν* has no tangential component if *ϵ*≠1. Thus, *ν* is an eigenvector of *D*^{2}*q* on ∂*A* and *D*^{2}*q* can only differ on each side of ∂*A* by a multiple of *ν*⊗*ν*. Moreover, recalling that *θ*=(1−*ϵ*)/(1+*ϵ*), the above equations tell us exactly what the eigenvalues are. Let us write *D*^{2}*qν*=*λ*^{+}*ν* in the *A* side, and *D*^{2}*qν*=*λ*^{−}*ν* in the *Q*\*A* side. We have(3.14)(3.15)which mean that *λ*^{+}=1 and *λ*^{−}=−1. This finishes the proof. ▪

The proof is slightly simpler if we assume *Q*\*A*=*A*+*e*, for *e*=(1/2, …, 1/2). In this case *A*_{eff}=*B*_{eff}, the upper bound on the sum becomes an upper bound for each trace, and (3.3) implies (3.4).

All that the proof requires is that *a*(*x*)+*b*(*x*) is a constant, so that we can apply the identity

We could consider periodic functions *a*(*x*) and *b*(*x*) that do not take only two values as long as *a*+*b* is a constant. The effective conductivities *A*_{eff} and *B*_{eff} are defined accordingly by (2.1) and (2.2). After redoing the computation in the proof of proposition 3.1, we would obtain some upper bound for tr *A*_{eff}+tr *B*_{eff} depending on the possible choices of *a* and *b*.

Note that the function *q* is completely defined by (3.6). The extra equation (3.7) can be understood as a free boundary condition; thus it imposes a restriction on the form of an optimal set *A*.

From theorem 3.3, we can draw an interesting conclusion. As expected, the optimality of a set *A* does not depend on the value of *ϵ*. This was pointed out by Torquato *et al*. (2003) based on their numerical computations.

It is also simple to show that a laminar composite is optimal. Indeed, if we choose *A*={(*x*, *y*, *z*): 1/4<*x*<3/4}, then the function *q* depends only on one dimension andwhich clearly satisfies the condition (3.7). Thus, the laminar structure achieves the optimal value for tr *A*_{eff}+tr *B*_{eff}, though it is not cubically symmetric and *A*_{eff} and *B*_{eff} are not scalar matrices.

Note that the laminar structure achieves the optimal value in any dimension. However, in dimension two, the upper bound can never be achieved by an isotropic composite due to the well-known phase exchange identity in two dimensions: *a*_{eff}.*b*_{eff}=*ϵ*.

## 4. Further analysis of the free boundary problem

The free boundary problem given by (3.6) and (3.7) has never been studied to our knowledge. At the present time, we are not able to show even the existence of other solutions besides the one-dimensional laminar one shown above. The main relevant question for the theory of composites would be whether or not a cubically symmetric solution exists. The computations by Torquato *et al*. (2002, 2003) suggest that there might be a solution whose free boundary coincides with or at least is very close to the Schwartz *P* surface. However, we have not been able to find any link between the free boundary problem and the mean curvature of the interface. In this section, we study the potential implications of the boundary ∂*A* being a minimal surface. We obtain indirect consequences that are easy to test numerically. We report some numerical results that strongly suggest that the structure with the Schwartz *P* minimal surface as the interface cannot achieve the upper bound (1.1).

Let us assume that we have a cubically symmetric optimal configuration corresponding to a set *A*. Equivalently, from theorem 3.3, we would have a function *q* solving the following overdetermined problems:(4.1)and(4.2)where *M*(*x*) is a symmetric matrix for every point *x*∈∂*A* such that *M*(*x*)*ν*=0 and *ν* is the normal vector to ∂*A* at *x*.

*Assume that ∂A is a smooth connected periodic surface with mean curvature zero for which the overdetermined problems* *(4.1)* *and* *(4.2)* *have a solution q*. *Then, q and q*_{ν} *are constants on ∂A, where q*_{ν} *denotes the normal derivative of q*.

As it is well known, the smoothness assumption is a consequence of the mean curvature zero up to dimension seven.

Given any smooth surface *S*, we can write the Laplacian of any *C*^{2} function *u* in the following way:where Δ_{S}*u* is the intrinsic Laplacian of *u* on *S*; *μ* is the sum of the principal curvatures of *S* (which is (*n*−1) times the mean curvature); and *u*_{νν} is the second derivative of *u* in the normal direction.

In case *S* has mean curvature zero, the formula simplifies to

The function *q* would be smooth on each side and up to ∂*A*. We can apply the above formula to *q* on each side of ∂*A*. For example, on the ‘*A*’ side we have Δ*q*=−1 from (4.1) and also *q*_{νν}=−1 from (4.2). Then, . Therefore, we have Δ_{S}*q*=0. If we do the same reasoning from the ‘*Q*\*A*’ side, we would also obtain Δ_{S}*q*=0 and no extra information.

Since ∂*A* is assumed to be connected, any harmonic function on it must be constant. Therefore, *q* is constant on ∂*A*. This proves the first part of lemma 4.1.

In order to prove the second part, we consider *q*_{ν}=*Dq*.*ν* and take a derivative in a direction *τ* tangential to ∂*A*. We haveFrom (4.2), *τ*^{t}*D*^{2}*qν*=0. Since ∂_{τ}*ν* is a tangential vector (this is true for any smooth surface) and *q* is constant on ∂*A*, *Dq*.∂_{τ}*ν*=0. Therefore, ∂_{τ}(*q*_{ν})=0 and *q*_{ν} must be constant on ∂*A*. This finishes the proof. ▪

Note that in the problems (4.1) and (4.2), we can add an arbitrary constant to the function *q*. Thus, we can make *q*=0 on ∂*A* in lemma 4.1. Using the maximum principle on equation (4.1), we see that the function *q* would be negative in the interior of *A* and positive outside of *A*. Therefore, if a set *A* was optimal and ∂*A* had mean curvature zero, the function *q* would also be a solution of the following free boundary problem:(4.3)which is completed by the compatibility condition that *q* must be *C*^{1} across the free boundary {*q*=0}. A similar problem has attracted some attention recently that is obtained by exchanging −1 and 1 in equation (4.3) (see Henrik Shahgholian (in press) and references therein). However, our opposite choice of signs greatly changes the nature of the problem. A more similar equation was considered by Andersson & Weiss (2006). The function *q* is a critical point of the nonconvex functional(4.4)

On the other hand, if we consider the positive part of *q* only, we arrive at a better known free boundary problem(4.5)

This is usually referred to as the *one-phase* problem and has been studied extensively within the free boundary community (e.g. Aguilera *et al*. 1986; Caffarelli 1998).

A strange consequence of a zero mean curvature free boundary would be that the three equations (4.1) and (4.2), (4.3) and (4.5) would have the same solution *q*.

Recalling that *U*=*θDq*, we observe that one implication of lemma 4.1 is that |*U*| would have to be constant on ∂*A*, if ∂*A* had zero mean curvature. We ran a finite element code to solve the cell problem in [−1,1]^{3}. We tried first with a uniform mesh and then with more accurate adaptive mesh. We used different nodal approximations to the *P* surface (Schwarz & Gompper 1999; Gandy *et al*. 2001). We tested the computed effective conductivities and indeed it was very close to optimal. But the computed values of |*U*| showed a significant difference at the points (1/2, 1/2, 1/2) and (1, 1/2, 0) which lie on the interface (approx. 0.32 and 0.39). These values did not vary much when we changed from the uniform mesh to the adaptive mesh, or when we added more degrees of freedom. This suggests that the optimal structure cannot be separated by the *P* surface.

The main idea for our test is that the effective conductivity may not be very sensitive to small changes in the domain, and thus hard to test numerically. However, the value of the normal derivative of *q* on ∂*A* is more sensitive to changes in the shape of *A* and easier to compute. In this way, what we test is an indirect consequence of a zero mean curvature boundary. The results that we obtained show strong evidence that if ∂*A* is the Schwartz *P* surface, then the upper bound (1.1) is not attained.

Another test that we performed was to compute a solution to equation (4.3) numerically. We set up an elementary finite difference scheme where the value of *q* at each grid point is given by the average of the six neighbours plus *h*^{2}/6 if *q* is positive, or minus *h*^{2}/6 if *q* is negative at that point (*h* is the distance between two neighbouring points on the grid). Then we solved that equation on the grid by an iterative method using a suitable relaxation parameter. We tested the problem in the cube [−1,1]^{3} with various grid sizes (from 20 to 250) and in every case we obtained that the normal derivative of *q* at the free boundary point (1/2, 1/2, 1/2) is approximately 0.47 and the normal derivative at the free boundary point (1, 1/2, 0) is 0.37. Even though we have not done any error analysis for this finite element scheme, the values we obtain when varying the grid size are all very similar. This seems to be a significant difference that also suggests that the minimal surface structure may not be optimal.

At the present time, the source code of our tests can be found on the website of the author (www.cims.nyu.edu/∼silvestr/source).

The second test that we described would not suffice to give a definitive answer to the question. The functional (4.4) is not convex, thus we cannot prove uniqueness for equation (4.3). We cannot rule out the possibility that our numerical approximation is picking a different solution that is also cubically symmetric; however, this seems to us unlikely. For that reason, the first test provides stronger evidence.

If the Schwartz *P* surface is not the optimal configuration for the sum of two effective conductivities, then the question is why does it approximate the upper bound with such a high level of precision. According to Torquato *et al*. (2003), if the conductivities of the two phases are 0.1 and 1, then the effective conductivities range between 0.38 and 0.427 for all possible structures. This is a small range. A surface that has a similar shape to the optimal (assuming there is actually a surface that achieves the upper bound) may have a very similar effective conductivity.

The high accuracy with which the Schwartz *P* surface configuration approximates the upper bound suggests that the free boundary problems (4.1) and (4.2) may have a cubically symmetric solution with a free boundary near the Schwartz *P* surface. Our computations suggest, however, that the free boundary does not coincide exactly with this surface. In any case, equations (4.1) and (4.2) define a new type of free boundary problem that would be worth studying and understanding better.

## Acknowledgments

I would like to thank Luca Heltai for help with the implementation of an adaptive finite element code for one of the numerical tests. I would also like to thank very specially Robert Kohn for telling me about the problem and for several useful discussions related to this paper.

## Footnotes

- Received January 13, 2007.
- Accepted May 13, 2007.

- © 2007 The Royal Society