The pressure field of an acoustic Bessel beam of zeroth order is expressed [11], as follows, p i = A J 0 e j k z z e − j ω t , J 0 = J 0 k R R , R = x 2 + y 2 (5) where J ₀ denotes the regular Bessel function of zeroth order, and A is the magnitude of the wave. According to the representation of a Bessel beam using Durnin rings [42, 43], J 0 k R R = 1 2 π ∫ 0 2 π e j k R x ⁡ cos ⁡ φ + k R y ⁡ cos ⁡ φ d φ , (6) this type of incident beam can be generated from the superposition of plane travelling waves coming towards a given axis at an incidence angle of β, as illustrated in Figure 1B. The Bessel beam is of the travelling type along the beam axis, i.e., z-direction, while fixed pressure and velocity nodes are (standing wave profile) along the R-direction. This special acoustic beam is axisymmetric, which implies the same acoustic pressure for points that are at the same radius from the axis and ∂ _θ p = 0. The pressure nodes are indicated in Figure 1A, and a zone in the vicinity of the beam axis is shaded as a region of interest for investigating the radiation force and torque exerted off-axis objects. The velocity field of this Bessel beam and its spatial derivatives are as in the following, ∇ p i ⋅ e z = j A k ⁡ cos ⁡ β J 0 e j k z z e − j ω t , ∇ p i ⋅ e R = − A k ⁡ sin ⁡ β J 1 e j k z z e − j ω t , v i ⋅ e z = A ρ f c f cos ⁡ β J 0 e j k z z e − j ω t , v i ⋅ e R = j A ρ f c f sin ⁡ β J 1 e j k z z e − j ω t , ∇ v i : e z e z = j A k ρ f c f cos 2 ⁡ β J 0 e j k z z e − j ω t , ∇ v i : e R e R = j A k ρ f c f sin 2 ⁡ β J 0 − J 1 k R R e j k z z e − j ω t , ∇ v i : e R e z = ∇ v i : e z e R = − A k 2 ρ f c f sin ⁡ 2 ⁡ β J 1 e j k z z e − j ω t , (7) where : operator denotes inner product of two second-order tensors, ∇ = ∂ R e R + 1 R ∂ θ e θ + ∂ z e z , and the dyadic products of the unit basis vectors are denoted by e _z e _R , e _z e _z , and e _R e _R . From Eq. 7, it can be seen that all the derivative fields are also axisymmetric, i.e., ∇v _i : e _R e _θ = ∇v _i : e _θ e _z = ∇v _i : e _θ e _θ = v _i ⋅e _θ = ∇ p _i ⋅e _θ = 0. These field expressions of the Bessel beam of the zeroth order in Eqs 5, 7 are required to formulate the dependence of the radiation force and torque on cylindrical coordinates (R, θ, z) and cone angle β.

