Acoustic equipmentThree single-element flat transducers of various dimensions and frequencies were used in this study. These commercially available flat transducers (American Piezo Co., USA) have active elements of ~0.9–1.1 mm thicknesses, with OD = 50, 35 and 25 mm with center frequency, f0, of 2.28, 1.86 and 2.24 MHz, respectively. Each transducer was encapsulated with aluminum housing, as provided by the manufacturer. The efficiency of the transducers in converting the electrical power to acoustic power was around 49%, as reported by the manufacturer. To drive the transducers the power generator device TPO-102 (Sonic Concept Inc., USA) with 4–210 Watt power range and 50 Ω output impedance was used with a built-in sin function generator. A step-down electrical matching unit (model #90-4496, American Piezo Co., USA) to connect the transducers with 20 Ω impedance to the power generator with 50 Ω was employed to ensure the best electrical matching between the transducer and the power generator. To minimize the distortion of the passing wavefront and reflection as well as maintaining the acoustic power transmission, a thin shell-type barrier separates the water medium from the built chamber. Moreover, to maintain a good energy transfer and balanced coupling between the transmission hologram and the transducer, we used mechanical clamping (bolting) accompanied by the epoxy paste matching layer, ensuring the best combination for acoustic matching. The conductive thermal paste is homogeneously distributed to the surface of the phase plate by spin coating (WS-400B/Lite) at ~3500 RPM. In order to acquire the pressure patterns a needle hydrophone with 0.2 mm diameter (model no. D1602) along with the hydrophone booster amplifier (Precision Acoustics Ltd., UK) is used.The speed of sound in deionized (DI) water at the laboratory temperature as ~1480 m/s and the transmission phase plate material as ~2430 m/s was measured using through transmission technique56,57 with needle hydrophone while the speed of sound in PDMS with different mixing ratios was obtained via high-speed images to be used for the further finite element simulations. Mostly we used 10:1 mixing ratio of the base PDMS to the curing agent, and the speed of sound obtained was ~1040 m/s. Furthermore, we performed the pressure mapping with an in-house setup including an acrylic test tank of 1 × 0.75 × 0.75 m dimension filled with DI water with a needle hydrophone mounted on a computer-controlled positioner. An oscilloscope (Tektronix DPO2400) was used to acquire the needle hydrophone signal which was then connected to MATLAB R2020b to store and post-process the signals. The positioner movement and the data acquisition from the needle hydrophone were synchronized with the computer and triggered by the input burst pulse to the transducer.Printing material preparationSYLGARD 184 kit (Krayden, Canada) is used as the printing material in this paper. The PDMS system (monomer and curing agent) was mixed with the defined mixing ratios. The mixture was degassed using a vacuum for 45 min. The PDMS system was colored by adding the 2% w/w oil-based dye (Winsor & Newton, UK) to the base monomer.Holography-based SCL experiment and high-speed imagingA 1 mM solution of luminol (3-aminophthalhydrazide, Sigma-Aldrich, Canada) is prepared, and the pH is adjusted to 12 using NaOH (Sigma-Aldrich, Canada). The pH is continuously monitored in real time using a pH 315i meter (Wissenschaftlich-Technische Werkstätten, WTW, GmbH, Germany). Subsequently, 0.5 M sodium carbonate (Na2CO3) is added to the solution. The sodium carbonate is produced by heating sodium bicarbonate (NaHCO3, commonly known as baking soda) to 100 °C.Luminol solution is poured inside a build volume that is used for printing material container and placed in front of the flat transducer with a hologram attached with a distance of the designed object plane (20 mm). To reduce the ambient light for capturing SCL images while sonication, we performed the experiment inside a dark room. Single lens EOS 500D (Canon, Japan) camera is used for capturing images of sonochemical reactions. The exposure time of the DSLR camera was set to a maximum 30 s to ensure capturing the blue SCL light from desired geometries.High-speed imaging performed with Fastcam SA-Z (Photron Inc, USA) high-speed camera. The high-speed imaging experiment presented in Fig. 3 was recorded with 6000 frames/s and a 6× magnification lens. To improve the brightness of the footages, Halogen light source is employed through all the experiment.Volume deposition rate in HDSP vs. DSPThe polymerization rate within our HDSP system can be correlated with acoustic power and pressure, although with certain complexities inherent to the process. Importantly, in HDSP, the polymerization rate is influenced not only by the delivered acoustic power but also by the DC of the transmitted acoustic power applied to the resin.To quantify this relationship in our system, we conducted a series of experiments in which we systematically varied both the acoustic power and the DC delivered to the resin. We have documented these variations and their impact on the printing time required to generate a solid wall by extruding a holographically generated line. This is detailed in Supplementary Fig. 3c, which illustrates the combination of power and DC with the printing time.We have conducted a more detailed analysis and included additional data in Supplementary Fig. 4. This figure illustrates our parametric study on printing a wall with dimensions of 15 × 1 × 20 mm3, as depicted schematically in Supplementary Fig. 4d.We have expanded our analysis to include a new representation that specifically illustrates the relationship between printing power and printing time when the DC is varied, as shown in Supplementary Fig. 4a. We observed that increasing the DC results in extended printing times. Additionally, we introduce the concept of interaction strength, σ (W), calculated by multiplying the acoustic power (W) and the DC (%). Supplementary Fig. 4b demonstrates how this interaction influences the printing time, allowing us to derive corresponding polymerization rates under these conditions.This analysis enables us to establish a quantifiable correlation between acoustic power, DC, and the polymerization rate, in terms of a volumetric deposition rate (VDR) function, given in Eq. 1. Our findings indicate that the total power required to print geometries of varying circumferential lengths differs, even though the intensity needed for printing each voxel remains constant. Consequently, to enhance the versatility of the VDR function presentation, we incorporated the intensity used in the parametric study, calculated as power divided by the circumferential area of the voxels targeted by acoustic holography for solidification. Supplementary Fig. 4c shows the plot of this numerically derived VDR function, modeled using a second-order polynomial (with R2 = 0.9564):$$\begin{array}{c}{{{\rm{VDR}}}}\,(I,{{{\rm{DC}}}})={a}_{0}+{a}_{10}I+{a}_{01}{{{\rm{DC}}}}+{a}_{20}{I}^{2}\\ \quad \quad \quad+{a}_{11}I.{{{\rm{DC}}}}+{a}_{02}{{{\rm{DC}}}}^{2}\\ {{{\rm{where}}}}:\,\left\{{a}_{0}=1.24\times {10}^{4}\,,\,{a}_{10}=-2.88\times {10}^{4}\,,\,{a}_{01}=-725.89,\right.\\ \left.{a}_{20}=1.96\times {10}^{4}\,,\,{a}_{11}=982.74\,,\,{a}_{02}=11.53\right\}\end{array}$$
(1)
where \(I=\frac{P}{{A}_{{{\rm{ROI}}}}}\) is the acoustic intensity on the image with area of AROI and delivered acoustic power P. The VDR function is based on two key input parameters, intensity of the printing, I (W/mm2) and DC. The results are clearly depicted in Supplementary Fig. 4c, which shows nonlinear traits. The quadratic nature of the numerically obtained VDR function is consistent with the quadratic-like trends observed in Supplementary Fig. 4a, b.The time saved by HDSP due to stationary source, and not moving the focal point as in DSP, translates directly into faster overall printing speeds. The ability to solidify a large area at once, rather than point-by-point, is a clear advantage in scenarios where production speed is a critical factor. We have experimentally compared the printing of a wall of 15 × 1 × 20 mm3 obtained by the two methods. Since we know the total volume of the printed part, and by selecting a modest power/feed rate setting for both methods of DSP and HDSP, the time taken for complete printing can be measured. The printing time for HDSP using 30 W and DC 50% (electrical input) was ~30 s while for DSP with 240 mm/min feed rate, the printing took 12.56 min. Comparing the printing time for HDSP vs. DSP, we can conclude that since the HDSP has more voxels to simultaneously print, it requires much higher energy compared to the DSP, and that’s why it is faster. The printing time of HDSP is calculated using VDR and the printing time for DSP is path length divided by the feed rate of the source. Therefore, printing time can be experimentally obtained for the HDSP vs. DSP. Supplementary Table 1 shows this comparison in detail.3D scanning and accuracy measurement of hologramsWe used Formlabs 2.0 and 3.0 to fabricate the holograms. High-temperature as well as Clear resin was used with acceptable speed of sound measured between ~2650 and ~2430 m/s, respectively. The advantage of incorporating high-temperature resin is that it can withstand the temperature increase caused by the high-intensity passing wave, without deforming and damaging the whole structure of the holograms.To assess the accuracy of the fabricated thickness pattern on the transmission holograms, we used 3D measuring laser Confocal Microscopy (Olympus, Lext OLS4100, Japan) with measuring accuracy down to 80 µm. The point clouds were then imported to CATIA V5 followed by deviation analysis with the original theoretical hologram 3D model.Hologram calculation and COMSOL 3D simulationsThe retrieved phase pattern required for the acoustic holography using transmission phase hologram was calculated with the IAS method35, performed in Matlab2020b. The phase pattern was then converted to the thickness map of each hologram which was then fabricated with Formlabs 2.0 SLA printer. Next section briefly describes a protocol for generating a required acoustic field and printing for HDSP. To validate the accuracy of the projected field by holograms, 3D finite element simulations were conducted in COMSOL Multiphysics 6.0, acoustic module through a high-performance computing center with 32core CPU and 500GB RAM.HDSP process flowThe working flow of the HDSP to print an object (Supplementary Fig. 1) begins by retrieving the phases necessary for creating the complex holographic field. This is done by preparing binarized image(s) of the cross-section of the objects to be printed, denoted as \({R}^{j}\) (j indicating the index of each image) placed at \(z={z}_{t}^{j}\). These images are used to build the complex pressure field as:$${p}_{t}^{j}(x,\, y,\, z)={R}^{j}(x,\, y,\, z){e}^{i{\phi }_{t}^{j}(x,\, y,\, z)}$$
(2)
where \({\varphi }_{t}^{j}\) is initially set to a random value. Then angular spectrum method takes the complex pressure field in the spectral domain of \({{{\bf{k}}}}=({k}_{x},\, {k}_{y},\, {k}_{z})\) which can be found by 2D Fourier transform as Eq. 3:$${P}_{t}^{j}({k}_{x},\, {k}_{y},\, {z}_{t}^{j})={\iint }_{\infty }{p}_{t}^{j}(x,\, y,\, {z}_{t}^{j}){e}^{-j({k}_{x}x+{k}_{y}y)}{{{\rm{d}}}}x{{{\rm{d}}}}y$$
(3)
and employs the propagator function H,$$H({k}_{x},\, {k}_{y},\, \pm {z}_{t}^{j})={e}^{i{k}_{z}(\pm {z}_{t}^{j})}$$
(4)
to propagate it back to the hologram plane at \(z=0\):$${P}_{h}^{j}({k}_{x},\, {k}_{y},\, 0)={P}_{t}^{j}({k}_{x},\, {k}_{y},\, {z}_{t}^{j})H({k}_{x},\, {k}_{y},\, -{z}_{t}^{j})$$
(5)
where \({k}_{z}=\sqrt{{k}^{2}-{k}_{x}^{2}-{k}_{y}^{2}}\), with k being wavenumber in the medium. In Eq. 3, a positive sign corresponds to the forward propagation toward the target plane, while a negative sign corresponds to the backward propagation toward the hologram plane. Finally, the 2D inverse Fourier transforms the complex pressure from the spectral domain to the spatial domain as:$${p}_{h}^{j}(x,\, y,\, 0)=\frac{1}{4{\pi }^{2}}{\iint }_{\infty }{P}_{h}^{j}({k}_{x},\, {k}_{y},\, 0){{{\rm{d}}}}{k}_{x}{{{\rm{d}}}}{k}_{y}$$
(6)
The phase map of the backpropagated complex field on the hologram plane (Eq. 6) can be retrieved as \({\phi }_{h}^{j}\). Next, the retrieved \({\phi }_{h}^{j}\) and imposed constrains at the source amplitude are used to continue the iteration and reconstruct the target object cross-section image, \({A}^{j}\). The algorithm iteratively propagates the complex field back and forth between the object plane at \(z={z}_{t}^{j}\) and the source plane at \(z=0\), while imposing the source and desired target amplitude, so that the reconstructed pressure image \({A}^{j}\) converges to the reference pressure image\({R}^{j}\). In each iteration step, the reconstructed image quality check will ensure the minimum error with the true image of the object. The evaluating criteria for measuring the quality of the reconstructed image can be mean squared error (MSE), NMSE, SNR and Correlation (see section “Image quality analysis”). Finally, the total phase map \({\varPhi }_{h}(x,\, y)\) on the hologram plane necessary for the complex holographic field can be obtained by summation of all the retrieved phase maps, \({\phi }_{h}^{j}\). After obtaining the phase map, the required complex acoustic field can be fine-tuned in its two constituents, the required power and the obtained phases. The acoustic power required for igniting the solidification can be adjusted by the pressure amplification, PA, which can be obtained theoretically for each hologram as \({{{\rm{P{A}}}}}^{j}\). To generate a complex geometry, the input power \({P}_{{{\rm{in}}}}\) must be adjusted in accordance with the experimentally determined threshold power required for generating a single voxel, \({P}_{{{\rm{thrsh,\, voxel}}}}\). This input power can be found using Eq. 7, which was derived from experimental data.$${P}_{in}={P}_{thrsh,\, {{{\rm{voxel}}}}}\times \left(\frac{P{A}_{voxel}}{P{A}^{j}}\right)$$
(7)
Next, the phase map \({\varPhi }_{h}\) can be converted to thickness map and fabricated using \(h={\varPhi }_{h}/({k}_{m}-{k}_{h})\), where \({k}_{h}\) and \({k}_{m}\) correspond to wavenumber in hologram material and medium, respectively. In this step, to prevent the formation of sudden artifacts, we bound the \({\varPhi }_{h}\) values to be within 0–2π. Finally, the complex acoustic field will be ready to ignite the printing within the region of interest. Supplementary Fig. 1 illustrates the process flow of the HDSP.Material characterizationTo evaluate the effecting parameters on the printing process of HDSP, further characterization along with input parametric analysis was performed. A line object created by HDSP and extruded 10 mm in z-direction to form a wall, was chosen for this experiment. From the parametric analysis of the printed walls, it was observed that the dominating input parameters affecting the printing process are, input power, DC% and feed rate. Supplementary Fig. 3c presents the relation between power, printing time, DC% and microstructure width form within the printed walls. First of all, increasing the power and DC% separately led to decreasing printing time, since more energy will be delivered to the printing ROI, while causing the denser bubble cloud formation in some locations within the printing object due to increased cavitation pressure. However, decreasing the DC% was proved to diminish the presence of cavitation bubble clouds and results in more transparent objects. As such comparison, the image on the right side of Supplementary Fig. 3b shows a part with bubble cloud pillars formed within the structure printed with 50% DC and a few seconds printing time, while the image on the top depicts a fully transparent wall with the same power and 14% DC but relatively longer printing time. The DC% is a relative term expressing the ratio of burst signal duration to the interval time between consecutive bursts. Supplementary Fig. 3d shows a box plot depicting the difference in measured width of the bubble cloud pillars formed under conditions of 2.5 and 5 ms intervals, both with an average power of 6 W. The regression line represents the trend observed in the data. Hence, altering the interval of the input signal could be another affecting parameter to control the cavitation bubble cloud formation inside the printing part.Image quality analysisIn order to have a precise comparison of the generated images of the printing object and the original image, a frequency range of 1.5–4.5 MHz with OD of transducer between 25–100 mm was considered. First, the original binary image of the desired object to be printed with a width of 20 mm at the base point of evaluation with f0 = 1.5 MHz and OD = 2 5 mm was selected with each pixel value one-fifth of the wavelength with 126 × 126 size of image. By conserving the object width for each case through padding and resizing the image matrix, a reconstructed image of the desired object can be achieved. Next, the obtained images were resized back to the 126 × 126 pixels using interpolation without reducing the content followed by normalization, where can be compared with the base point image.The equation for the Correlation between the two images can be found by Eq. 858:$$\begin{array}{c}{{{\rm{Correlation}}}}=\\ \frac{{\sum}_{j}{\sum}_{i}({R}_{i,\, j}-\bar{R})({A}_{i,\, j}-\bar{A})}{\left({\left({\sum}_{j}{\sum}_{i}{({R}_{i,\, j}-\bar{R})}^{2}\right)\left({\sum}_{j}{\sum}_{i}{({A}_{i,\, j}-\bar{A})}^{2}\right)}\right)^{0.5}}\end{array}$$
(8)
where R is the reference image and A is the reconstructed image. The bar represents the mean value of the image. The NMSE formula can be represented as Eq. 942:$${{{\rm{NMSE}}}}=\frac{{||R-A||}_{2}^{2}}{{||R||}_{2}^{2}}$$
(9)
Other than evaluating parameters discussed in Fig. 6d, e, PSNR and SSIM59 assist to compare the generated image and actual image of the printing object, more accurately. Hence, according to the PSNR equation given in Eq. 10 as:$${{{\rm{PSNR}}}}= 20{\log }_{10}\frac{{A}_{\max }}{\sqrt{{{\rm{MSE}}}}}\,,\, \\ {{{\rm{MSE}}}}= \frac{1}{m}\sum \sum ({A}_{i,\, j}-{R}_{i,\, j})^2$$
(10)
Supplementary Fig. 5a, b illustrates the PSNR and normalized SSIM, demonstrating the direct relation of the frequency, f0 and source aperture, OD, on image quality for the 20 mm image to be generated at 20 mm distance target plane. By increasing the frequency and transducer OD, the PSNR of the obtained images of the object improves and structural similarity between the desired image and reconstructed image increases. Moreover, by keeping the source aperture and frequency constant, the distance of the desired object plane from the source can highly affect the resolution of the object. Supplementary Fig. 5c blue shows the dependency of the object plane location to the SNR for the 50 mm transducer with 2.28 MHz frequency. As can be seen for ztarget range between 16 and 20 mm, we can expect better SNR. Also, the PSNR of the reconstructed image vs. target plane position can be seen in Supplementary Fig. 5c red for the 50 mm transducer and 2.28 MHz frequency. The better SNR and PSNR can be achieved with a 16–20 mm distance of the object plane from the source. Hence, improved peak amplitude and less noise and side lobes can be achieved by wise selection of the effective parameters.Robot’s trajectory in HDSPThe robot’s end effector, where the printing platform is attached to, trajectory and velocity need to be calculated based on the nature of the printing mechanism in HDSP using a robotic arm (Supplementary Fig. 7a). The printing location in HDSP, which is the image plane, is stationary while the printing platform (end effector) moves in the working space carrying already printed object. Therefore, in order to preserve constant printing speed and printing direction at the image plane, an inverse kinematic approach is developed for HDSP to generate the trajectory and velocity of the platform for a given part geometry.Assuming the part geometry follows a curve C0(u), a virtual curve as shown in Supplementary Fig. 7b, where u is the curve parameter. This curve is defined in the printing coordinate system’s (PCS) origin and discretized by length Δs into N elements. Ck(u) represents a family of curves where kth element of this curve passes through the PCS’s origin and k = 0, …, N. In order to keep the printing direction constant (along the z-axis of PCS), the tangent of Ck(u) at z = 0, \({{{{\bf{V}}}}}_{0}^{k}\) should be coincident with the z-axis. \({{{{\bf{V}}}}}_{0}^{k}\) can be calculated as$${{{{\bf{V}}}}}_{0}^{k}={\frac{d{{{{\bf{C}}}}}^{k}(u)/du}{|d{{{{\bf{C}}}}}^{k}(u)/du|}\Bigg|}_{z=0}\cdot \varDelta s$$
(11)
However, \({{{{\bf{V}}}}}_{0}^{k}\) creates an angle θk with z-axis as shown in Supplementary Fig. 7b, c. In order to compensate for θk and adjust the orientation of Ck(u) at PCS origin, Ck(u) is transformed by the following transformation$${{{{\bf{C}}}}}^{k}(u)=\left[\begin{array}{ccc}1 & 0 & 0\\ 0 & \cos (-{\theta }^{k}) & -\,\sin (-{\theta }^{k})\\ 0 & \sin (-{\theta }^{k}) & \cos (-{\theta }^{k})\end{array}\right]\cdot {{{{\bf{C}}}}}^{k}(u)+{{{{\bf{V}}}}}_{0}^{k}$$
(12)
The trajectory of the end effector, Pk, can be found using Eq. 13 for each kth element as$${{{{\bf{P}}}}}^{k}={{{{{\bf{C}}}}}^{k}\left(u\right)\Big|}_{u=0},$$
(13)
and the vector defining the orientation of the end effector, \({{{{\bf{n}}}}}_{p}^{k}\), can be written as$${{{{\bf{n}}}}}_{p}^{k}={\frac{d{{{{\bf{C}}}}}^{k}(u)/du}{|d{{{{\bf{C}}}}}^{k}(u)/du|}\bigg|}_{u=0}.$$
(14)
In order to maintain a constant printing speed, Vpr, at the location of the image plane, the end effector transitional velocity, \({V}_{e}^{k}\), can be calculated as$${V}_{{{{\rm{e}}}}}^{k}=\Big|{{{{\bf{P}}}}}^{k}-{{{{\bf{P}}}}}^{k-1}\Big|/\left(\Big|{{{{\bf{V}}}}}_{0}^{k}\Big|/{V}_{{{{\rm{pr}}}}}\right).$$
(15)
Therefore, for any given C0(u), the trajectory, orientation and the transitional velocity of the end effector can be calculated using Eqs. 13–15, respectively. Supplementary Fig. 7b–f shows an example of a calculated trajectory at multiple steps during the printing of the object.
Holographic direct sound printing | Nature Communications
