Complex diffraction limited acoustic fields can be generated from a single element transducer using inexpensive 3-D printable acoustic kinoforms. This is extremely promising for a number of applications. However, the lack of ability to vary the field limits the potential use of this technology. In this work, this limitation is circumvented using multi-frequency acoustic kinoforms for which different acoustic fields are encoded onto different driving frequencies. An optimisation approach based on random downhill binary search is introduced for the design of the multi-frequency kinoforms. This is applied to two test cases to demonstrate the technique: a kinoform designed to generate the numerals “1,” “2,” and “3” in the same plane but at different driving frequencies, and a kinoform designed to generate 3 sets of eight foci lying on a circle with a driving-frequency-dependent radius. Field measurements from these samples confirmed that multi-frequency acoustic kinoforms can be designed that switch between different arbitrary, pre-designed, acoustic field patterns in the target plane by changing the driving frequency.

The precise manipulation of acoustic fields in 3-D is vital to a wide range of applications in physical acoustics.^{1–3} Recently, a new approach for the generation of complex diffraction limited acoustic fields in 3-D was introduced.^{1} This works by mapping the continuous wave output of a planar single element transducer onto a pre-determined phase distribution which then diffracts to form a target acoustic field. The mapping is accomplished using a 3-D printed kinoform attached to the front face of the transducer. Differences in sound speed between the kinoform and the coupled medium result in a phase change proportional to the thickness of the structure at each point.

The kinoforms are extremely cheap to fabricate, forego the complicated driving electronics required for 2-D arrays, and scale easily to MHz frequencies. However, one drawback is that the acoustic field generated by each kinoform is fixed. This limitation can potentially be overcome using multi-frequency kinoforms for which different target patterns are encoded onto different frequencies. The field at a particular target depth could then be varied simply by altering the driving frequency of the transducer. The goal of this work was to demonstrate the use of multi-frequency kinoforms in acoustics.

Within optical holography, multi-frequency kinoforms are of interest for a variety of applications, particularly colour displays. This has led to a range of algorithms for their design. Two simple approaches are spatial and depth division.^{4} In the former, the aperture of the kinoform is sub-divided into separate regions for each frequency.^{5} For the latter, a kinoform is calculated for a single frequency with target images at different depths using a multi-plane iterative Fourier transform algorithm (IFTA). The focal distance of kinoforms in the nearfield scales inversely with wavelength, so different patterns are shifted into the target plane as the frequency is changed.^{6} Other works have used approaches that search directly for a thickness at each position that maximises a quality-metric based on the target field for each frequency. Several methods have been used to solve this optimisation problem including direct search,^{7,8} optimal-rotation-angle (ORA),^{9–11} and the Yang-Gu algorithm for phase retrieval.^{12}

To design the kinoforms in this work, the direct search approach outlined in Kim *et al.* was adapted.^{7} The implementation was partly based on our earlier work, in which surface profiles for the generation of arbitrary optoacoustic fields were designed with a random downhill binary search.^{13}

For a kinoform with a spatially varying thickness described by *h*(*x*, *y*) attached directly to a planar transducer driven at a frequency *f*, the acoustic field in a plane *p*(*x*, *y*) directly after propagation through the kinoform can be approximated by $p(x,y)=A(x,y)ei\Delta \varphi (x,y)$. Here, *A*(*x*, *y*) is the surface pressure of the transducer and Δ*ϕ*(*x*, *y*) is the phase offset introduced by propagation through the kinoform. This phase offset is related to the thickness by

where *k _{h}* and

*k*are the wavenumbers in the kinoform and coupled medium, respectively, and

_{c}*h*is the maximum thickness of the kinoform. This model approximates the kinoform as an infinitely thin phase screen and neglects transmission losses.

_{m}Assuming that the transducer is placed in the *z* = 0 plane, the acoustic field *p*(*x*, *y*, *z*) for *z* > 0 can then be calculated using the Rayleigh-Sommerfeld integral^{14}