Acoustic radiation force and torque are exerted on an object within an incident field as a result of radiation stresses ⟨ σ ⟩ that can mathematically be expressed as follows, 〈 σ 〉 = − 1 2 κ f 〈 p 2 〉 − 1 2 ρ f 〈 v 2 〉 I − ρ f 〈 v v 〉 , (8) where κ f = 1 / ρ f c f 2 is the mean fluid compressibility, ⟨⟩ operator denotes a time-averaging over one wave period [2]. The first term on the right hand side of Eq. 8 is called radiation pressure that is generated by the radiated momentum of the scattered pressure field, and the second term represents the Reynolds stresses indicating the contribution from the object’s surface oscillations to the radiation momentum [1, 2, 5, 8, 9, 44]. Acoustic radiation force and torque are obtained from radiation stresses ⟨ σ ⟩, as follows, F = ∫ Γ 〈 σ 〉 · n d Γ , T = ∫ Γ x × 〈 σ 〉 · n d Γ , (9) where Γ denotes a surface enclosing the object and n is the outward normal vector of Γ surface. By using the far-field approach and integrating the stresses on any fictitious surface enclosing the object [2, 8, 45], the radiation force and torque exerted on an object with arbitrary shape and size in the Rayleigh scattering limit ka ≪ 1 is expressed [5], as follows, F = F d + F c , F d = − 〈 α p p 2 ρ f ∇ p i 2 + j ω α v v v i ⋅ ∇ v i ⟩ r = 0 , F c = 〈 1 ρ f α p v ⋅ p i ∇ v i − v i ∇ p i ⟩ r = 0 . T = T d + T c , T d = j ω α v v v i × v i r = 0 , T c = j ω p i α v p × v i r = 0 , (10) where d and c subscripts indicate the force and torque associated with direct polarization and Willis coupling, respectively. By substituting Eqs 5, 7 into Eq. 10 and using the reciprocity properties of the polarizability tensor [5, 22, 35, 38, 46], the analytical expressions of the force and torque are obtained, as in the following, F d ⋅ e z = − 2 E i k ρ f 1 κ f Im α p p cos ⁡ β J 0 2 + k c f Re α v v z z cos 3 ⁡ β J 0 2 + k c f Re α v v R R sin 2 ⁡ β ⁡ cos ⁡ β J 1 2 , F d ⋅ e R = − 2 E i k ρ f − 1 κ f Re α p p sin ⁡ β J 0 J 1 + k c f Re α v v R z sin 2 ⁡ β ⁡ cos ⁡ β J 0 2 + J 1 2 − J 0 J 1 k R R + k c f Im α v v z z cos 2 ⁡ β ⁡ sin ⁡ β J 0 J 1 − k c f Im α v v R R sin 3 ⁡ β J 0 J 1 − J 1 2 k R R , (11) F c ⋅ e z = − 2 E i k ρ f 2 c f Re α p v R sin ⁡ β ⁡ cos ⁡ β J 0 J 1 , F c ⋅ e R = − 2 E i k ρ f c f Im α p v R sin 2 ⁡ β J 0 2 + J 1 2 − J 0 J 1 k R R , (12) T d ⋅ e z = − 2 E i k ρ f − c f Re α v v θ z sin ⁡ β ⁡ cos ⁡ β J 0 J 1 + c f Im α v v θ R sin 2 ⁡ β J 1 2 , T d ⋅ e θ = − 2 E i k ρ f c f Re α v v z z sin ⁡ β ⁡ cos ⁡ β J 0 J 1 − c f Im α v v z R sin 2 ⁡ β J 1 2 + c f Im α v v R z cos 2 ⁡ β J 0 2 + c f Re α v v R R sin ⁡ β ⁡ cos ⁡ β J 0 J 1 , T d ⋅ e R = − 2 E i k ρ f − c f Im α v v θ z cos 2 ⁡ β J 0 2 − c f Re α v v θ R sin ⁡ β ⁡ cos ⁡ β J 0 J 1 , (13) T c ⋅ e z = − 2 E i k ρ f − 1 κ f Re α v p θ sin ⁡ β J 0 J 1 , T c ⋅ e θ = − 2 E i k ρ f 1 κ f Re α v p z sin ⁡ β J 0 J 1 + 1 κ f Im α v p R cos ⁡ β J 0 2 , T c ⋅ e R = − 2 E i k ρ f 1 κ f Im α v p θ cos ⁡ β J 0 2 , (14) where E i = A 2 / 4 ρ f c f 2 denotes the energy density of the incident Bessel beam. Both partial forces F _d and F _c have components only in z- and R-directions, which implies objects being pushed or pulled in axial and radial directions in an axisymmetric acoustic field. For the partial torques T _d and T _c , in addition to the component in the θ-direction, there exist two more components the z- and R-plane, which means a full 3D torque is induced by an axisymmetric acoustic field. This is in contrast to the case of a standing, plane wave, where only in-plane torque components are generated [5, 22]. Compared to a plane wave, this torque property of the Bessel beam shows that a 3D rotational manipulation an object with arbitrary shape can be realized even using 2D propagation type, axisymmetric beams. The real and imaginary parts of all polarizability coefficients, from the cylindrical description of the polarizability tensor in Eq. 4, are required to obtain the three components of the radiation torque; while only the z- and R-dependent coefficients are required to calculate both partial forces. This indicates the significance of accounting for the shape of the object to find the radiation force and torque in a Bessel beam of zeroth order.

This case corresponds to zero radial variation of the wavefront k _R = 0, which means the Bessel beam degenerates to a travelling plane wave. By substituting sin β = 0, cos β = 1, J ₀ = 1, J ₁ = 0, and J ₁/k _R R ≈ 1/2, the partial radiation forces and torques in Eqs 11–14 become F d ⋅ e z = − 2 E i k ρ f 1 κ f Im α p p + k c f Re α v v z z , F d ⋅ e R = 0 , (15) F c ⋅ e z = 0 , F c ⋅ e R = 0 , (16) T d ⋅ e z = 0 T d ⋅ e θ = − 2 E i k ρ f c f Im α v v R z , T d ⋅ e R = − 2 E i k ρ f − c f Im α v v θ z , (17) T c ⋅ e z = 0 , T c ⋅ e θ = − 2 E i k ρ f 1 κ f Im α v p R T c ⋅ e R = − 2 E i k ρ f 1 κ f Im α v p θ . (18) Equations 15–18 are identical to expressions in Eq. [26] of Ref. [5], verifying our formulation of the radiation force and torque under a Bessel beam at the limit case of β = 0. It is also found that Willis coupling effects show zero contribution to the force and only generates torque in the in-plane directions of a travelling plane wave. The forces and torques generated by a travelling plane wave are constant and independent of the location of the objects, meaning that the radiation force and torque fields are uniform and only depend on the polarizability response of the object.

