## Abstract

Composites with subwavelength features exhibit effective properties that depend on microstructure morphology and materials, which can be adjusted to obtain enhanced characteristics. We detail the systematic design of electromagnetic metamaterials composed of dielectric inclusions in a ferroelectric matrix that, under an applied voltage, present an optimized effective tunability higher than the bulk due to a nonlinear local electric field enhancement. The effect of volume fraction, losses, and biasing field on homogenized properties is investigated and the analysis of the photonic band diagram is carried out, providing the frequency dependence of the anisotropic effective index and tunability. Such metaceramics can be used in microwave antennas and components with higher reconfigurability and reduced power consumption.

Published by The Optical Society under the terms of the Creative Commons Attribution 4.0 License. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.

## 1. Introduction

Multifunctional and reconfigurable systems are increasingly demanded for a wide range of applications in Electromagnetics and Photonics [1–3]. In parallel, the development of artificial media and metamaterials [4,5], allowing unprecedented control of electromagnetic (EM) waves by carefully engineered subwavelength structures, has stimulated the research interest in tunable materials. Indeed, the ability to introduce tunable components has been established as a straightforward route to EM and metamaterials reconfigurability, albeit the added fabrication complexity in comparison with static devices. As a matter of fact, it is envisioned that reconfigurable metamaterial structures have the potential to mitigate some of the most significant limitations of static metamaterials such as narrow bandwidth operations, high losses, and tolerance sensitivity. Often, in addition to the necessary control logic, adding reconfigurability to a metamaterial structure requires the introduction of extra components and materials to a standard design, hence potentially affecting the system’s performance, power consumption and weight. Consequently, there are still significant theoretical and technical obstacles to the production and deployment of reconfigurable metamaterials within commercial products, particularly regarding versatility, operational bandwidth, performance, complexity of realization, robustness, speed of tuning, and costs. Those engineering challenges are driving efforts of academia and industry towards new technological solutions and experimental testing of reconfigurable metamaterials, which is currently a very active area of research. Experimental demonstrations across frequency ranges have been reported and include microwave devices [6,7], lenses [8,9] tunable band diagrams [10,11], reconfigurable orbital angular momentum [12], topological insulators [13], infrared emitters [14] negative effective index media [15–17], hyperbolic metamaterials [18] and optically controlled metasurfaces [19] using different stimuli and tuning mechanisms such as mechanically adjustable elements [20–24], microfluidics [25,26], graphene [27] or liquid crystals [28,29].

In microwave engineering, ferroelectric materials play a crucial role in applications needing reconfigurability, such as antenna beam steering, phase shifters, filters, and tunable power splitters, [30]. Both thin films and bulk ceramics are used for frequency agile components [31,32] and metamaterials [15,33]. Ferroelectric materials possess a strong nonlinear dependence of their permittivity $\varepsilon$ on an applied electric field $E$, with key requirements for antenna and microwave applications being large tunability and low losses. The permittivity values are usually high even at microwave frequencies, and therefore it has been considered to mix ferroelectric materials with low-index and low-loss non-tunable dielectrics or to synthesize porous ceramics in order to reduce both permittivity values and losses. The effective parameters of those composites have been investigated [34–37] using analytical effective medium approaches for low filling fraction of dielectric and have been successfully compared with numerical simulations and experiments. It has been shown that the permittivity can be reduced while losses are much less sensitive to the dielectric phase addition, and in some situations lead to a small increase of the tunability of the mixtures. Recently, refined theoretical and numerical models have been proposed [38–40], and the predicted behaviour compared well to experimental data [41]. Those methods account for the nonlinear coupling of electrostatic biasing and field dependent permittivity. It was shown that the permittivity can be reduced while maintaining high tunability and low losses, due to local field enhancement in the ferroelectric phase.

We propose here to design the microstructure of metamaterials composed of a ferroelectric and low index non tunable dielectric (see Fig. 1) through topology optimization [42]. Our aim is to find structures with desired effective properties: high tunability and minimal losses. Using a nonlinear coupled model that takes into consideration the subwavelength field enhancement and ferroelectric response, we systematically obtain the material distribution in the unit cell and study the effect of volume fraction and biasing field strength. Furthermore, we conduct a multi-objective optimization for maximal tunability and minimal losses and present the Pareto front, revealing a trade-off between the two characteristics. We finally study the photonic bands of a particular example of optimized geometry, showing an adjustable bandgap, derive the frequency dependence of the effective permittivity and compare the results with the homogenization approach.

## 2. Methods