Here, Ω defines the area over which the transducer output is non-zero, *d*Ω is a planar area element (in the *z* = 0 plane), and R is the distance between the coordinates (*x*, *y*, *z*) and $(x\u2032,y\u2032,0)$.

The inverse problem is to find a kinoform *h*(*x*, *y*) that for a set of design frequencies *f _{n}*, $n\u2208\mathbb{N}$, maps the field of the transducer onto a set of frequency-specific target pressure distributions

*p*(

_{n}*x*,

*y*,

*z*). Here,

_{d}*z*is the depth of the target plane and

_{d}*n*is the index of the driving frequency. To solve this, first the surface coordinates of the kinoform are discretised such that a square transducer occupying an L × L region of the (

*x*,

*y*) plane is represented by a set of N × N square pixels with sides of length $\Delta =LN$. For other transducer shapes (i.e., circular), the set of squares is restricted to those that lie within the area over which the transducer output is non-zero. The surface height is also discretised such that the height at each coordinate (

*x*,

_{i}*y*) is given by $hij=m\Delta z,\u2009m\u2208\mathbb{N}$ and

_{i}*m*Δ

_{z}<

*h*. Here,

_{m}*ij*is the coordinate index.

^{13}

For this work, rather than define continuous target distributions *p _{n}*(

*x*,

*y*,

*z*), sets of target points

_{d}*r*= {(

_{n}*x*,

_{m}*y*,

_{m}*z*),

_{d}*m*= 1,…,

*M*}

_{n}were defined for each frequency. The optimisation then searched for a map

*h*that maximised the pressure across the different sets of target points by maximising the following cost function:

_{ij}Here, *p _{n}* is the value of the pressure at a set of target points

*r*as defined above, evaluated for a driving frequency

_{n}*f*,

_{n}*S*is a smoothness term evaluated using a second-order accurate finite difference approximation of the Laplacian,

^{13}and

*α*and

*β*are the weighting parameters. The first term $E(|pn|)$ is the average pressure generated over the target points for each design frequency when the transducer is driven at that design frequency. The second term is the standard deviation of the same parameter. The first term maximises the average amplitude over the target pattern for each frequency, while the second term penalises variations in the pressure across the target patterns. The third term is a penalty on sharp isolated variations in the kinoform height. The forward model does not account for reflections from the kinoform structure, and hence, this penalty mitigates rapid variations that could lead to more reflections.

The kinoform was initialised with a uniform height $hij0=hm/2$. From kinoform $hijk$, the optimisation proceeded by selecting a pixel at random (from those not yet updated) and assigning it the height that maximised the cost function in Eq. (3). This was then repeated until all pixels had been updated, resulting in the updated kinoform $hijk+1$. The optimisation was terminated when the number of pixels changed from the previous iteration fell below 0.5% of the total pixels on the surface.

If the smoothness term S was weighted inappropriately, then the optimisation converged to kinoforms that were either flat or had large height differences between adjacent pixels. To find a suitable weighting, the optimisation was run once to convergence with *β* = 0. The output was then used to initialise a second instance of the optimisation with *β* chosen such that $E(|pn|)=\u22120.5S$ for the initial kinoform state. The variable *α* was assigned a value of 0.3.

The approximation of the kinoform as a thin phase element in the forward model can be inaccurate for substantial element thicknesses and large height variations.^{15–17} As a result, significant differences can occur between the prediction of the optimisation and the actual field generated by the kinoform. Ideally, a full-wave model such as the k-Wave toolbox^{18} would be used to more accurately model the field. However, at present, such models are too slow to incorporate into the optimisation. In practice, to ensure the optimisation converged to valid results, the maximum thickness *h _{m}* was constrained.

Two test cases were fabricated to demonstrate the realisation of multi-frequency acoustic kinoforms. The first was designed to generate the numerals “1,” “2,” and “3” at three distinct frequencies. The second was designed to generate 3 sets of 8 discrete foci distributed evenly on circles of radii 10.5, 7, and 3.5 mm. These were selected to demonstrate that both patterns of pressure and discrete foci can be created with this approach.