This case corresponds to the long wavelength limit in the radial R direction, and k _R = k sin β ≈ kβ. This case is of interest since it represents relatively weaker propagation of acoustic energy in the R direction, i.e., the standing wave profile and fixed pressure and velocity nodes. By using the asymptotic approximation of the Bessel function in this limit, i.e., J ₀ ≈ 1 and J ₁ ≈ k _R R/2, sin β ≈ β, cos β ≈ 1 and β ² ≪ β, Eqs 11–14 can be approximated, as follows, F d ⋅ e z = − 2 E i k ρ f 1 κ f Im α p p + k c f Re α v v z z , F d ⋅ e R = − 2 E i k ρ f − k 2 κ f Re α p p β 2 R + k 2 c f 2 Im α v v z z β 2 R + k c f 2 Re α v v R z β 2 , (19) F c ⋅ e z = − 2 E i k ρ f k c f Re α p v R β 2 R , F c ⋅ e R = − 2 E i k ρ f c f 2 Im α p v R β 2 , (20) T d ⋅ e z = − 2 E i k ρ f − k c f 2 Re α v v θ z β 2 R , T d ⋅ e θ = − 2 E i k ρ f k c f 2 Re α v v z z + α v v R R β 2 R + c f Im α v v R z , T d ⋅ e R = − 2 E i k ρ f − c f Im α v v θ z − k c f 2 Re α v v θ R β 2 R , (21) T c ⋅ e z = − 2 E i k ρ f − k 2 κ f Re α v p θ β 2 R , T c ⋅ e θ = − 2 E i k ρ f k 2 κ f Re α v p z β 2 R + Im α v p R , T c ⋅ e R = − 2 E i k ρ f 1 κ f Im α v p θ . (22) For a given cone angle β, we can see from Eqs 19–22 that the radiation force and torque spatially depend only on their radial coordinate R, changing linearly, while these fields are independent of the axial coordinate z, due to the travelling wave nature of propagation in this direction. Compared to the limit case of β = 0, it was found that the additional terms proportional to β ² appear in the force and torque expressions T _c ⋅e _R , as result of the weak propagation in the R-direction, except for the radial component of the Willis coupling torque T _c ⋅e _R that remains unchanged. These analytical expressions reveal the sensitivity of the radiation force and torque, both direct and Willis coupling parts, to the change in the wave front that is the weak Bessel-type undulations. This provides β angle as another degree of freedom for manipulation of objects with arbitrary shape, using this type of radially weak Bessel beams.

There is a velocity node, corresponding to maximum pressure, on the axis of the zeroth-order Bessel beam. This is of interest since it is a relatively good example of focused acoustic beams that can be produced using a meta-lens [47] or a phased transducer array [48]. When the centroid of an object is located on the beam axis, the expressions of acoustic radiation force and torque in Eqs 11–14 are further simplified, as follows, F d ⋅ e z = − 2 E i k ρ f 1 κ f Im α p p cos ⁡ β + k c f Re α v v z z cos 3 ⁡ β , F d ⋅ e R = − 2 E i k ρ f k c f 2 Re α v v R z sin 2 ⁡ β ⁡ cos ⁡ β , (23) F c ⋅ e z = 0 F c ⋅ e R = − 2 E i k ρ f c f 2 Im α p v R sin 2 ⁡ β , (24) T d ⋅ e z = 0 T d ⋅ e θ = − 2 E i k ρ f c f Im α v v R z cos 2 ⁡ β , T d ⋅ e R = − 2 E i k ρ f − c f Im α v v θ z cos 2 ⁡ β , (25) T c ⋅ e z = 0 T c ⋅ e θ = − 2 E i k ρ f 1 κ f Im α v p R cos ⁡ β , T c ⋅ e R = − 2 E i k ρ f 1 κ f Im α v p θ cos ⁡ β . (26) From Eq. 23, we found that an object on the beam axis is subjected to only an axial force, i.e., in the z-direction as a result of the direct polarization. ² Since the range of values of the cone angle is 0 ≤ β < π/2, cos β > 0 and the direction of the axial force is determined from the polarization coefficients, i.e., Im(α _pp ) and Re ( α v v z z ) , implying that some objects can experience pull-in effects from a Bessel beam of zeroth order at the ka ≪ 1 limit. The axial components of the Willis-coupling partial force F _c and both direct and Willis coupling torques T _d and T _c are zero. Equations 25, 26 show that the dependence of the direct and Willis coupling-induced torques on cone angle β are cos² β and cos β, respectively, implying that partial torque fields T _d and T _c are independent of each other and the shape asymmetry needs to be accounted for separately. Similarly, the partial force fields F _d and F _c are independent, as can be seen from Eqs 23, 24. Finally, the most significant contribution of the Willis coupling is the radial component of the radiation force, coming from F _c ⋅e _R , meaning that objects with shape asymmetry are pushed away from the beam axis in the direction that Im ( α p v R ) is non-zero and maximum, if the radiation torque is accounted for. This depends on the orientation of the objects and its shape asymmetry with respect to the beam axis. For axisymmetric shapes with one degree of asymmetry along the axis of revolution, e.g., shapes of a Helmholtz resonator, unhealthy red blood cells [49], or a meta-atom with controlled Willis coupling [22], this radial force is zero if the axis of symmetry of the object is aligned with the beam axis. It is inferred that, in general, objects tends to be off-axis in such acoustic beams; hence, the axisymmetric wavefront provides the same level of control compared to a 2D arbitrary wavefront.