The study is restricted to two-dimensional metamaterials with square unit cell $\Omega$ of size $d<<\lambda$ and periodic boundary conditions, where $\lambda$ is the wavelength under consideration. Since we target the microwave frequency range ($f=3.8\,$GHz, $\lambda \simeq 8\,$cm), the maximal value for the periodicity is roughly $\lambda /10=8\,$mm. On the other hand, we consider a bulk ferroelectric so that nano-scale effects such as domain walls motion [43], and micro-scale grains ($\ell \simeq 10\,{\mu {\rm m}}$, [44]) are not resolved in the present model but are incorporated in the dielectric permittivity. This places a minimal length scale for the validity of the study at around $10 \ell \simeq 100\,{\mu {\rm m}}$. The materials we employ are barium strontium titanate (BST) and a low index ($\varepsilon _\textrm {d}=3$) dielectric such as common thermoplastic polymer acrylonitrile butadiene styrene (ABS).

#### 2.1 Coupling nonlinear ferroelectric permittivity and electrostatics

The tunability of barium strontium titanate $\textrm {Ba}_{x}\textrm {Sr}_{1-x}\textrm {TiO}_3$ with a barium ratio of $x = 0.6$ samples fabricated using a conventional sintering method was measured in the static regime (DC) and at microwave frequencies (MW) [44] and fitted to a Vendik model [45] to represent the electric field dependence of the permittivity:

To calculate the total electric field in the metamaterial, one has to solve for the potential $V$ satisfying Gauss’s law for a given permittivity distribution $\varepsilon$:

biased by a constant uniform electric field $\boldsymbol E_\textrm {B} ={E_\textrm {B}} \boldsymbol {x}$. The electric field $\boldsymbol E=-\mathbf{\nabla} V$ derived from the solution of Eq. (2) depends on the permittivity distribution, which itself depends on the electric field. The coupled system formed of Eqs. (1) and (2) is solved iteratively using a Picard method until convergence on the norm of the electric field $||\boldsymbol E_{i+1} - \boldsymbol E_i||<10^{-5}$. Although initially constant, the permittivity in the ferroelectric phase is spatially varying due to the non-uniform distribution of the total electric field.#### 2.2 Homogenized properties

The aim of homogenization approaches is to retrieve a simplified, but almost equally precise, description of the material response by averaging out material properties at the subwavelength level. Research in this area has been deeply linked to and influenced by the development of metamaterials and their application in physics, particularly for electromagnetic phenomena. Analytical models for the effective permittivity previously employed in the literature such as Maxwell-Garnett or Bruggeman theories are limited to a few canonical shapes of the inclusions, and cannot handle arbitrary geometries and media with spatially varying properties, which has to be accounted for in our model because of the field induced local permittivity change. In this study, the effective permittivity of the metamaterials is calculated using a two scale convergence homogenization technique [46,47]. We will focus on 2D geometries and TM polarization (since the TE polarization effective permittivity is trivially the arithmetic mean of the basic components permittivity), resulting in a scalar wave equation for the $z$ component of the magnetic field. The main idea is to select two scales in the study: a microscopic one (the size of the unit cell) and a mesoscopic one (the size of the bulk), controlled by a real parameter $\nu >0$. Denoting $\boldsymbol {r} = (x, y)^\textrm {T}$ the position vector, the approach consists in introducing the ansatz ${{\mathcal {H}}}(\boldsymbol {r})={{\mathcal {H}}}_0(\boldsymbol {r},\boldsymbol {r}/\nu ) + \nu {{\mathcal {H}}}_1(\boldsymbol {r},\boldsymbol {r}/\nu ) + \nu ^{2} {{\mathcal {H}}}_2(\boldsymbol {r},\boldsymbol {r}/\nu ) + \cdots$ for the magnetic field ${{\mathcal {H}}}$ solution of time harmonic Maxwell’s equation and performing an asymptotic analysis as $\nu \rightarrow 0$. One then needs to find the solutions $\psi _j$ of two annex problems $\mathcal {P}_j$, $j={x, y}$:

The homogenized permittivity $\tilde {{\varepsilon }}$ is then obtained as:

where $\left \langle {.} \right \rangle$ denotes the mean value over the unit cell, $\boldsymbol I$ is the $2\times 2$ identity matrix and $\boldsymbol {\phi }_{ij} = \langle \varepsilon \mathbf{\nabla} \psi _i \rangle _j$ are correction terms. This analysis allows us to obtain the effective parameters at higher frequencies, in contrast with capacitance-based models valid in the electrostatic regime. Contrarily to most homogenization procedures that are based on a quasi-static approximation, in the two scale convergence method, the frequency is fixed and the characteristic size of the system (the periodicity of the composites) tends to zero, making the study of the frequency dependence of the effective parameters possible.#### 2.3 Topology optimization

Inverse design, where a specific target is searched according to engineering constraints, is an actively researched topic in the nanophotonics and metamaterial community [48–50]. Topology optimization has led to the design of a large range of devices such as invisibility cloaks [51,52], illusion devices [53], metasurfaces [54,55] and metamaterials with tailored properties [56–59], showing non intuitive material distributions, allowing to manipulate light matter interaction with unprecedented capability. In the topology optimization procedure, the permittivity distribution is parametrized by a continuous scalar density function $\rho \in [0,1]$. To avoid small features and pathological “chessboard” patterns, we use a filtered density obtained by solving the following Helmholtz equation [60] with periodic boundary conditions on $\Omega$:

with $R$ a real positive parameter characterizing the filter radius ($R=0.05d$ in the rest of the paper). In order to enforce a binary design, we apply a threshold projection [61,62]:## 3. Results

Numerical results are obtained using `Python` and open source libraries: Eqs. (1), (2), (3) and (5) are solved by the finite element method with the `fenics` library [64], the sensitivity of the objective $g$ with respect to the design variables $p$ are computed with the adjoint method [65] using automatic differentiation through `dolfin-adjoint` [66] and the method of moving asymptotes (MMA, [67]) is employed for the gradient based optimization with the `nlopt` package [68].

#### 3.1 Maximizing tunability

The first objective is to maximize the tunability of the composites along the direction of the applied electric field defined as $\tilde {\eta }(E)=\tilde {\varepsilon }_{xx}(0)/\tilde {\varepsilon }_{xx}(E)$:

*i.e.*the proportion of BST in the metamaterial). Indeed, there is an infinity of possible material distribution with the same volume fraction, but they would lead to different effective material properties. Our aim is to find the distribution that maximize our optimization objective. Those constraints are imposed as an inequality because the MMA algorithm we use cannot enforce equality constraints. To target approximately a given volume fraction $f_\textrm {cons}$, we set $f_\textrm {min} = f_\textrm {cons}-\delta _f$ and $f_\textrm {max} = f_\textrm {cons}+\delta _f$ and $\delta _f$ is decreased linearly from 0.1 to 0.01 during the course of the optimization. The optimization is initialized with a density $\rho _0 = 1-\textrm {exp}(-(x^{2}+y^{2})/r_f^{2})$, with $r_f = \sqrt {(1 - f) / \pi \log {2})}$. Additionally, our numerical experiments have shown that field enhancement is suppressed when the inclusion are connected along the $y$-axis, we thus constrain further the material to be ferroelectric (

*i.e.*$\rho =1$) at the top and bottom of the unit cell in a strip of height $0.05d$.

Results are displayed on Fig. (2(a)), where we plot the effective tunability enhancement of structures optimized for different volume constraints and a biasing field of 1 kV/mm. Tunability higher than bulk BST is achieved in each case except for the smallest volume fraction where it is 4% weaker. A maximum enhancement of 34% is achieved for a structure containing 83% of ferroelectric, which indicates a trade-off between composites containing more tunable material and the volume of dielectric phase which induces a field concentration within the unit cell. Note that the obtained effective properties are anisotropic due to to the asymmetry of material distribution, but we are only focusing on the properties parallel to the applied electric field. As $f$ increases, both components of the effective permittivity tensor increase (see Table (1)), with strong anisotropy for lower values of filling fraction, which correspond to structures connected along the $x$-axis. As the target volume fraction increases, the topology evolves from layered structure along $y$ to cross like shapes, a “bow tie” looking dielectric structure for the highest enhancement ($f=83$%) to a needle like inclusion with long side along $y$ for $f=83$%. These non-trivial topologies show the generality and versatility of the approach which does not require any geometrical parametrization, generating freeform microstructures. The field maps of the electric field and corresponding permittivity in subfigures (2(b-e)) show that the field enhancement is maximal between the dielectric parts of the microstructure, which leads to a localized tuning in the ceramic phase.

The previous optimization was performed for a fixed biasing field of 1 kV/mm. In order to assess the response of the metamaterials, we computed the effective tunability as a function of $E$ (see Fig. (3)). The maximal values of tunability enhancement are obtained as expected around 1 kV/mm. We note that the permittivity model has an inflexion point around $E=1.18\,$kV/mm so that a larger local field variation leads to a stronger change in permittivity. The improved tunability is maintained, albeit smaller, for different biasing intensities except at higher fields. Further optimization studies not reported here with lower and higher field strength led to similar material topologies, indicating the robustness of the design process.

#### 3.2 Bi-objective optimization: high tunability and low loss metamaterials

Our next aim is to find metamaterial structures that will at the same time enhance the effective tunability while reducing the effective losses. We now use a complex permittivity value $\varepsilon _\textrm {f}^\textrm {C}$ with the same dependence on the electric field (1) but with constant value of the loss tangent $\tan \delta _\textrm {f}=10^{-3}$, that is $\varepsilon _\textrm {f}^\textrm {C}(E)=\varepsilon _\textrm {f}(E)\,(1-\tan \delta _\textrm {f})$. We define the loss reduction factor as:

where the homogenized loss tangent $xx$ component is $\tan \tilde {\delta }_{xx}(E) = -\mathrm {Im}\, \tilde {\varepsilon }_{xx}(E) / \mathrm {Re}\,\tilde {\varepsilon }_{xx}(E)$. We then have to solve the following maximization problem:The Pareto front of this bi-objective problem is given in Fig. (4), where each point represent a different value for $p$. One can notice that there is a trade-off between those two concurrent goals and observe two regimes: the first with enhanced tunability (18% to 25%) and moderate loss reduction (18% to 24%), and the second where the tunability is weaker than bulk BST (-11% to -7%) but with a stronger decrease in losses (33% to 43%). Depending on the application requirements, one may choose between different topologies to mitigate the advantages and drawbacks of the two material properties. Note that the Pareto front is effectively discontinuous, and is jumping in between the two regimes of high tunability/moderate loss and low tunability/low losses. This is unrelated to the parameter $p$ which is just a numerical coefficient used in the optimization procedure.

We now look at the response of the optimized composites as the biasing field is varied. A metric for the performance of ceramics that takes into account tunability and losses is the commutation quality factor (CQF) defined as:

For the structures associated with the second region of the Pareto set ($p = 0.10$ and $p = 0.25$), the tunability decreases significantly as the field increases (see Fig. (5(a))), while losses are strongly reduced, and this reduction actually increases with biasing strength (see Fig. (5(b))). We show in Fig. (5(c)) the normalized CQF defined as $\kappa = \tilde {K}/K_\textrm {f}$, where $\tilde {K}$ and $K_\textrm {f}$ correspond to the CQF for the homogenized metamaterial and bulk BST respectively. For $p=0.10$, the QCF is reduced while it is of the same order than than the bulk for the case $p=0.25$, for all electric fields. This highlights the limits of such a metric where it is hard to distinguish independently the contributions of losses and tunability: indeed one can have, as the later case, a material with reduced tunability and losses with similar QCF. For the metaferroelectrics belonging to the other optimized region ($p = 0.35$ and $p = 0.95$), a maximum of tunability enhancement is observed as expected around 1KV/mm and decreasing for higher fields, similarly to the lossless case. Losses are reduced significantly and the stronger the bias, the smaller the loss tangent. Finally, in both cases the associated QCF is much larger than the bulk (around 4 times at 1KV/mm) and decreases as $E$ increases. The enhanced performances reported here are clearly a trade-off in between those two conflicting objectives (losses and tunability), but other considerations may come into play in the choice of and adequate metamaterial such as the volume fraction or the effective anisotropy.

## 4. Photonic bands properties

We now investigate the spectral properties of one example of optimized composite (Section (3.1) for $f=0.47$). The eigenvalues $k_n=\omega _n/c$, of Maxwell’s operator for TM polarization, *i.e* the solutions of:

`GetDP`[69] with its interface to the

`SLEPc`library [70]. The effective permittivity can be obtained from the band diagram according to:

Figure (6(a)) shows the photonic band diagram associated with the optimized lossless composite for $f=0.47$ obtained in Section (3.1) (see inset in Fig. (2) for the unit cell topology) and observe a tunable band gap between the first and second bands. When unbiased, the photonic band gap is centred at $\omega =0.0829\omega _0$ with a width of $w=0.0173\omega _0$, and when the electric field is applied, its centre is blueshifted to $0.0848\omega _0$ (~2%) and slightly broadened ($w=0.0175\omega _0$, ~0.6%). Another much narrower bandgap between the second and third bands is present as well, with centre $\omega =0.114\omega _0$ (resp. $0.118\omega _0$) and width $0.0017\omega _0$ (resp. $0.0019\omega _0$) for the unbiased and biased case respectively. Note that while in the first case the two band gap overlap when the electric field is on or off, in the later case the two forbidden propagation bands are disjoint since the bandwidth is smaller.

The isofrequency contours for the first band are displayed on (6(b)). When no field is applied (full lines), we have elliptical contours with long axis along $k_y$, meaning a diagonal anisotropic effective index with higher values along $x$, which is consistent with the results of the previous homogenization study. This is valid for $\omega \rightarrow 0$, and when the frequency increases so does the effective index. Applying a voltage on the photonic crystal deforms the isofrequency contours, increasing the semi-axis of the ellipse along $k_y$ but almost keeping the same along $k_x$: in other words, the tunability is high in the $x$ direction but practically negligible in the orthogonal direction, in agreement with the results from Section (3.1).

To corroborate those qualitative results, we extracted the effective permittivity tensor from the first band and plotted the results in Fig. (7(a)). Firstly, the values agree very well with those obtained with two scale homogenization (see the markers at $\omega =0$), and we observe an increase of the components of the effective permittivity. The permittivity along $x$ is enhanced further as we get closer to the edge of the gap, and hence, due to the shift of the bandgap from applying a voltage, the tunability along $x$ is increased, reaching a value of around 4 at $\omega =0.041\omega _0$ (see Fig. (7(b))). This effect, even if it is narrowband, provides a way to even further enhance the tunability in this kind of metastructures. Note that the design was optimized to enhance the homogenized tunability, not to maximize bandgaps / tunability frequency dependence, which is out of the scope of this paper.