Both kinoforms were designed for a 3.1 cm PZT piston transducer (Olympus, Japan) with a centre frequency of 2.25 MHz and a −6 dB bandwidth between 1.9 and 3.1 MHz. The image depth *z _{d}* in all cases was 2.5 cm, and both kinoforms were discretised using a spacing Δ and Δ

_{z}of 100

*μ*m. For calculating the wavenumber

*k*, the kinoform sound speed was assumed to be 2495 ms

_{h}^{−1}.

^{13}The maximum thickness

*h*was set to 2 mm to provide a 2

_{m}*π*modulation for the lowest design frequency. This was calculated by inverting Eq. (1) and using 1500 ms

^{−1}as the sound speed for the medium.

The design frequencies for the three patterns were set to 1.9, 2.5, and 3.1 MHz to maximise the spacing within the usable bandwidth. This was necessary to minimise crosstalk, the unwanted generation in the field at one design frequency of distributions encoded onto other design frequencies. This occurs as, for low *h _{m}*, the phase offsets introduced for different wavelengths are similar, leading to similar diffraction patterns. Fortunately, due to the chromatic variation in focal distance for near-field kinoforms, the depth at which the crosstalk is generated differs from the target depth. This depth scales approximately with the ratio of the design and driving frequencies, hence the need to maximise the frequency spacing to separate the crosstalk and target planes.

^{9,17,19}

The calculated kinoforms were printed using a high-resolution PolyJet printer (Objet350 Connex, Stratasys, Eden Prairie, MN, USA) using VeroBlack as the substrate. A 1 mm thick backing was added to both kinoform designs for structural rigidity. Additionally, a ring was added to the outside to enable the kinoform to slot onto the front face of the transducer. A rendering of the kinoform design and the fabricated kinoform for the 3 pattern sample can be seen in Figs. 1(a) and 1(b).

The acoustic fields generated by the kinoforms were measured in a 40 × 40 × 60 cm test tank with a two axis computer controlled positioning system (Precision Acoustics, Dorchester, UK). Measurements were made using a calibrated 0.2 mm needle hydrophone (Precision Acoustics, Dorchester, UK). The transducer was suspended inside the tank and aligned such that its front face was parallel to the measurement plane. The needle hydrophone was then suspended ∼5 cm from the transducer aligned with its beam axis.

The transducer was driven using a signal generator (33522A, Agilent Technologies, Santa Clara, CA, USA) connected via a 75W power amplifier (A075, E&I, Rochester, NY, USA). The driving signal in each case was a sinusoid at one of the design frequencies. The number of cycles was varied between 35 and 56 to ensure that the signal covered the minimum and maximum travel times from the kinoform to the measurement plane. The driving voltage varied from 40 to 60 V peak-to-peak for different scans. The field was recorded over a 6 × 6 cm area centred on the middle of the kinoform. This was done for both kinoforms and for each design frequency. The step size was varied from 0.24–0.3 mm to ensure that a *λ*/2 sampling criterion was met for each driving frequency. Time domain signals were recorded at each position using a digital oscilloscope with a sampling rate of 400 MHz and 32 averages. To calculate the 3-D wavefield from each measurement set, the k-Wave toolbox was used to back-propagate them. Prior to projection, each dataset was upsampled by a factor of two and bandpass filtered between 100 kHz and 4 MHz to reduce noise.

The output field of each kinoform for each driving frequency was also simulated to compare with the experimental measurements. These simulations were carried out using the k-Wave toolbox^{18} on a 360 × 360 × 360 grid with a 100 *μ*m point spacing. The kinoform designs were inserted into the simulation as volumes with a sound speed of 2495 ms^{−1} and a density of 1190 kg m^{−3}.^{13} The transducer was modelled by a circular source mask of diameter 3.1 cm driven at a single design frequency *f _{n}*. In the background medium, the acoustic properties were set to those of water. The maximum pressure generated in the target plane

*z*for each design frequency in both simulation and experiment can be seen in Figs. 2(a)–2(f) for the continuous pattern kinoform and in Figs. 3(a)–3(f) for the discrete foci kinoform.