We consider this special off-axis case such that k _R R ≪ 1, as shown by the yellow-shaded strip in Figure 1A. This condition applies to long wavelength limit ka ≪ 1, and we assume the radial distance is close to zero (vicinity of beam axis), while 0 ≤ β < π/2. By using the asymptotic approximation of Bessel functions J ₀ and J ₁, the radiation force and torque in Eqs 11–14 are approximated, as follows, F d ⋅ e z = − 2 E i k ρ f 1 κ f Im α p p cos ⁡ β + k c f Re α v v z z cos 3 ⁡ β + k 3 c f 4 Re α v v R R sin 4 ⁡ β ⁡ cos ⁡ β R 2 , F d ⋅ e R = − 2 E i k ρ f − k 2 κ f Re α p p sin 2 ⁡ β R + k 2 c f 2 Im α v v z z cos 2 ⁡ β ⁡ sin 2 ⁡ β R + k c f 2 Re α v v R z sin 2 ⁡ β ⁡ cos ⁡ β 1 + k 2 2 sin 2 ⁡ β R 2 − k c f 4 Im α v v R R sin 4 ⁡ β R , (27) F c ⋅ e z = − 2 E i k ρ f k c f Re α p v R sin 2 ⁡ β ⁡ cos ⁡ β R , F c ⋅ e R = − 2 E i k ρ f c f Im α p v R sin 2 ⁡ β 1 2 + k 2 ⁡ sin 2 ⁡ β 4 R 2 , T d ⋅ e z = − 2 E i k ρ f − k c f 2 Re α v v θ z sin 2 ⁡ β ⁡ cos ⁡ β R + k 2 c f 4 Im α v v θ R sin 4 ⁡ β R 2 , (28) T d ⋅ e θ = − 2 E i k ρ f k c f 2 Re α v v z z sin 2 ⁡ β ⁡ cos ⁡ β R − k 2 c f 4 Im α v v z R sin 4 ⁡ β R 2 + c f Im α v v R z cos 2 ⁡ β + k c f 2 Re α v v R R sin 2 ⁡ β ⁡ cos ⁡ β R , T d ⋅ e R = − 2 E i k ρ f − c f Im α v v θ z cos 2 ⁡ β − k c f 2 Re α v v θ R sin 2 ⁡ β ⁡ cos ⁡ β R , (29) T c ⋅ e z = − 2 E i k ρ f − k 2 κ f Re α v p θ sin 2 ⁡ β R , T c ⋅ e θ = − 2 E i k ρ f k 2 κ f Re α v p z sin 2 ⁡ β R + Im α v p R cos ⁡ β , T c ⋅ e R = − 2 E i k ρ f 1 κ f Im α v p θ cos ⁡ β . (30) From Eq. 30, it was found that the T _c ⋅e _R is constant and independent of the R-coordinate. Other components of the partial radiation force and torques change with the radial coordinate R, e.g., proportional to R and R ². There are also some constant terms in F _d ⋅e _z , F _c ⋅e _R , T _d ⋅e _θ , T _d ⋅e _R , T _c ⋅e _θ , and T _c ⋅e _R , which provide the base value of these force and torque fields in the k _R R ≪ 1 region irrespective of the object’s radial location.