## 5. Conclusion

We presented an inverse design methodology to obtain microstructures of ferroelectric/dielectric composites with desired and enhanced properties. This inverse homogenization approach allows us to find the material distribution within the periodic unit cell that gives a maximized effective tunability for a given volume fraction of ceramic. Tunability enhancement as high as 34% compared to bulk BST for a filling fraction of 87% are obtained. Taking dissipation into account, we have shown that one has to mitigate two competing objectives, namely high tunability and low losses. In addition, we have studied the spectral properties of one optimized photonic crystal, showing a tunable bandgaps, frequency dependent enhanced tunability and confirming the homogenization approach. Extension to 3D metamaterial topologies, fabrication and experimental validation of the designed composites will be the subject of future studies. The method proposed here provide guidelines for the systematic design of metamaterials with desired properties, allowing further control of electromagnetic waves in antennas and radio frequency applications. The proposed designs can provide a range of permittivity values suitable for instance for the implementation of graded index lenses employed in beam steering antennas. At the same time, the optimized tunability would increase the reconfigurability of such devices, avoiding the need for complex feed network. On the other hand, a co-design of and integrated biasing network of electrodes will be needed. The general approach followed here can be applied to other type of tunable systems in different frequency ranges where local field enhancement and amplification is relevant, including for instance liquid crystals based devices, ferromagnetic metamaterials and field-enhanced carrier dynamics in doped semiconductors. Finally, since ferroelectric material behaviour is inherently governed by multiple physical processes, it would be of great interest to optimize in parallel the effective electro-mechanical and piezoelectric properties of metaferroelectrics, in order for example to mitigate the brittleness of ceramic materials and obtain composites with high tunability, low losses and strong mechanical compliance.

## Funding

Engineering and Physical Sciences Research Council (EP/N010493/1, EP/P005578/1, EP/R035393/1).

## Acknowledgments

The authors would like to thanks Henry Giddens and Hangfeng Zhang for performing the measurements of ferroelectric permittivity used in this paper.

## Disclosures

The authors declare no conflicts of interest.

## Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

## References

**1. **G. Oliveri, D. H. Werner, and A. Massa, “Reconfigurable electromagnetics through metamaterials—a review,” Proc. IEEE **103**(7), 1034–1056 (2015). [CrossRef]

**2. **J. P. Turpin, J. A. Bossard, K. L. Morgan, D. H. Werner, and P. L. Werner, “Reconfigurable and tunable metamaterials: a review of the theory and applications,” Int. J. Antennas Propag. **2014**, 1–18 (2014). [CrossRef]

**3. **J. Y. Ou, E. Plum, L. Jiang, and N. I. Zheludev, “Reconfigurable photonic metamaterials,” Nano Lett. **11**(5), 2142–2144 (2011). [CrossRef]

**4. **H. Chen, C. Chan, and P. Sheng, “Transformation optics and metamaterials,” Nat. Mater. **9**(5), 387–396 (2010). [CrossRef]

**5. **N. Zheludev, “The road ahead for metamaterials,” Science **328**(5978), 582–583 (2010). [CrossRef]

**6. **S. Gevorgian, E. Carlsson, E. Wikborg, and E. Kollberg, “Tunable microwave devices based on bulk and thin film ferroelectrics,” Integr. Ferroelectr. **22**(1-4), 245–257 (1998). [CrossRef]

**7. **R. Jakoby, P. Scheele, S. Muller, and C. Weil, “Nonlinear dielectrics for tunable microwave components,” in 15th International Conference on Microwaves, Radar and Wireless Communications (IEEE Cat. No.04EX824), vol. 2 (2004), pp. 369–378 Vol.2.

**8. **K. Chen, Y. Feng, F. Monticone, J. Zhao, B. Zhu, T. Jiang, L. Zhang, Y. Kim, X. Ding, S. Zhang, A. Alù, and C.-W. Qiu, “A reconfigurable active Huygens’ metalens,” Adv. Mater. **29**(17), 1606422 (2017). [CrossRef]

**9. **M.-C. Tseng, F. Fan, C.-Y. Lee, A. Murauski, V. Chigrinov, and H.-S. Kwok, “Tunable lens by spatially varying liquid crystal pretilt angles,” J. Appl. Phys. **109**(8), 083109 (2011). [CrossRef]

**10. **J. Lo, J. Sokoloff, T. Callegari, and J. P. Boeuf, “Reconfigurable electromagnetic band gap device using plasma as a localized tunable defect,” Appl. Phys. Lett. **96**(25), 251501 (2010). [CrossRef]

**11. **R. Getz, D. M. Kochmann, and G. Shmuel, “Voltage-controlled complete stopbands in two-dimensional soft dielectrics,” Int. J. Solids Struct. **113-114**, 24–36 (2017). [CrossRef]