_{d}For both kinoforms, the target distributions are clearly realised in the experimental data. This demonstrates that using multi-frequency kinoforms, different complex pressure distributions can be generated by switching the driving frequency. Additionally, the variation in pressure across the pattern and in the background is similar between the simulations and experiments. This confirms that a full-wave model can be used to accurately predict the fields prior to fabrication, suggesting that better results could be obtained by using a full-wave model in the optimisation. However, this would come at significant computational cost [for the kinoforms in this work, the forward model was run approximately 30 million times.] The background for both datasets also shows some crosstalk due to the constraints on the frequency spacing applied by the transducer bandwidth.

The variation in crosstalk with frequency spacing was quantified in more detail numerically. A set of kinoforms was calculated, each designed to generate the numerals “1” and “2” at a depth of 2.5 cm. The design frequency for the “1” was fixed at 2 MHz, and the design frequency for the “2” was increased in increments of 100 kHz from 2 to 6 MHz. Each was calculated for a 3.1 cm aperture using a spacing of 100 *μ*m and a *h _{m}* value of 2 mm. The field generated by each kinoform at its design frequencies was then simulated using the k-Wave toolbox. To evaluate the crosstalk, the images of the “1” and “2” were used to define a foreground and background for both patterns. The crosstalk was calculated as the ratio of the foreground pressure to background pressure for the “1” in the field at the design frequency of the “2” and vice-versa. This was evaluated at all depths for each kinoform and the plane with the highest value selected, along with the value at the target depth. These were each averaged across both design frequencies to give a single number.

The resulting variation in the peak crosstalk and in-plane crosstalk as a function of frequency spacing can be seen in Fig. 4. For comparison, the contrast (ratio of foreground to background pressure) measured across the target distribution for each kinoform is also shown. Figure 4 shows that as the design frequency ratio increases, the crosstalk decreases, falling to 33% of the target contrast for a ratio of 3. This is because the phase offsets introduced for the design frequencies become less similar as the spacing increases. The in-plane crosstalk decreases more rapidly, converging to a value of ∼2 by a frequency ratio of 1.5. This occurs as the chromatic variation in focal distance increasingly separates the crosstalk plane from the target depth, eliminating the in-plane distortion.

The ability of multi-frequency kinoforms to encode multiple distributions of pressure is also affected by several other parameters. As additional patterns are added, the solution quality drops. This is because the kinoforms need to focus over a greater number of target points using the same space-bandwidth product or number of pixels.^{1,7} Similarly, increasing the image complexity or area decreases the contrast. The complexity of features that can be realised is also limited by diffraction for each wavelength. Conversely, the solution quality improves as the aperture size and modulation depth or *h _{m}* increase. However, as discussed above, the accuracy of the forward model decreases as

*h*increases.

_{m}In addition to the relative spacing of the design frequencies, the choice of frequencies also has an impact. For higher frequencies, the performance is limited by acoustic attenuation and the resolution of the 3-D printing process. At lower frequencies, the ability to realise an arbitrary pressure distribution is constrained by the diffraction limit. Figure 5 illustrates the effect of frequency and frequency spacing on the field more clearly. It shows the field of the continuous pattern kinoform used in the experiment compared to that of a kinoform designed to generate the same patterns using design frequencies of 2.5, 4, and 6 MHz. The reduction in crosstalk and significant improvement in the fidelity are apparent.

This work has shown that multi-frequency kinoforms can be implemented in acoustics for the formation of different target distributions at a common target depth. This could have applications in several areas of acoustics where the manipulation of acoustic pressure in 3-D is required, such as particle manipulation and haptic displays. This work could be extended in several ways. Different cost functions could be incorporated into the optimisation to control the phase of the target field in addition to the amplitude. The forward model could also be altered to incorporate a full-wave model for certain problems, for example, axisymmetric distributions where the number of optimisation variables is sufficiently small that the number of simulations required becomes feasible.

This work was funded by the Engineering and Physical Sciences Research Council and the UCL EPSRC Centre for Doctoral Training in Medical Imaging.