To investigate the influence of non-planar wavefront of the zeroth order Bessel beam on the acoustic radiation force and torque, we consider an object with asymmetric shape such that major asymmetry occurs along the z-direction for θ = 0 and the polarizability coefficients become α p p Ω s = 7.2 × 1 0 − 3 + 0.0 j , α v p = − j ω ρ f α p v α p v α p p = 8.4 × 1 0 − 4 + 7.2 × 1 0 − 4 j 3.4 × 1 0 − 10 − 4.8 × 1 0 − 5 j 0.0 , − , 7.4 × 1 0 − 2 j , α v v α p p = 1.0 , × , 1 0 − 4 − 2.39 j 3.6 × 1 0 − 9 − 2.89 , × , 1 0 − 5 j 0.0 , + , 1.3 × 1 0 − 3 j 1.0 , × , 1 0 − 4 − 2.39 j 5.0 , × , 1 0 − 8 − 2.4 × 1 0 − 4 j Sym. 1.0 , × , 1 0 − 4 − 2.43 j . (31) where Ω _s denote the volume of the object, and these values satisfy the acoustic energy balance, i.e., they are less than the maximum admissible polarizability [38]. Although these are just examples of the values of polarizability coefficients, they correspond to a sphere of size a with a blind circular hole which leads to the asymmetry along the hole axis [5], i.e., z-direction in this study. The non-zero coefficients with relatively smaller values correspond to the numerical discretization of this object for the calculation of monopole and dipole moments using Boundary Element Method, as outlined in Ref. 5. Nevertheless, these small values indicate the level of asymmetry in x- and y-directions, compared to the intended one in the z-direction.

The ratio of partial forces and torques in cases II to IV against the force and torque of a plane travelling wave, i.e., case I, is considered for comparison. The magnitude of the force and torque for case I are expressed, as follows, Q = F d + F c ⋅ F d + F c , Z = T d + T c ⋅ T d + T c . (32) First, we investigate the case of weak Bessel beam, i.e., case II and Eqs 19–22, by choosing β = 5°, corresponding to k _R /k ≈ 0.1. The three components of the partial forces and torques are presented in Figure 2, for relatively large range of 0 < k _R R < π.

The direct partial force F _d shows a much larger radial component than the axial z-component, which is almost the same as the axial force component of a plane travelling wave, as shown in Figure 2A. Similar differences between components are observed in the Willis-coupling partial force F _c , as shown in Figure 2B; however, the contribution of F _c to the total force is negligible compared to F _d . It is noted that the only non-zero force component for the reference case of a plane travelling wave is the axial component of the direct partial force F _d ⋅e _z . In all cases I to IV, the force component in the θ-direction is always zero due to the axisymmetric wavefront of the Bessel beam. Comparing the partial torques in Figures 2C,D, we observed that the component in the θ-direction is the largest and the Willis coupling contribution is much larger than the direct partial torque. The R-component of the partial torques is of the same of order of magnitude, but the z- and θ-components of Willis coupling torque T _c are much larger than those of T _d . The force and torque components changes almost linearly with respect to the radial R-coordinate, except for R → 0 corresponding to region in the vicinity of the beam axis, i.e., R = 0. These results indicate that the effects of asymmetry of the object’s shape manifest as a relatively large torque component in the θ-direction, for Bessel beam of small β angle representing a weak non-planar travelling wave.

Next, the case of objects at the vicinity of the beam axis, i.e., case IV, is shown in Figure 3. The same object with polarizability tensor given in Eq. 31 is placed one radius a away from the beam axis, i.e., R/a = 1 corresponding to kR ≈ 0.03. The force and torque results are shown for cone angle β in the range of 0–90°. As shown in Figure 3A, the direct partial force F _d has a larger component in the radial R-direction. The component in the z-direction undergoes a sign reversal at β ≈ 78°, from opposite (pull) to same (push) as the wave propagation direction along the z-axis. The Willis coupling partial force F _c , shown in Figure 3B, is smaller than the F _d by several orders of magnitude and can be neglected. Comparing the partial torques T _d and T _c , shown in Figures 3C,D, respectively, it is observed that the asymmetry in shape results in relatively larger torque components in the θ-direction, and the magnitude increases as the cone angle β approaches π/2, corresponding to the limit of zero axial propagation. Finally, the significant difference between the magnitudes of the force and torque components correspond to the values of the polarizability coefficients for this study, which comes from a discretized geometry of a sphere with a blind circular hole. Nevertheless, they can be considered as an indication of the force and torque acting upon an object with asymmetric shape from several directions. Finally, our results show that the cone angle β can be used as a design parameter for axial manipulation using zeroth order Bessel beam.