**12. **Q. Liu, Y. Zhao, M. Ding, W. Yao, X. Fan, and D. Shen, “Wavelength- and OAM-tunable vortex laser with a reflective volume Bragg grating,” Opt. Express **25**(19), 23312 (2017). [CrossRef]

**13. **X. Cheng, C. Jouvaud, X. Ni, S. H. Mousavi, A. Z. Genack, and A. B. Khanikaev, “Robust reconfigurable electromagnetic pathways within a photonic topological insulator,” Nat. Mater. **15**(5), 542–548 (2016). [CrossRef]

**14. **X. Liu and W. J. Padilla, “Reconfigurable room temperature metamaterial infrared emitter,” Optica **4**(4), 430–433 (2017). [CrossRef]

**15. **H. Zhao, L. Kang, J. Zhou, Q. Zhao, L. Li, L. Peng, and Y. Bai, “Experimental demonstration of tunable negative phase velocity and negative refraction in a ferromagnetic/ferroelectric composite metamaterial,” Appl. Phys. Lett. **93**(20), 201106 (2008). [CrossRef]

**16. **H. Němec, P. Kužel, F. Kadlec, C. Kadlec, R. Yahiaoui, and P. Mounaix, “Tunable terahertz metamaterials with negative permeability,” Phys. Rev. B **79**(24), 241108 (2009). [CrossRef]

**17. **F. Castles, J. A. J. Fells, D. Isakov, S. M. Morris, A. A. R. Watt, and P. S. Grant, “Active metamaterials with negative static electric susceptibility,” Adv. Mater. **32**(9), 1904863 (2020). [CrossRef]

**18. **J. A. Roberts, S.-J. Yu, P.-H. Ho, S. Schoeche, A. L. Falk, and J. A. Fan, “Tunable hyperbolic metamaterials based on self-assembled carbon nanotubes,” Nano Lett. **19**(5), 3131–3137 (2019). [CrossRef]

**19. **Q. Wang, E. T. F. Rogers, B. Gholipour, C.-M. Wang, G. Yuan, J. Teng, and N. I. Zheludev, “Optically reconfigurable metasurfaces and photonic devices based on phase change materials,” Nat. Photonics **10**(1), 60–65 (2016). [CrossRef]

**20. **E. Hajiesmaili and D. R. Clarke, “Reconfigurable shape-morphing dielectric elastomers using spatially varying electric fields,” Nat. Commun. **10**(1), 183–187 (2019). [CrossRef]

**21. **Z. Han, K. Kohno, H. Fujita, K. Hirakawa, and H. Toshiyoshi, “MEMS reconfigurable metamaterial for terahertz switchable filter and modulator,” Opt. Express **22**(18), 21326–21339 (2014). [CrossRef]

**22. **K. Bai, X. Cheng, Z. Xue, H. Song, L. Sang, F. Zhang, F. Liu, X. Luo, W. Huang, Y. Huang, and Y. Zhang, “Geometrically reconfigurable 3D mesostructures and electromagnetic devices through a rational bottom-up design strategy,” Sci. Adv. **6**(30), eabb7417 (2020). [CrossRef]

**23. **Z. Wang, L. Jing, K. Yao, Y. Yang, B. Zheng, C. M. Soukoulis, H. Chen, and Y. Liu, “Origami-based reconfigurable metamaterials for tunable chirality,” Adv. Mater. **29**(27), 1700412 (2017). [CrossRef]

**24. **J. Zendejas, J. Gianvittorio, Y. Rahmat-Samii, and J. Judy, “Magnetic MEMS reconfigurable frequency-selective surfaces,” J. Microelectromech. Syst. **15**(3), 613–623 (2006). [CrossRef]

**25. **T. S. Kasirga, Y. N. Ertas, and M. Bayindir, “Microfluidics for reconfigurable electromagnetic metamaterials,” Appl. Phys. Lett. **95**(21), 214102 (2009). [CrossRef]

**26. **M. Wang, C. Trlica, M. R. Khan, M. D. Dickey, and J. J. Adams, “A reconfigurable liquid metal antenna driven by electrochemically controlled capillarity,” J. Appl. Phys. **117**(19), 194901 (2015). [CrossRef]

**27. **W. Y. Kim, H.-D. Kim, T.-T. Kim, H.-S. Park, K. Lee, H. J. Choi, S. H. Lee, J. Son, N. Park, and B. Min, “Graphene–ferroelectric metadevices for nonvolatile memory and reconfigurable logic-gate operations,” Nat. Commun. **7**(1), 10429 (2016). [CrossRef]

**28. **D. H. Werner, D.-H. Kwon, I.-C. Khoo, A. V. Kildishev, and V. M. Shalaev, “Liquid crystal clad near-infrared metamaterials with tunable negative-zero-positive refractive indices,” Opt. Express **15**(6), 3342 (2007). [CrossRef]

**29. **S. Slussarenko, A. Murauski, T. Du, V. Chigrinov, L. Marrucci, and E. Santamato, “Tunable liquid crystal q-plates with arbitrary topological charge,” Opt. Express **19**(5), 4085 (2011). [CrossRef]

**30. **A. Tagantsev, V. Sherman, K. Astafiev, J. Venkatesh, and N. Setter, “Ferroelectric materials for microwave tunable applications,” J. Electroceram. **11**(1/2), 5–66 (2003). [CrossRef]

**31. **O. Vendik, E. Hollmann, A. Kozyrev, and A. Prudan, “Ferroelectric tuning of planar and bulk microwave devices,” J. Supercond. **12**(2), 325–338 (1999). [CrossRef]

**32. **M. J. Lancaster, J. Powell, and A. Porch, “Thin-film ferroelectric microwave devices,” Supercond. Sci. Technol. **11**(11), 1323–1334 (1998). [CrossRef]

**33. **T. H. Hand and S. A. Cummer, “Frequency tunable electromagnetic metamaterial using ferroelectric loaded split rings,” J. Appl. Phys. **103**(6), 066105 (2008). [CrossRef]

**34. **K. Astafiev, V. Sherman, A. Tagantsev, and N. Setter, “Can the addition of a dielectric improve the figure of merit of a tunable material?” J. Eur. Ceram. Soc. **23**(14), 2381–2386 (2003). [CrossRef]

**35. **L. Jylha and A. H. Sihvola, “Tunability of granular ferroelectric dielectric composites,” Prog. Electromagn. Res. **78**, 189–207 (2008). [CrossRef]

**36. **V. O. Sherman, A. K. Tagantsev, N. Setter, D. Iddles, and T. Price, “Ferroelectric-dielectric tunable composites,” J. Appl. Phys. **99**(7), 074104 (2006). [CrossRef]

**37. **V. Sherman, A. Tagantsev, and N. Setter, “Tunability and loss of the ferroelectric-dielectric composites,” in 14th IEEE International Symposium on Applications of Ferroelectrics, 2004. ISAF-04. 2004, (IEEE, 2004), pp. 33–38.

**38. **L. Padurariu, L. Curecheriu, C. Galassi, and L. Mitoseriu, “Tailoring non-linear dielectric properties by local field engineering in anisotropic porous ferroelectric structures,” Appl. Phys. Lett. **100**(25), 252905 (2012). [CrossRef]

**39. **A. Cazacu, L. Curecheriu, A. Neagu, L. Padurariu, A. Cernescu, I. Lisiecki, and L. Mitoseriu, “Tunable gold-chitosan nanocomposites by local field engineering,” Appl. Phys. Lett. **102**(22), 222903 (2013). [CrossRef]

**40. **B. Vial and Y. Hao, “Enhanced tunability in ferroelectric composites through local field enhancement and the effect of disorder,” J. Appl. Phys. **126**(4), 044102 (2019). [CrossRef]

**41. **L. Padurariu, L. Curecheriu, V. Buscaglia, and L. Mitoseriu, “Field-dependent permittivity in nanostructured BaTiO_{3} ceramics: modeling and experimental verification,” Phys. Rev. B **85**(22), 224111 (2012). [CrossRef]

**42. **M. P. Bendsoe and O. Sigmund, * Topology Optimization: Theory, Methods, and Applications* (Springer Science & Business Media, 2013).

**43. **J. Guyonnet, “Domain walls in ferroelectric materials,” in Ferroelectric Domain Walls (Springer International Publishing, 2014), pp. 7–24.

**44. **H. Zhang, H. Giddens, Y. Yue, X. Xu, V. Araullo-Peters, V. Koval, M. Palma, I. Abrahams, H. Yan, and Y. Hao, “Polar nano-clusters in nominally paraelectric ceramics demonstrating high microwave tunability for wireless communication,” J. Eur. Ceram. Soc. **40**(12), 3996–4003 (2020). [CrossRef]

**45. **O. G. Vendik, L. T. Ter-Martirosyan, and S. P. Zubko, “Microwave losses in incipient ferroelectrics as functions of the temperature and the biasing field,” J. Appl. Phys. **84**(2), 993–998 (1998). [CrossRef]

**46. **G. Allaire, “Homogenization and two-scale convergence,” SIAM J. Math. Anal. **23**(6), 1482–1518 (1992). [CrossRef]

**47. **S. Guenneau and F. Zolla, “Homogenization of three-dimensional finite photonic crystals,” Prog. Electromagn. Res. **27**, 91–127 (2000). [CrossRef]

**48. **M. Otomori, J. Andkjaer, O. Sigmund, K. Izui, and S. Nishiwaki, “Inverse design of dielectric materials by topology optimization,” Prog. Electromagn. Res. **127**, 93–120 (2012). [CrossRef]

**49. **S. Molesky, Z. Lin, A. Y. Piggott, W. Jin, J. Vucković, and A. W. Rodriguez, “Inverse design in nanophotonics,” Nat. Photonics **12**(11), 659–670 (2018). [CrossRef]

**50. **N. M. Estakhri, B. Edwards, and N. Engheta, “Inverse-designed metastructures that solve equations,” Science **363**(6433), 1333–1338 (2019). [CrossRef]

**51. **J. Andkjær and O. Sigmund, “Topology optimized low-contrast all-dielectric optical cloak,” Appl. Phys. Lett. **98**(2), 021112 (2011). [CrossRef]

**52. **B. Vial and Y. Hao, “Topology optimized all-dielectric cloak: Design, performances and modal picture of the invisibility effect,” Opt. Express **23**(18), 23551 (2015). [CrossRef]

**53. **B. Vial, M. M. Torrico, and Y. Hao, “Optimized microwave illusion device,” Sci. Rep. **7**(1), 3929 (2017). [CrossRef]

**54. **Z. Lin, V. Liu, R. Pestourie, and S. G. Johnson, “Topology optimization of freeform large-area metasurfaces,” Opt. Express **27**(11), 15765–15775 (2019). [CrossRef]

**55. **Z. Lin and S. G. Johnson, “Overlapping domains for topology optimization of large-area metasurfaces,” Opt. Express **27**(22), 32445 (2019). [CrossRef]

**56. **A. R. Diaz and O. Sigmund, “A topology optimization method for design of negative permeability metamaterials,” Struct Multidisc Optim **41**(2), 163–177 (2010). [CrossRef]

**57. **O. Sigmund, “Systematic design of metamaterials by topology optimization,” in IUTAM Symposium on Modelling Nanomaterials and Nanosystems, R. Pyrz and J. C. Rauhe, eds. (Springer, 2009), IUTAM Bookseries, pp. 151–159.

**58. **S. Zhou, W. Li, Y. Chen, G. Sun, and Q. Li, “Topology optimization for negative permeability metamaterials using level-set algorithm,” Acta Mater. **59**(7), 2624–2636 (2011). [CrossRef]

**59. **S. Nishi, T. Yamada, K. Izui, S. Nishiwaki, and K. Terada, “Isogeometric topology optimization of anisotropic metamaterials for controlling high-frequency electromagnetic wave,” Int. J. Numer. Meth. Engng. **121**(6), 1218–1247 (2020). [CrossRef]

**60. **B. S. Lazarov and O. Sigmund, “Filters in topology optimization based on Helmholtz-type differential equations,” Int. J. Numer. Meth. Engng. **86**(6), 765–781 (2011). [CrossRef]

**61. **F. Wang, B. S. Lazarov, and O. Sigmund, “On projection methods, convergence and robust formulations in topology optimization,” Struct Multidisc Optim **43**(6), 767–784 (2011). [CrossRef]

**62. **A. Kawamoto, T. Matsumori, S. Yamasaki, T. Nomura, T. Kondoh, and S. Nishiwaki, “Heaviside projection based topology optimization by a PDE-filtered scalar function,” Struct Multidisc Optim **44**(1), 19–24 (2011). [CrossRef]

**63. **M. Bendsøe and O. Sigmund, “Material interpolation schemes in topology optimization,” Arch. Appl. Mech. Ing. Arch. **69**(9-10), 635–654 (1999). [CrossRef]

**64. **M. Alnæs, J. Blechta, J. Hake, A. Johansson, B. Kehlet, A. Logg, C. Richardson, J. Ring, M. E. Rognes, and G. N. Wells, “The FEniCS Project Version 1.5,” Arch. Numer. Softw. **3**(100), 20553 (2015). [CrossRef]

**65. **M. B. Giles and N. A. Pierce, “An introduction to the adjoint approach to design,” Flow, Turbul. Combust. **65**(3/4), 393–415 (2000). [CrossRef]

**66. **S. K. Mitusch, S. W. Funke, and J. S. Dokken, “Dolfin-adjoint 2018.1: automated adjoints for FEniCS and Firedrake,” J. Open Source Softw. **4**(38), 1292 (2019). [CrossRef]

**67. **K. Svanberg, “A class of globally convergent optimization methods based on conservative convex separable approximations,” SIAM J. Optim. **12**(2), 555–573 (2002). [CrossRef]

**68. **S. G. Johnson, “The NLopt nonlinear-optimization package,” (2014).

**69. **P. Dular, C. Geuzaine, F. Henrotte, and W. Legros, “A general environment for the treatment of discrete problems and its application to the finite element method,” IEEE Trans. Magn. **34**(5), 3395–3398 (1998). [CrossRef]

**70. **V. Hernandez, J. E. Roman, and V. Vidal, “SLEPc: a scalable and flexible toolkit for the solution of eigenvalue problems,” ACM Trans. Math. Softw. **31**(3), 351–362 (2005). [CrossRef]