Driving techniques for phased array systems.

Abstract

A method for driving a phased array system is provided. [Solution] A method for generating acoustic fields designed to produce mid-air haptics, parametric audio, acoustic levitation, and acoustic imaging, which mitigates the effects of air variations, produces trap points in the air, defines phased array optimization with respect to vectors for producing more consistent haptic effects, and defines one or more control points or regions in space through which the controlled acoustic field is passed. The method uses reduced representations to construct the acoustic field basis functions, performs efficient evaluation of complex-valued functions for massive throughput, generates Krylov subregions of the matrix, and controls the acoustic field with control points and / or regions different from those used to create the acoustic field basis functions.

Description

[Technical field]

[0001] Related Applications This application claims the benefit of the following eight U.S. provisional patent applications, all of which are incorporated herein by reference in their entireties:

[0002] 1. Application No. 62 / 433,544, Filing Date: December 13, 2016

[0003] 2. Application No. 62 / 436,457, Filing date December 20, 2016

[0004] 3. Application No. 62 / 437,791, Filing date December 22, 2016

[0005] 4. Application No. 62 / 438,512, filed December 23, 2016

[0006] 5. Application No. 62 / 594,574, Filing date: December 5, 2017

[0007] 6. Application No. 62 / 594,687, Filing date December 5, 2017

[0008] 7. Application No. 62 / 594,836, Filing date December 5, 2017

[0009] 8. Application No. 62 / 594,874, Filing date: December 5, 2017

[0010] Disclosure Fields This disclosure relates generally to improved techniques for driving phased array systems, primarily for generating acoustic fields designed to produce aerial haptics, parametric audio, acoustic levitation, and acoustic imaging. [Background technology]

[0011] background Continuous distributions of sound energy, which we call "acoustic fields," can be used for a variety of applications, including air-based haptic feedback.

[0012] In this acoustic field, one or more control points can be defined. These control points can be amplitude modulated with a signal, thereby generating vibrotactile feedback in the air. An alternative way to create feedback is to create control points that are not amplitude modulated and move them in space to generate a spatiotemporal modulation that can be felt.

[0013] The minimum size and maximum power of a control point created by the device depends on the array layout, shape and relative position of the control point or points generated. This may make it impossible to match the behavior of a control point in one location with its behavior in another. Further complications may occur when control points are dynamically moved through space. The resulting change in the haptic function of the points experienced by the user is undesirable and often unintuitive, resulting in confusion and loss of experiential immersion. A system that can mitigate the effect of air change and reverse it to provide a consistent haptic experience would be commercially valuable.

[0014] Furthermore, it has been shown that a general nonlinear optimization system can be used to generate trapping points in air. However, the mechanism used for this purpose has drawbacks including not evaluating potential fields suitable for levitation on phased array systems, using an unnecessarily expensive and time-consuming optimization process, and being limited to phase manipulation only.

[0015] Addressing these issues is important because the traps created by existing systems are not as strong as they could be and cannot be scaled accurately due to these shortcomings. This makes it difficult to use alternative systems that utilize acoustic levitation to provide additional visual and tactile systems to facilitate a multimodal approach.

[0016] Furthermore, an acoustic field can be controlled in space by defining one or more control points. Each point can be assigned a value equal to the desired amplitude at the control point. A physical set of transducers can then be controlled to create a sound field that exhibits the desired amplitude at the control points.

[0017] Due to the conditions of existing systems, the natural directional component of acoustic phased arrays has received little attention. For example, radiofrequency therapeutic medical systems, where the objective is to focus energy in a small area, are well approximated by a scalar pressure (the change in the number of air molecules involved) to determine the optimal solution for injecting energy, and therefore do not require a vector solution. Focusing on the behavior of the phased array vectors (which represent the direction of motion of the air molecules) can be ignored in most cases, so achieving this efficiently is not yet a priority.

[0018] However, to create haptic feedback in mid-air, it is necessary to exchange vibrations with human skin. This depends entirely on the interface with the wavefront and the momentum carried by the wave, which is a consequence of its relationship with the vector behavior of the wave. Not only this, but as systems are focused such that the angles of incidence between the transducer elements are increasingly different, the assumed linear relationship between momentum and pressure breaks down. It is impossible to ignore the vector components in this system, and scalar pressure is insufficient to consistently find a suitable solution to create mid-air haptic effects. For this reason, making devices aware of this momentum effect is commercially important.

[0019] Additionally, an acoustic field can be controlled by defining one or more control points in space. Each point can be assigned a value equal to the desired amplitude at the control point. A physical set of transducers can then be controlled to produce an acoustic field exhibiting the desired amplitude at the control point.

[0020] It is shown that the calculation of the eigenvectors corresponds to the solution of some very specific purely quadratic constrained problems in matrix algebra. Considering a particular problem in aerial haptics, if the eigenvectors can form a matrix that solves for the control coefficients of the transducer, the solution can form the basis of a device implementation that represents a simplification with respect to existing technology.

[0021] Furthermore, in order to use reduced representations to describe the basis functions, one must first develop reduced representations that allow simple (low entropy) descriptions of functions that are useful in terms of the details of the acoustic field.

[0022] The present disclosure will depend on what is theoretically accessible from the user's point of view.

[0023] Once such an accessible low-entropy representation is found, there will also be a way to convert it into its equivalent transducer basis (a transducer basis set describing the basis of one transducer acoustic field of a potentially large number of basis functions in a later solution). This will allow potential users to design a variety of different structures within the acoustic field, building bespoke systems for maximal manipulation of the energy flow in the field leading to the formation of focal points, Bessel beams, bottle beams, bending beams, etc., many different effects at once. These would extend the concept of control points to control regions encapsulating richer acoustic field control methods.

[0024] It has also been discovered that when considering momentum transfer in touch, there are several problems that must be addressed in order to use particle velocity as a proxy for wave direction.

[0025] Additionally, performing efficient parallel complex-valued multiplications with the same second operand was necessary for the massive throughput required to achieve a device with performance metrics consistent with those outlined in the requirements document. In particular, computing one element in parallel of the sum that led to the computation of a row of the required optimization matrix was sufficient to meet the requirements, but only if the subsystem was replicated a sufficient number of times. This pressing need for efficiency led to a different look at old algorithms, which gave rise to new approaches for evaluating complex transcendental functions on top of simple multiplication.

[0026] In this wider world, the resulting algorithms are applicable wherever an element of processing hardware must interact efficiently or in real time with complex-valued data with non-negligible coefficients. This includes, but is not limited to, communications including GPS, wireless, radio and satellite transmission and reception systems, audio processing, phased array systems including airborne haptics, time-of-flight imaging and digital interferometry, and computation of scattering effects. This could lead to a wider range of opportunities for the present invention.

[0027] Volder's algorithm (Volder, J. - The CORDIC computational method - 1959), also known as CORDIC, can multiply unit complex values, but cannot include coefficients in the second operand of the calculation without going beyond the calculation of the basic CORDIC algorithm.

[0028] Another algorithm invented in 1994, implemented by the author for this reason several years ago, known as BKM after its authors Bajard, Kla and Muller, is even more capable of achieving the goal of complete complex-complex multiplication (Bajard, J, Kla, S. and Muller, J. - BKM: A New Hardware Algorithm for Complex Elementary Functions - 1994). However, this algorithm has a complex implementation and therefore has not received as much attention as it would otherwise. This algorithm is able to perform complete complex-complex multiplication by reducing the complex-valued operands to their complex logarithms, adding them and then exponentiating them.

[0029] Using BKMs and possibly new algorithms, there are more efficient ways to achieve this complex-valued multiplication and complex-valued division, both serial and parallel.

[0030] A 1999 addendum to the BKM algorithm addressed the overly complex implementation of complex logarithms and simplified the implementation somewhat. (Bajard, J. and Imbert, L. - Evaluation of Complex Initial Functions: A New Version of BKM - 1999). This BKM approach from the 1994 and 1999 addendums forms the algorithm referred to in this document as simply BKM.

[0031] While trying to optimize the algorithm for incorporation into actual hardware, the idea of ​​a new technique emerged that prioritizes ease of integration into our, or more generally existing, systems over evaluating mathematical functions as described. This new technique, when developed into proof-of-concept code, turns out to have many novel properties that make it more promising for general adoption of hardware implementations than the original algorithm. In addition, if the exact numerical values ​​demanded for the original technique are required, further scaling may be applied to search for these results, while the algorithm as a whole remains competitive. In either case, the result is a technique with substantial advantages over the state of the art.

[0032] It is known that quadratic problems with unit norm can be solved by constructing a matrix whose eigenproblems are equivalent to the quadratic problem. A quadratic problem matrix greater than rank 1 can be generated by summing many rank 1 problem matrices and / or expanding the null space of an existing problem matrix. By generating a Krylov subspace of this matrix, an optimized subset of the implicit basis functions can be generated. By using the Krylov subspace as a basis set to linear systems, a correct or approximate evaluation of the basis set can be achieved using a large number of test functions.

[0033] Because how many eigenvectors are used and from where they are selected, the power of each subsequent eigenvector decreases, allowing for graceful degradation of a complex set of many basis functions.

[0034] Furthermore, the complex-valued linear system for defining the resulting acoustic field shape is JPEG0007599472000001.jpg710 basis functions are taken and recombined into a resulting field that describes the desired configuration of sound pressure at the carrier frequency.

[0035] Up to this point, we have used control point basis functions to operate the control points. Described here is a method that generalizes to apply a different set of control test functions to the initial control point or region basis functions.

[0036] This is important because it breaks the connection between the behavior of the initial "palette" in the acoustic field and the objective function used by linear systems: for example, it can be used to create a basis set of Bessel functions, which can then be used to actuate the control points. [Brief description of the drawings]

[0037] BRIEF DESCRIPTION OF THE DRAWINGS The accompanying drawings, together with the following detailed description, which refer to identical or functionally similar elements through separate views, are incorporated in and form a part of the specification, and serve to further illustrate embodiments of the concepts that comprise the claimed invention(s), and explain various principles and advantages of these embodiments.

[0038] [Figure 1] FIG. 1 shows a haptic waveform of arbitrary dynamic range that can be processed through a nonlinear remapping function to generate a new waveform of controlled dynamic range. [Diagram 2] FIG. 2 is a haptic waveform of arbitrary dynamic range that can be processed through a nonlinear remapping function to generate a new waveform of controlled dynamic range.

[0039] [Diagram 3] FIG. 3 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 4] FIG. 4 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Diagram 5] FIG. 5 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 6] FIG. 6 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 7] FIG. 7 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 8] FIG. 8 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 9] FIG. 9 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 10] FIG. 10 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 11] FIG. 11 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques. [Figure 12] FIG. 12 shows the effective amplitude modulation waveforms resulting from different configurations of control points applying different range control techniques.

[0040] [Figure 13] FIG. 13 is a projection of the sound pressure resulting from the optimization of the Gorkov related function. [Figure 14] FIG. 14 shows the projection of the sound pressure resulting from the optimization of the Gorkov related function. [Figure 15] FIG. 15 is a projection of the sound pressure resulting from the optimization of the Gorkov related function. [Figure 16] FIG. 16 shows the projection of the sound pressure resulting from the optimization of the Gorkov related function.

[0041] [Figure 17] FIG. 17 is a graph of the same wave with the same sound pressure level, but in a scenario where the effect of the vector components would cause a maximum haptic effect, but no haptic effect.

[0042] [Figure 18] FIG. 18 is a simulation of measurements of the sound field produced by the phased array in dB SPL, dB SIL and dB SVL.

[0043] [Figure 19] FIG. 19 shows a control point plot produced by solving a quadratic problem.

[0044] [Figure 20]FIG. 20 shows a cut through the sound field at the edge of the near field region parallel to the flat linear transducer array.

[0045] [Figure 21] FIG. 21 shows a graph of the Bessel function.

[0046] [Figure 22] FIG. 22 shows a cross-sectional view of an approximation to a higher order Bessel beam and a Mathieu beam. [Figure 23] FIG. 23 shows cross-sectional views of approximations to higher order Bessel beams and Mathieu beams.

[0047] [Figure 24] FIG. 24 shows a bent beam, a bottle beam, and a partial cross section of a bottle beam. [Diagram 25] FIG. 25 shows a bent beam, a bottle beam, and a partial cross section of a bottle beam. [Figure 26] FIG. 26 shows a bent beam, a bottle beam, and a partial cross section of a bottle beam.

[0048] [Figure 27] FIG. 27 shows an elongated wavelength beam, a snake beam, and a helical beam. [Figure 28] FIG. 28 shows an elongated wavelength beam, a snake beam, and a helical beam. [Figure 29] FIG. 29 shows an elongated wavelength beam, a snake beam, and a helical beam.

[0049] [Diagram 30] FIG. 30 shows the effect of removing transducers that contribute to the wave from the opposite direction to the desired wave. [Diagram 31] 31 shows the effect of eliminating transducers that contribute to the wave from the direction opposite to the desired wave.

[0050] [Diagram 32]FIG. 32 shows the convergence regions of the logarithm method and the power method without remapping. [Diagram 33] FIG. 33 shows the convergence regions of the logarithm method and the power method without remapping. [Diagram 34] FIG. 34 shows the convergence regions of the logarithm method and the power method without remapping.

[0051] [Diagram 35] FIG. 35 shows the convergence graphs for the logarithmic method and the power method.

[0052] Those skilled in the art will appreciate that elements in the figures are illustrated for simplicity and clarity and have not necessarily been drawn to scale - for example, the dimensions of some of the elements in the figures may be exaggerated relative to other elements to help to improve understanding of embodiments of the present invention.

[0053] Components of the apparatus and methods are represented, where appropriate, by conventional symbols in the drawings, showing only the specific details relevant to understanding the embodiments of the invention so as not to obscure the disclosure, as would be apparent to one of ordinary skill in the art having the benefit of the description herein. DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

[0054] Detailed Description The techniques, devices and methods described above may be used in combination or separately in an air-based haptic system.

[0055] 1. Dynamic Range Reduction for Airborne Haptic Devices

[0056] I. Psychophysical Modeling

[0057] The behavior of the system and the haptic output must be tailored to the way humans perceive touch. For example, as the size of a control point increases, the physical amplitude it creates does not change, but the haptic perception increases. Cross-modal effects, such as those from audio, can help improve haptic feedback if they are correlated. This leads to discrepancies between physically defined control points, necessitating the use of a psychophysical modeling step to achieve parity of perception.

[0058] II. Remapping using a psychohaptic model

[0059] Due to the physical limitations of the playback device, the exact waveform may not be reproducible as specified. The waveform may be constructed within a psychologically defined perceptual space. Alternatively, it may consist of measurements of strength and force on a substantially absolute physical scale, which may be a haptic texture, or a record of any constructed waveform potentially represented with further metadata. This waveform represents the motion of the control point given a "perfect" airborne haptic device. This is then perceptually evaluated using a psychophysical model that can describe and quantify the perceptual aspects of the reproduction of the control point. While reproducing this haptic effect, the control point may be in motion. The motion is included as part of the perceptual evaluation of the haptic effect. This can then be reinterpreted into the space of outputs available to an instance of the airborne haptic system, constrained by the user position, relative geometry, and available power limitations, and the system is configured to create a substantially equivalent control point.

[0060] The psycho-haptic model gives the device the ability to determine how the modulation waveform of the control points can be modified to induce substantially the same haptic response as the complete device. The process of achieving this remapping can also include a psycho-acoustic modeling step to prevent the resulting haptic waveform from being noisy, or to select a more audibly acceptable modulation that is close in terms of the perceived haptic space but less noisy. It is known that many waveforms can be made more haptically powerful by increasing the amount of audible noise. In some situations, this may be loud enough to be intolerable, while in others this may be desirable. The psycho-acoustic model can then be used in conjunction to bias the output to be quieter or more powerful, with the audio being placed in a haptic trade-off setting.

[0061] III. Mental Tactile Unit

[0062] When defining waveforms in a psychological space that are intended to be felt as equivalent in amplitude, units must be designed that take into account the psychological effects of control points generated in the air. The most notable of these is the integration of haptics into the hand. The focal area contained in each control point again varies in size and shape depending on the properties of the physical device, and this must be countered in the configuration of the psychotactile unit. Such a psychotactile unit would modify the physical control point amplitude it evaluates depending on the size and shape of the focal area, countering the effects of integration of touch sensations over a large area. This approach allows the same number of psychotactile units to code equivalent sensations, preserving haptic equivalence and physical instantiation between devices.

[0063] IV. Dynamic Waveform Mapping

[0064] To arrive at the best mapping to the target device, the complete waveform and its trajectory may be considered. However, because interactivity and therefore low latency are prerequisites, only parts of the waveform or individual samples are available. The entire waveform can be remapped to the capabilities of the target device, but to do this the waveform needs to be known in advance, so it can be seen as a pre-processing step of the known waveform. An alternative method is to remap the physical properties of the amplitude waveform samples with amplitude, since remapping of amplitude samples can be achieved in real time with interactive waveforms and is more suitable for providing haptic feedback in real time.

[0065] A variety of functions are considered good candidates to provide the remapping for the control point amplitude samples. These include sigmoidal functions along with log-linear curves, which can convincingly reproduce the perceptual sensory input on devices with low dynamic ranges in other modalities. While these curves generally do not require the specification of an upper limit in the mapping function, as shown in Figures 1 and 2, given that there is a finite amount of power available, the dynamic range can be used more effectively. Figure 1 shows a haptic waveform of any dynamic range 100, which can be processed through a nonlinear remapping function 110 to generate a new waveform of a controlled dynamic range 120. This can be used to take a haptic waveform of any dynamic range and convert it to a smaller range. Since the amplitude of the waveform can be infinite, this does not fill the entire dynamic range of the output, but rather results in unused output range. In Figure 2, similar to Figure 1, a nonlinear function 210 remaps the amplitude waveform 200 to the amplitude available to the target device 220. Knowing the spread of the wave allows the remapping function to be precisely tailored to the device-specific output range. This allows for the development of a full range of output devices.

[0066] For illustrative purposes, the following examples use minimalist linear scaling to demonstrate how mapping in air-based haptic devices avoids undesirable behavior.

[0067] The graph 300 shown in Figure 3 indicates that the first step in achieving a remapping of a waveform to a reproducible perception is to feed the data into a solver function that describes the waveform and its limits 300. In Figure 3, the solid line represents JPEG0007599472000002.jpg1112, and the dotted line is JPEG0007599472000003.jpg1119. It is a simple waveform input remapping function, where the y-axis is in device internal units.

[0068] This is achieved by defining the maximum amplitude that the waveform must be played to radiate without clipping, which at each sample time is: This is expressed as JPEG0007599472000004.jpg1219. This is the amplitude of the input waveform. This is fed into the solver as input along with JPEG0007599472000005.jpg1111.

[0069] Shown in Figure 4 is a graph 400 that illustrates the results of attempting to recreate this waveform without applying any modifications to the input data. The result is that the waveform is clipped due to insufficient power available at the control points to recreate the required amplitude. In Figure 4, the straight line JPEG0007599472000006.jpg1227, which is the amplitude This is the result of JPEG0007599472000007.jpg1226, which is the result of trying to directly output the user-specified amplitude in Figure 1, regardless of the dynamic range of the device.

[0070] The physical configuration of the active transducers in the system is JPEG0007599472000008.jpg88 position The desired control point vibration in JPEG0007599472000009.jpg78 is JPEG0007599472000010.jpg78, the function vibration is used to solve the function for the air haptic system. JPEG0007599472000011.jpg1229 can be defined and, if a given device is physically able to provide this at a given time, it will result in either the desired carrier amplitude at the control point, or the closest amplitude that can be reproduced at the location of the control point.

[0071] Evaluate JPEG0007599472000012.jpg1170 to find the maximum range to which the input amplitude can be mapped at this point. This does not need to be exact, since the control points do not move much between samples, just the iterations of the control points from the previous iteration. The value JPEG0007599472000013.jpg822 can be used. This value can also be filtered to remove rounding errors and other sources of noise. In this case, the remapping function takes a single sample of the waveform and divides it into the interval Section from JPEG0007599472000014.jpg1128 5 shows the graph 500 of the remapping function 500. From JPEG0007599472000016.jpg1440 Produces a direct map to JPEG0007599472000017.jpg935. From JPEG0007599472000018.jpg1110 From JPEG0007599472000019.jpg1166, value By remapping to the range obtained by JPEG0007599472000020.jpg1110, clipping of the waveform can be prevented. The result is JPEG0007599472000021.jpg1237.

[0072] V. Remapping multiple control points simultaneously

[0073] To use the system for multiple control points at once, such as the two sine wave example in Figure 6, the ranges must be evaluated simultaneously to provide assurance that the solver can correctly represent the waveform. As graph 600 shows, two control points can be evaluated simultaneously at different points and with different amplitudes. JPEG0007599472000022.jpg1110 and The two sensors are generated by a 3D scanner, one moving close to the transducer array and the other moving laterally relative to the transducer array. They are set to oscillate in amplitude with two sine waves of different frequencies.

[0074] Violating the power limits of the device results in underpowered and clipped waveforms, and multiple control points exceeding the device limits promote crosstalk between the points, disrupting their haptic effects as compared to a single control point. Graph 700 shown in Figure 7 shows an example of two moving control points. Figure 7 illustrates an attempt to recreate the two control points using a transducer array, resulting in complex and unpredictable behavior as different resonances arise and are replaced between the two points.

[0075] Shown in FIG. 8 is a graph 800 that uses a secondary solver evaluation to determine the power limitations that exist at both points, leading to the shown waveform, which can be written as follows: JPEG0007599472000024.jpg12110Here, one control point moves towards the device potentially becoming stronger, while the other control point moves away potentially becoming weaker across the surface of the device. As can be seen in Figure 8, the dynamic range limitations produced by the solver solution ensure that their waveforms are represented as intended and without crosstalk. Thus, Figure 8 shows two points whose amplitude envelopes have been reduced to a level guaranteed by the solver to be reproduced without error, or with a controlled and known level of error. They are of equal power, since this is the configuration that minimizes the square of the amplitude error at the ideal amplitude level.

[0076] VI. Redistribution of power across multiple simultaneous control points

[0077] The amplitude range can be set to different values ​​and dynamically set through time to manage power levels. In this way, manipulating the range determines the proportion of the array output that is concentrated at a particular control point. This is shown in Figures 9 and 10 with a simple manipulation of power from one control point to another. There are also more complex out-of-phase relationships that allow for dynamics of movement between control points as shown in Figures 11 and 12.

[0078] Shown in Figure 9 is a graph 900 with input waveforms and ranges for two input control points. These have different waveform ranges, which are also used as hints to allow each point to have more or less power available from the device. Shown in Figure 10 is a graph 1000 with the output obtained from the inputs in Figure 9. In this configuration, the power generated by the array is transferred from one control point to the other while they are moving away from and towards the array, respectively. This allows more power to be used than if the control points were expected to have equal amplitudes, as predicted in Figures 6, 7 and 8.

[0079] Shown in FIG. 11 is a graph 1100 with an alternating function where one or a group of control points can be effectively modeled by matching the range to a sine wave. This ensures that each waveform reproduces according to the power requirements of each side of the modulation system and that each side of the system has access to the full array power. This results in time slicing without the possibility of any output clipping occurring. Shown in FIG. 12 is a graph 1200 with the output obtained from the input of FIG. 11 but without the grouping described. This reduces the output dynamic range as the power budget is unnecessarily shared between all the control points.

[0080] 2. Efficient generation of trap points for acoustic levitation

[0081] I. Optimal objective for suspended spherical particles

[0082] The optimization problem required to levitate a particle centers around the use of the Gorkov energy potential. This is a scalar function that evaluates the energy of each possible position of the particle in an acoustic field. Finding a minimum or low energy state at a spatial location corresponds to finding a potentially stable trapping point. However, there are acoustic effects that the Gorkov potential does not take into account, such as acoustic streaming, where acoustic waves induce bulk motion in a fluid. For this reason, it is also preferable to ensure that trapping does not occur in regions of fluctuating pressure that are often determined to be the feasible minimum of the energy potential, and are included in the definition of the potential.

[0083] Because the direction and magnitude of external forces, such as gravity and acoustic streaming effects, cannot be derived from the Gorkov potential alone, it is often preferable to calculate the Laplacian of the Gorkov scalar potential. This allows further control over which axes to favor over others in order to generate a more stable configuration of the inward pushing forces.

[0084] II. Evaluation of the Golkovs potential

[0085] The Gorkov potential field, an energy field that allows spherical particles to be suspended in the presence of infinitesimal energy, is described as (Bruus 2012) as: JPEG0007599472000025.jpg3898 where, JPEG0007599472000026.jpg79 is the speed of sound through an acoustic medium, JPEG0007599472000027.jpg88 is the speed of sound through a spherical particle, JPEG0007599472000028.jpg99 is the density of the acoustic medium, JPEG0007599472000029.jpg89 is the density of a spherical particle, and the time-averaged pressure The left side, including JPEG0007599472000030.jpg79, is the representation of potential energy, and the particle velocity vector The right hand side involving JPEG0007599472000031.jpg79 is conservatively related to the potential energy. Because we make the assumption that the Gorkov potential exists in the active or far field, these pressure and velocity terms can be exchanged with in-phase time harmonic terms. This suggests that the time harmonic Gorkov potential field can be rewritten as: JPEG0007599472000032.jpg13121

[0086] Extending this to many sources gives us: JPEG0007599472000033.jpg23159

[0087] The problem with this is that the cross terms effectively mean that it becomes: JPEG0007599472000034.jpg39155The result is a quadratic, squared increase in the number of transducers.

[0088] III. Potential optimization by eigenproblems

[0089] The potential definition, neglecting time harmonic effects and considering the pressure and particle velocity input to each transducer, can be written as: JPEG0007599472000035.jpg1398JPEG0007599472000036.jpg58106JPEG0007599472000037.jpg1028Here, JPEG0007599472000038.jpg78 is the index of the desired trap, and the extra JPEG0007599472000039.jpg1146 is the transducer JPEG0007599472000040.jpg88 and Trap point from JPEG0007599472000041.jpg76 JPEG0007599472000042.jpg86. The phase and amplitude change as the transducer moves to a non-zero level. Using JPEG0007599472000043.jpg923, the optimal complex transducer activation coefficient to achieve a potential well has the maximum eigenvalue JPEG0007599472000044.jpg is given by 88 eigenvectors.

[0090] If we select JPEG0007599472000045.jpg77 in this way, the objective function JPEG0007599472000046.jpg918 is maximized, which means Due to the negative sign in the definition of JPEG0007599472000047.jpg914, a minimum is obtained in the Gorgov potential field because JPEG0007599472000048.jpg821 is the definition of an eigenvector, It is defined that the image is JPEG0007599472000049.jpg833, and by definition, JPEG0007599472000050.jpg87 is an eigenvalue of the matrix, Since JPEG0007599472000051.jpg822 is maximized, this should yield the maximum value, if it exists. A standard statement of the problem is as follows: JPEG0007599472000052.jpg2161

[0091] However, as previously discussed in the literature, minimizing the energy of the Gorkov field produces a focal point, not actually a trap point. Generating a trap point is about minimizing the pressure fluctuations at the trap. This is non-trivial, as the definition of pressure fluctuations is complex-valued and cannot be meaningfully minimized, but instead represents a quantity of fluctuation that must be zeroed out.

[0092] IV. Creating a null space within the eigenvalue problem to create zero pressure

[0093] The zero pressure condition can be constructed as a set of constraints that evaluate to zero pressure at a given point in the acoustic field. Pressure can only evaluate to a non-zero value at a given location if there exists some part of the solution vector in the high dimensional phase space that is parallel to the known optimal pressure focusing solution. Trap points JPEG0007599472000053.jpg87, Vector We have JPEG0007599472000054.jpg1212, which looks like this: JPEG0007599472000055.jpg12107

[0094] Hence the trap point The constraint for zero pressure in JPEG0007599472000056.jpg89 can be written as a constraint equation: JPEG0007599472000057.jpg1025

[0095] It doesn't matter what the amplitude is for this, the constraint needs to be unit length, so: JPEG0007599472000058.jpg2332

[0096] The standard description of the problem is as follows: JPEG0007599472000059.jpg4577

[0097] Deriving the eigenvalue system from Golub (1973) and restating it for complex values ​​leads to the following argument:

[0098] Given the above system with JPEG0007599472000060.jpg719, the characteristic equation of the system can be written as: JPEG0007599472000061.jpg1092 where, JPEG0007599472000062.jpg77 and JPEG0007599472000063.jpg77 is the Lagrange multiplier. This is JPEG0007599472000064.jpg88 We can differentiate by the conjugate of the vector and get the system: JPEG0007599472000065.jpg943

[0099] Multiplying this by JPEG0007599472000066.jpg910 gives us: JPEG0007599472000067.jpg1561

[0100] JPEG0007599472000068.jpg78, assuming it is unit length, JPEG0007599472000069.jpg921, and in the previous problem definition JPEG0007599472000070.jpg921. Substituting both of these gives us: JPEG0007599472000071.jpg930

[0101] Substituting this into the previous derivative, we get: JPEG0007599472000072.jpg843

[0102] This is as follows: JPEG0007599472000073.jpg1043

[0103] This is a description of a specific system, and JPEG0007599472000074.jpg827 is not symmetric, so computing the eigenvectors is difficult. However, JPEG0007599472000075.jpg818 is idempotent As shown in JPEG0007599472000076.jpg949, the eigenvalues ​​do not change with the ordering, so we get: JPEG0007599472000077.jpg961

[0104] It doesn't hurt to extend this to more constraint vectors, or multiply the left hand side further: JPEG0007599472000078.jpg11110After many multiplications we finally get: JPEG0007599472000079.jpg1462 here JPEG0007599472000080.jpg89 represents the Moore-Penrose pseudoinverse.

[0105] In summary, the solution to the initial problem is This can be found by finding the eigenvector with the maximum eigenvalue of JPEG0007599472000081.jpg846. The multiplication of JPEG0007599472000082.jpg819 can be viewed as adding a rank-1 modification to the null space of the solution matrix, which adds a basis vector.

[0106] Finally, it would be useful to add support for computing the Laplacian of the Gorkov potential in a simple way. The Gorkov potential solution found here produces trap points, but does not allow control over which direction the force is applied from, which may even be directions where there is no force. To achieve this, it is necessary to optimize for the most divergent points in the gradient of the latent field, and then reweight the gradient of this field, which is generally speaking the Laplacian of the Gorkov field.

[0107] V. Modifications for generating the Laplacian of a Gorkov field

[0108] The modification to minimize the Laplacian of the Gorkov field is JPEG0007599472000083.jpg89 matrix only needs to be adapted. Discarding the effect of wavefront curvature from each transducer and considering them as plane waves only removes singularities that complicate the formulation and makes the particle velocity constant. The resulting plane waves from each transducer can be written in the following equation: JPEG0007599472000084.jpg1185

[0109] The matrices can be expressed as equivalent systems, except for the Laplacian case: All constraints for JPEG0007599472000085.jpg19150JPEG0007599472000086.jpg912 are It will be migrated to JPEG0007599472000087.jpg913. Taking the Laplacian of this formulation gives the modified matrix formulation: JPEG0007599472000088.jpg19134

[0110] Using this value instead of the corresponding element of the JPEG0007599472000089.jpg88 matrix gives the solution in terms of optimizing the Laplacian of the Gorkov potential. JPEG0007599472000090.jpg1331 can be written as: JPEG0007599472000091.jpg13136

[0111] This is also JPEG0007599472000092.jpg93 It can also be split into contributions from one direction: JPEG0007599472000093.jpg11122

[0112] In this method, weighting can be applied to generate constraint axes to bias the forces generated from each direction.

[0113] VI. Implementation on device

[0114] The system herein represents an apparatus that can be used to calculate one or more levitation traps based on a simple eigensystem solver, which can potentially be implemented in device firmware to create embodiments that can use acoustic levitation to restrain, trap, and transport items.

[0115] This is shown in Figures 13-16.

[0116] Figure 13 shows the Gorkov Laplacian with uniform weighting. The left side 1310 is the xy plane parallel to the transducer array. The right side 1320 is the xz plane perpendicular to and slicing through the planar transducer array.

[0117] 14 shows the Gorkov Laplacian with the z direction heavily weighted. The right hand side 1420 is the xz plane that is perpendicular to and slices through the planar transducer array.

[0118] Figure 15 shows the Gorkov maximization of the zero pressure condition. The left side 1510 is the xy plane parallel to the transducer array, which describes the helical trap, indicated by the non-uniform shading around the ring of uniform intensity. The right side 1520 is the xz plane perpendicular to and slicing through the planar transducer array below.

[0119] In FIG. 16, additional eigenmodes 1610, 1620, 1630 correspond to trap solutions that also solve the core problem. These are less effective because they do not maximize the objective function. This is because the modified Gorkov potential, which is the value of the objective function, is small for these modes, orthogonal to the principal eigenvector. The objective function value achieved is reflected in the eigenvalues ​​associated with each eigenvector. The eigenvectors associated with small eigenvalues ​​represent these alternative, less effective solutions.

[0120] 3. Phased array device focusing vector effects applied to airborne haptics

[0121] I. active sound intensity

[0122] active sound intensity JPEG0007599472000094.jpg910 is the pressure of the acoustic medium JPEG0007599472000095.jpg88 and particle velocity vector JPEG0007599472000096.jpg77 can be used to define an acoustic wave as follows: JPEG0007599472000097.jpg1136Here, the activation strength is a purely real vector. In the near field of the source, JPEG0007599472000098.jpg1334, but in the far field, the reaction field dominates. JPEG0007599472000099.jpg819. For this reason, Only JPEG0007599472000100.jpg98 is considered, ensuring realism.

[0123] Consider the following monotonic plane wave source: JPEG0007599472000101.jpg1536 which gives: JPEG0007599472000102.jpg1269

[0124] Where: JPEG0007599472000103.jpg87 is the pressure amplitude, JPEG0007599472000104.jpg88 is the particle velocity amplitude in each direction.

[0125] Since the sound intensity, and therefore the energy flux, transmitted by the wave along a given axis is of interest, the dot product to measure the scalar contribution is: JPEG0007599472000105.jpg1582

[0126] Due to the form of the sound wave equation: JPEG0007599472000106.jpg2745The principle of superposition applies to both pressure vectors and particle velocity vectors. JPEG0007599472000107.jpgAssuming 76 sound sources, each approximated as a plane wave at the point of interest, we obtain: For efficiency, it is desirable to solve for the appropriate variables using a linear system. This means that it is not possible to solve for the sound intensity in a given direction, since it consists of two independent variables. Existing techniques use a complex-valued scalar pressure JPEG0007599472000109.jpg1120Solution: We then assume that the direction of each of the individual waves is unimportant, or that they are approximately parallel. Instead, we use the particle velocity vector JPEG0007599472000110.jpg1121, which is important for mid-air haptic feedback applications.

[0127] This produces a solution that reflects the directionality of the waves and the resulting energy transfer, giving a complex-valued particle velocity vector This leads to two potential techniques to extend the current solution mechanisms to take advantage of JPEG0007599472000111.jpg917.

[0128] II. Solution for Complex Values

[0129] The solution for complex values ​​involves using an acoustic model for the values. This model can be queried given the transducer position along with the control point positions to determine the resulting phasor at that point in the field. This is a complex value, JPEG0007599472000112.jpg1044 here JPEG0007599472000113.jpg88 is the index of a sample point in the sound field, JPEG0007599472000114.jpg77 corresponds to the index number of the active transducer: JPEG0007599472000115.jpg3292 where, JPEG0007599472000116.jpg77 is JPEG0007599472000117.jpg829 is the intended phaser in the air at location JPEG0007599472000118.jpg67 is the transducer Transducer activation coefficients from JPEG0007599472000119.jpg929, initial phase and amplitude required to generate the behavior in air It is identified as JPEG0007599472000120.jpg77.

[0130] however, Since JPEG0007599472000121.jpg87 is complex-valued and contains both phase and amplitude specifications, the phase offset can be modified at will. This can be exploited to find the phase value that opens up the largest possible subspace of amplitudes. This is achieved by finding the principal eigenvector y of the system. JPEG0007599472000122.jpg38144

[0131] JPEG0007599472000123.jpg78 is a sound field For each control point in JPEG0007599472000124.jpg829, represent the desired amplitude at each point in the acoustic field, JPEG0007599472000125.jpg77 is a per-column normalization applied to cancel the inner products to 1 in each diagonal of the matrix. The eigenvectors of this matrix The normalized components of JPEG0007599472000126.jpg77 are the initial amplitude JPEG0007599472000127.jpg76 is a phasor vector that when multiplied with promotes reinforcement between each of the pressure points in the field. Corresponds to JPEG0007599472000128.jpg77.

[0132] The minimum norm solution of an underdetermined problem can be solved using the following form: JPEG0007599472000129.jpg958

[0133] The key to this solution is the matrix The inverse of JPEG0007599472000130.jpg921, which is the first step in solving the system via Cholesky decomposition, can be written as: JPEG0007599472000131.jpg32100This is the matrix of the eigensystem. This is also meant to help derive JPEG0007599472000132.jpg77. But you can also manipulate the system in the following ways: This results in a very small linear system and subsequent solution size, which can be solved by any technique (QR decomposition, Cholesky decomposition, etc.) and kept this small until released. This C matrix and this Neither of the JPEG0007599472000134.jpg819 vectors have an index into the transducer until they are computed: JPEG0007599472000135.jpg1050

[0134] There may be many transducers, so at the last possible JPEG0007599472000136.jpg88 is multiplied to multiply the waveform emitted from each transducer, JPEG0007599472000137.jpg725 It is often worth converting back to transducer activation.

[0135] III. Resolving one or more anisotropic control points

[0136] Creating control points in the acoustic field to have a specified normal vector or direction along each of them is one way to incorporate momentum requirements into a phased array solver. Giving each control point a normal vector means that rather than solving for the pressure at each control point, we effectively solve for the x, y and z particle velocities, and therefore the components of the wave momentum, at the points of interest. To achieve that, the linear system equations are expanded from above to include three components for each particle velocity, instead of the one for the pressure defined earlier. Those three components are: It is written as JPEG0007599472000138.jpg1254 and is defined as follows: JPEG0007599472000139.jpg3575In each case, JPEG0007599472000140.jpg1012 is the control point Transducer in JPEG0007599472000141.jpg86 JPEG0007599472000142.jpg76 is a vector showing the direction of propagation from the point across which momentum is measured. JPEG0007599472000143.jpg98 is the control point is the normal vector direction of JPEG0007599472000144.jpg86. The kernel linear system for the solution is: JPEG0007599472000145.jpg69102This shows the three components of the velocity of the particles in the acoustic medium at each control point. However, if we use a reduced representation approach to the solution of the system, i.e., if we specify the C matrix to reduce the amount of work required, the number of entries in this matrix increases by a factor of nine, making it a much more intimidating problem to compute. First, each element of the matrix C JPEG0007599472000146.jpg923 is the submatrix: JPEG0007599472000147.jpg2184

[0137] This is the eigenvalue matrix Using the existing approach for JPEG0007599472000148.jpg86 works fine and can be constructed using the same approach, but the amplitude of the user input is JPEG0007599472000149.jpg1054 and JPEG0007599472000150.jpg1168. The cross terms of the eigensystem matrix can be thought of as making a measurement on one axis at the location of one control point, while considering the effect of adding a waveform on another axis to yet another control point. JPEG0007599472000151.jpg918 now contains three times as many components, Multiplication by JPEG0007599472000152.jpg810 can return the system to the transducer activation coefficients needed to reproduce the desired momentum vector in the air above the device.

[0138] This solution is particularly interesting for the application of airborne haptic feedback, since it opens up the possibility of creating waves approaching the hand from different angles. This allows the acoustic waves to induce shear waves in the skin on the boundary, which opens up a whole new set of sensations that can be created using airborne haptic techniques, at the expense of extra necessary computations.

[0139] IV. Solving for particle velocities at one or more control points in an acoustic medium

[0140] Another technique to achieve the conversion from a solver that targets pressure values ​​to one that can target momentum is to generate a linear system solver that works perfectly with the acoustic particle velocities (as in the original pressure approach) and is also complex-valued for monotonic harmonic times. Moreover, it is a scalar: it determines the direction of the maximum particle velocity at the control point, and thus the resulting wave direction. This can be achieved by computing JPEG0007599472000153.jpg89 and then multiplying to make each matrix entry into a single complex-valued scalar. The dot product is This can be applied as a post-processing to the initial matrix construction, since it affects all sampling directions by the scale factor equally across 730. However, the resulting wave direction is not affected by the solver, since it needs to be known before the solver can work. When applied to building aerial haptic devices, this approach would be most suitable for small devices where airborne, but relatively close-range feedback is best, and is an extension to current solving mechanisms.

[0141] To achieve this, the acoustic model for each transducer is recast as a problem to determine particle velocities in an acoustic medium. Thus, the alpha vector is a system that needs to be solved for a given control point, as the complex-valued particle velocity measured along the direction of maximum motion. JPEG0007599472000155.jpg1150, and the system that needs to solve for a given control point is JPEG0007599472000156.jpg1956. But unfortunately, JPEG0007599472000157.jpg910 cannot be known a priori, but is necessary to construct the linear system matrix: JPEG0007599472000158.jpg3390

[0142] Thus, the acoustic medium can be modified to produce the values: JPEG0007599472000159.jpg3652Here, the formula is as follows. JPEG0007599472000160.jpg2456

[0143] These values ​​are, where: JPEG0007599472000161.jpg77To construct the dot product of the matrices, the following must be used: JPEG0007599472000162.jpg32100

[0144] However, there are two possible routes to achieve this. The necessary dot products are defined as follows: JPEG0007599472000163.jpg26123

[0145] This can be written as follows: JPEG0007599472000164.jpg31129, where JPEG0007599472000165.jpg89 is a dyadic product, which requires maintaining nine separate sums per matrix element (three for the diagonal) to later compute the normal vector. This has many similarities to that of the anisotropic control point approach. This is beneficial in that it provides an initial matrix configuration for both techniques, thus allowing a common implementation infrastructure to support both techniques. However, in many implementations, it may be preferable to maintain these nine separate summation variables by using a two-stage modeling system instead. In such a system, the first stage runs a simplified model that calculates only the velocity amplitude and wave direction, and then sums the control points. JPEG0007599472000167.jpg1499 finds three velocity directions for each. This is a simplified model since it is decoupled from the phase calculation. The summation of these gives the final strongest amplitude component, respectively, which can then be used to calculate the direction over which to sum in the second stage.

[0146] To summarize: JPEG0007599472000168.jpg5164Now in step 2, the model can be computed as follows: JPEG0007599472000169.jpg2456Inside the transducer model, before starting any other stages of the process. The remaining stages proceed in much the same way as before, only the direction in which the velocities act is different.

[0147] This approach enables the application of mid-air haptic feedback to more directly address momentum transfer to the hand, rather than pressure, improving the fidelity of sensations and therefore the measurement of perceptual behavior, which is essential for the expansion of mid-air haptics as a practical technology.

[0148] V. Attached drawings

[0149] Figure 17 shows that the same wave with the same sound pressure level (SPL) produces different effects on a surface depending on the angle of incidence. On the left 1710, a wave traveling parallel to the surface (right) strikes it with some energy E and reflects away with slightly less energy as a result of momentum exchange. On the right 1720, a wave traveling perpendicular to the same surface does not reflect and exchange energy, thus leaving the same amount of energy. As a result, it is not possible to produce a haptic effect.

[0150] In FIG. 18, shown in the top two rows 1810, 1820 are simulations of measurements of the sound field generated by the phased array at two different focal points made by the time-of-flight method in dB SPL, dB SIL and dB SVL, respectively, with light grey being +10 dB, white being 0 dB and dark grey being -10 dB. 0The bottom row 1830 (dB SPL-dB SIL), (dB SPL-dB SVL) and (dB SIL-dB SVL), respectively, with the same color scheme, reveals the differences between each measurement method. Since the SVL field measurement simulation is based on the velocity of the particles in the medium, the waves can be further canceled by the cancellation of the additional vector components. As a result, dB SVL is always lower than dB SIL and dB SIL is lower than dB SPL, but they are similar distances apart in log dB space due to the multiplication of their respective components and the effective geometric mean of the SIL measurements.

[0151] 4. A unique system solution for the quadratic problem in mid-air haptic systems

[0152] I. The Rayleigh-Ritz theorem for quadratic optimization

[0153] Matrix JPEG0007599472000170.jpg89 is a square positive semidefinite Hermitian JPEG0007599472000171.jpg718 matrix. Since the matrix is ​​a square Hermitian matrix, each eigenvalue must be real. Since it is positive semidefinite, these eigenvalues ​​are Must be positive as JPEG0007599472000172.jpg924.

[0154] Then, order the eigenvalues ​​as follows: JPEG0007599472000173.jpg972

[0155] Matrix JPEG0007599472000174.jpg1042 is shown to be the solution to the optimization problem: JPEG0007599472000175.jpg1557 and JPEG0007599472000176.jpg1659

[0156] This is the Rayleigh-Ritz theorem. Therefore, The JPEG0007599472000177.jpg77 vector can be found as the eigenvector corresponding to the detected eigenvalue.

[0157] II. Quadratic optimization for phased array systems.

[0158] Considering a phased array system, the following It is useful to define it as JPEG0007599472000178.jpg88: JPEG0007599472000179.jpg1195JPEG0007599472000180.jpgThis vector solution can be used to drive an array of source transducers using 66 monotonic sources. JPEG0007599472000181.jpg77 is the unit norm To make the eigenvector of JPEG0007599472000182.jpg823, we use the matrix JPEG0007599472000183.jpg88 Define the Hermitian matrix of JPEG0007599472000184.jpg717 as: JPEG0007599472000185.jpg4664

[0159] This shows the formation of the optimization. All we have to do is create the matrix JPEG0007599472000186.jpg88 requires a specific pure quadratic (matrix JPEG0007599472000187.jpg88 must be a symmetric or Hermitian matrix) to input the objective function.

[0160] Hermitian matrix The objective function that generates JPEG0007599472000188.jpg88 has five potentially useful purely quadratic forms, so we use an existing eigenvector solver to maximize them and get a solution that can be used to drive airborne ultrasound haptics:

[0161] 1. Maximum pressure optimization (all SPL): a. Proportional to the square of the pressure, maximize JPEG0007599472000189.jpg89; b. This is a commonly used objective function to drive phased array systems.

[0162] 2. Maximum particle velocity optimization (all SVL) is: a. proportional to the square of the particle velocity b. take into account the directionality of the waves, as this is necessary to obtain a consistent tactile perception.

[0163] 3. But the maximum particle velocity optimization along a given direction (available SVL): a. is proportional to the square of the particle velocity along the given direction JPEG0007599472000191.jpg830; b. This is due to the wave directionality and the (skin) surface Taking into account both the effect of the normals on JPEG0007599472000192.jpg77.

[0164] 4. Maximum sound intensity (energy flux density) optimization (all SILs): 4.a. Proportional to sound intensity JPEG0007599472000193.jpg2134:b. This takes into account the directionality of the waves, as this is necessary to obtain a consistent tactile perception.

[0165] 5. Optimization for maximum acoustic intensity (energy flux density) along a given direction (total SIL): a. Proportional to the acoustic intensity along a given direction b. Maximize the wave direction and the (skin) surface Taking into account both the effect of the normals on JPEG0007599472000195.jpg88.

[0166] III. Optimized for maximum pressure

[0167] A quadratic optimizer that can maximize the pressure at each control point will model the quantity: JPEG0007599472000196.jpg21113For each control point j, we sum the objective function: JPEG0007599472000197.jpg19112This means that the dominant eigenvector is solved. JPEG0007599472000198.jpg817 matrix elements as follows: JPEG0007599472000199.jpg1973 here JPEG0007599472000200.jpg810 is the positive and real amplitude of the wave pressure and control points from the initial release JPEG0007599472000201.jpg76, Transducer Complex phase offset on the wave traveling from JPEG0007599472000202.jpg76 JPEG0007599472000203.jpg819. Weighting JPEG0007599472000204.jpg710 uses a quadratic filter to provide relative amplitude control between control points. JPEG0007599472000205.jpg79 terms are used to reweight the input power The unit norm condition of 819 is provided by this, which can be combined with a quadratic weighting factor and a control loop to accurately address individual pressure amplitudes.

[0168] IV. Optimizing for Maximizing Particle Velocity

[0169] A quadratic optimizer that can maximize particle velocity would model the quantity: JPEG0007599472000207.jpg58123This is done for each control point j. This amounts to summing the following objective function: JPEG0007599472000208.jpg68120This means that the dominant eigenvector is solved. JPEG0007599472000209.jpg718Matrix elements are as follows: JPEG0007599472000210.jpg18125 here JPEG0007599472000211.jpg712 is the control point from the initial emission JPEG0007599472000212.jpg86, Transducer The x-direction part of the induced particle velocity amplitude on the wave traveling from JPEG0007599472000213.jpg76 (this is real and positive) and JPEG0007599472000214.jpg918 is a complex phase offset. Weighting JPEG0007599472000215.jpg79 uses a quadratic filter to provide relative particle velocity amplitude control between control points. JPEG0007599472000216.jpg814 terms are used to weight the input power. This is provided by the unit norm condition of 820. Combining this with quadratic weighting factors and a control loop allows accurate use of the individual velocity amplitudes.

[0170] V. Optimization to maximize particle velocity for a given set of directions

[0171] A quadratic optimizer that can maximize particle velocity for a given set of directions would model the quantity: JPEG0007599472000218.jpg57138This is done for each control point j. This amounts to summing the following objective function: JPEG0007599472000219.jpg67133This makes each element of the n×n matrix for which the dominant eigenvector is to be solved: JPEG0007599472000220.jpg18144, where JPEG0007599472000221.jpg914 is the amplitude of the induced particle velocity JPEG0007599472000222.jpg77 is the directional part (which is a real and positive number), JPEG0007599472000223.jpg717 is the control point from the initial emission JPEG0007599472000224.jpg77 Transducer is the complex phase offset on the wave traveling from JPEG0007599472000225.jpg77. JPEG0007599472000226.jpg78 uses a quadratic term to provide relative amplitude control between control points. JPEG0007599472000227.jpg830. Global amplitude scaling is used to reweight the input power The unit norm condition for JPEG0007599472000228.jpg820 is given by this. In combination with a quadratic weighting factor and a control loop, this can be used to precisely address individual amplitudes.

[0172] VI. active sound reinforcement

[0173] active sound intensity JPEG0007599472000229.jpg88 is the pressure of the acoustic medium JPEG0007599472000230.jpg77 and particle velocity vector JPEG0007599472000231.jpg67 can be used to define an acoustic wave as follows: JPEG0007599472000232.jpg1135Here, the activation strength is a purely real vector quantity. In the near field of the source, the activation field is instead However, beyond one or two wavelengths in the far field, JPEG0007599472000234.jpg817. For this reason, Only JPEG0007599472000235.jpg88 is considered and assumed to be guaranteed reality. Now consider a monotonic plane wave source: JPEG0007599472000236.jpg1436We get the following formula: JPEG0007599472000237.jpg1368 where, JPEG0007599472000238.jpg98 is the pressure amplitude, JPEG0007599472000239.jpg77 is the particle velocity amplitude in each direction.

[0174] Since the sound intensity, and therefore the energy flux, transmitted by a wave along a given axis is also important, a dot product is generated to measure the scalar contribution: JPEG0007599472000240.jpg1581

[0175] By forming the acoustic wave equation: JPEG0007599472000241.jpg2642 The principle of superposition applies to both the pressure vector and the particle velocity vector, of which we are interested, JPEG0007599472000242.jpg67Assuming the acoustic source is approximated as a plane wave, the equation is obtained as follows: JPEG0007599472000243.jpg39106

[0176] In the problem area, JPEG0007599472000244.jpg1028, It is well known that JPEG0007599472000245.jpg1019 is in phase, in which case the placement of the conjugates does not matter, and it is possible to reverse the order of the conjugates to the conjugate pressures instead of the particle velocities. (This can be verified by extending the above equations and using commutativity to reverse the order of the complex powers of each term.) By inference, this is a Hermitian matrix It can be written in a way that produces JPEG0007599472000246.jpg89; JPEG0007599472000247.jpg88126

[0177] VII. Optimized for maximum acoustic reinforcement

[0178] The first quantities that produce scalar values ​​are: JPEG0007599472000248.jpg38134 However, this is assumed to be the far field, which can also be written as: JPEG0007599472000249.jpg37145This is the fourth power of the pressure sum, but because of the structure of complex powers, two of these are the same and can be eliminated using the square root.

[0179] Thus, a quadratic optimizer that can maximize the acoustic intensity would model the quantity: JPEG0007599472000250.jpg18159Each control point For JPEG0007599472000251.jpg77, this would result in the sum of the following objective functions: JPEG0007599472000252.jpg48114This means that the dominant eigenvector is solved. JPEG0007599472000253.jpg818Matrix elements are as follows: JPEG0007599472000254.jpg19151, where JPEG0007599472000255.jpg1013 is the partial direction of the normal vector before the wave JPEG0007599472000256.jpg76, JPEG0007599472000257.jpg810 is the positive and real amplitude of the wave pressure on the wave traveling from transducer q at control point j from the initial emission, and JPEG0007599472000258.jpg818 is a complex phase offset. JPEG0007599472000259.jpg78 uses the square root of a function containing the fourth power of pressure to provide relative particle velocity amplitude control between control points. Global amplitude scaling is used to reweight the quadratic terms to keep the input power This is provided by the unit norm condition of 921. Combining this with a quadratic weighting factor and a control loop allows for accurate handling of individual velocity amplitudes.

[0180] VIII. Optimization for Maximizing Acoustic Enhancement in a Given Set of Directions

[0181] A quadratic optimizer that can maximize the acoustic intensity along a given direction will model the quantity: JPEG0007599472000262.jpg36159However, assuming this is in the distal field, this can be written as: JPEG0007599472000263.jpg37159This is for each control point for JPEG0007599472000264.jpg76. This results in a sum of the objective function: JPEG0007599472000265.jpg81126 This is the dominant The eigenvectors that generate each element of the JPEG0007599472000266.jpg718 matrix should be solved: JPEG0007599472000267.jpg30112 where, JPEG0007599472000268.jpg912 is the induced particle velocity amplitude JPEG0007599472000269.jpg66 is the part in the direction (this is a real positive number), JPEG0007599472000270.jpg710 is the positive and real amplitude of the wave pressure, JPEG0007599472000271.jpg818 is the control point from the initial radiance JPEG0007599472000272.jpg818 Transducer is the complex phase offset on the wave traveling from JPEG0007599472000273.jpg76. JPEG0007599472000274.jpg710 uses a quadratic JPEG0007599472000275.jpg1441 terms are used to reweight the input power The unit norm condition on JPEG0007599472000276.jpg819 is provided by this, which can be used in combination with quadratic weighting factors and a control loop to accurately handle the individual velocity amplitudes.

[0182] IX. Finding the weighting factor

[0183] The Rayleigh-Ritz theorem shows that any of the above matrices can be maximized by finding the dominant eigenvector and increasing the power up to the physical limit of the transducer. This can be reached through power multiplication. While defining each of the control point sums separately, the matrix JPEG0007599472000277.jpg78 is composed of vector When multiplying by JPEG0007599472000278.jpg66, the evaluation of the objective function can be performed as an additional side effect. This allows the maximization to find the amplitude that, when weighted, is appropriate for the solution of the system. can be used to find an appropriate value for JPEG0007599472000279.jpg79.

[0184] The power iteration process is described as follows: JPEG0007599472000280.jpg1535 where: JPEG0007599472000281.jpg1128

[0185] We can use a power iterative process to determine the objective function at each control point.

[0186] To define it: JPEG0007599472000282.jpg1939 and: JPEG0007599472000283.jpg1928At each iteration, instead of calculating: JPEG0007599472000284.jpg834

[0187] That is, the objective function at each control point can be evaluated as follows: JPEG0007599472000285.jpg944

[0188] Finally, as before We can do the calculation for JPEG0007599472000286.jpg89: JPEG0007599472000287.jpg1739The final eigenvector solution is as above. However, when calculated in this way, JPEG0007599472000288.jpg1021 is the control point This gives us the intermediate result JPEG0007599472000289.jpg77. To fix the level that represents the desired control point output, we use a set of The average objective function value is taken to allow the value of JPEG0007599472000290.jpg820 to be directly compared with the objective value given by the user. In this case, the weighting coefficients can be updated in the next iteration using a derivative-free method aimed at optimizing expensive functions, such as the Nelder-Mead algorithm.

[0189] After several iterations of derivative-free optimization combined with the power method, intermediate results JPEG0007599472000291.jpg78 should be used to drive the transducer. This therefore allows the system to track the desired value of the objective function over time while allowing for some error. The advantage of this system is that it is a very simple way to drive a transducer array. The biggest disadvantage is that it requires a lot of calculations, but the calculations are effectively just complex multiplications and are easily implemented with sufficient hardware.

[0190] X. Attached drawings

[0191] Shown in Figure 19 are: Top Left 1910: Control points created by solving a quadratic particle velocity problem using non-orthogonal vectors; JPEG0007599472000292.jpg712 Shown as a cut through the plane. Top right 1920: JPEG0007599472000293.jpg712 The same control points, viewed as a cut through the plane. Bottom left: 1930 By comparison, these are the control points created by solving a quadratic pressure problem. JPEG0007599472000294.jpg812 Shown as a cut through the plane. Bottom right 1940: JPEG0007599472000295.jpg711A cross section through the solution of the pressure problem across a plane.

[0192] 5. Control regions for phased arrays specified via a quadratic problem.

[0193] This disclosure describes a method for automatically generating basis functions using user-supplied design parameters. These are input to an efficient optimization over quadratic functions in a space of given transducer bases set online. This gives the user an accessible design space through online optimization that can produce a variety of different effects. A matrix can be constructed by evaluating a binomial product over the transducer effects using a linear operator called the "quadratic control region" since it is a control region for a quadratic problem. The principal eigenvector of this matrix may be selected as the best solution from a set of orthonormal solutions, which are the eigenvectors of this matrix. The optimization produces as output transducer activation coefficients - the raw complex values ​​needed to drive the individual transducer elements to reproduce a sound field or phenomenon with the desired characteristics.

[0194] These raw complex values ​​can then be used to further comprise the basis functions of a linear system for subsequent use, the solution of which allows the user to simultaneously produce all desired effects while reweighting the output power available to each.

[0195] I. Linear Control Region

[0196] The linear control region is a complex-valued weighting function is defined as an operator that performs an integral with a linear acoustic quantity such as pressure using the JPEG0007599472000297.jpg1397This is used to evaluate the acoustic field. It should be noted that when JPEG0007599472000298.jpg825 is paired with a similar phase distribution in the acoustic quantities, the complex conjugate is used to maximize the function evaluation. Any approximation to this integral (sampling, etc.) can be used, or in fact it can be used to approximate any linear operator over the sound pressure field. If it is also linear, a further integral can be constructed to optimize the control region using the vector particle velocity: JPEG0007599472000299.jpg1498 where, JPEG0007599472000300.jpg67 is the particle velocity of the medium, JPEG0007599472000301.jpg66 is a three-dimensional vector field, with each dimension being a complex number, which acts as a weighting function for the same purpose as above. In this method, the various patterns of acoustic quantities are evaluated linearly and the phase A particular phase of the result is maximized when it matches a "template" created by a scalar or vector function of JPEG0007599472000302.jpg66.

[0197] Also, It should also be noted that JPEG0007599472000303.jpg825 may be near zero (or composed of zero vectors) to create a single estimate of multiple seemingly disjoint regions.This formulation can be extended to any linear operator over a field of linear acoustic quantities, so that any integer approximation (e.g., sampling) can be used.

[0198] A phased array can be described by a set of complex-valued coefficients that control each transducer: To maximize the amount of audio encoded into JPEG0007599472000304.jpg78, we use the following input: Driving an array of acoustic transducers at a fixed frequency with a 66 vector can be expressed as a linear algebra problem, the intent of which is as follows: JPEG0007599472000306.jpg2150

[0199] The problem statement is eigenvectors, here This is also the definition of JPEG0007599472000307.jpg66, so this is a matrix This can be solved by taking the principal eigenvector of the matrix Constructing the JPEG0007599472000309.jpg78 is accomplished by first considering the linear acoustic quantities emitted by each transducer as basis functions in the global space of acoustic quantities of interest. JPEG0007599472000310.jpg831 is a method for achieving maximization of pressure or other acoustic quantities at a point, represented by its absolute position. However, by using the control area integrals described, it is possible to govern the acoustic quantity radiated into the field of each individual transducer element. JPEG0007599472000311.jpg922 is JPEG0007599472000312.jpg920, resulting in a column vector: JPEG0007599472000313.jpg2138

[0200] By invoking linearity, we can multiply the complex activation values ​​that drive each transducer to derive a vector: JPEG0007599472000314.jpg2148 where, JPEG0007599472000315.jpg67 is the number of transducers. Given that JPEG0007599472000316.jpg68 is the result of a linear operator acting on a linear field, it can be Hermitian transposed and multiplied by itself to produce a quadratic form, ultimately resulting in: JPEG0007599472000317.jpg957

[0201] Weighting Function Modify JPEG0007599472000318.jpg825 and add the matrix By solving the eigenvectors of JPEG0007599472000319.jpg78, a variety of different beams, surfaces and other acoustic phenomena can be generated and maximized to generate regions within the acoustics that are optimally actuated.

[0202] The null space generation technique discussed above ("Eigensystem Solution for a Second-Order Problem in Aerial Haptic Systems") is also used to extend this method of describing acoustic levitation, as well as to enable null steering and aerial texturing by introducing null points into the field of points of the describing structure. It should be noted that this process increases the rank of the matrix, so no simplifications are possible with its use.

[0203] Equivalently, for a number of control domains, JPEG0007599472000320.jpgMultiple optimizations can be performed by adding together 78 matrices, which can also include matrices representing control points. This produces a matrix with an increased rank, which is not possible in the present case, although some simplified methods are used.

[0204] By taking the results into account and iteratively reweighting each point, the system can be modified to converge to the required geometric configuration of acoustic quantities, which is typically necessary when amplitudes are not configured at the output.

[0205] This is achieved by summing the generated acoustic quantities and reweighting the system to more closely approximate the result in the next iteration. For example, the iterative reweighting for an integration point of a sub-split volume or a single pressure point might look like this: JPEG0007599472000321.jpg1587 where, JPEG0007599472000322.jpg814 is the next coarse iteration of the technique, JPEG0007599472000323.jpg76 is the current iteration, JPEG0007599472000324.jpg77 is the desired pressure (which may be complex), and the denominator refers to the total pressure obtained by the current iteration. The fraction may be rendered as real only, to affect only the amplitude and not the phase of the resulting field.

[0206] The Bessel beam is weighted in the direction of travel by a weighting function This is generated by generating a unit pressure while modifying the phase angle forward as the beam moves forward. This is sampled at points as long as the angle through which each successive point passes is small enough to sample the phase along the beam axis.

[0207] The sampling of the phase is done in terms of the angle it travels through, It is important to have more than two phase samples per JPEG810. This is because the increments along the path are In the case of JPEG0007599472000327.jpg66, the wave direction cannot be determined. This is also because, in many cases, for each wavelength of the monochromatic frequency distance that the phase covers, JPEG0007599472000328.jpg means that the phase must increase forward to complete a rotation of 810 radians. This is an approximation to a path integral of the form described in the definition of the control region.

[0208] As a result of maximizing the integral written in this way, the transducer will produce an optimal Bessel-like beam by its theoretical capabilities for a simple beam defined progressively in a straight line moving away from the array, and will produce optimal pressure in this configuration. Shown in Figure 20 is a cross section of the field at the edge of the near-field region parallel to a planar linear transducer array, revealing the structure of the complex-valued transducer activation coefficients. A simple progressive phase beam gives a Bessel function approximation (left 2000), while an undersampled helical beam gives a Bessel function with dipole behavior (right 2010), where the dipole creates two structures that activate in phase opposition to each other. Due to the effects of transducer resolution limits, finite array size, and other physical considerations that the optimization must avoid, the beam may not be an exact Bessel beam, but an optimization that takes into account the physical limitations of the real transducer array may work more effectively than a simple evaluation of the Bessel beam function. For large arrays, this closely follows the Bessel function, and FIG. 21 shows a graph 2100 of the Bessel function plotted with pressure when considered radially from the center of the array of FIG.

[0209] A further property of this system is that transducers of any position and orientation with any generated monochromatic field function may be used, which is a key advantage in Bessel beam generation over sampling methods, which by definition only work for planar arrays.

[0210] By creating multiple out-of-phase paths, phase shifts and phase singularities are implicitly created by winding a spiral wave around the path. This path is also curved, resulting in a high degree of configurability, while using a simple quadratic optimization built into the eigenpath approach. This is achieved by inducing beam rotation through the definition of separate phase windings on many parallel beamlines that are closer than the beam width but out of phase with each other. Higher order Bessel beams can then be formed by setting nearby phases to cancel each other at the centerline of the beam, as shown in Figure 22. In Figure 22, two or more samples of progressive phase around a circle (and thus the phase of the spiral arrangement) give rise to the formation of a Bessel function with angular momentum as seen on the left 2210. To identify the phase around the spiral, JPEG0007599472000329.jpgBy using a value greater than 810 radians, the right side (2220) JPEG0007599472000330.jpg78 beams can be generated using higher order Bessel beams.

[0211] These phase lines are an example of how disjoint domains can be defined as a single domain optimization, sampled with points or any other integral approximation, working to nonzero values ​​only within the volume to optimize. Higher order Mathieu-type elliptical helical beams can be generated by introducing phase beams, where the phase increases incrementally and at intervals around the ellipse. As shown in Figure 23, elliptical helical beams with similar properties to higher order Bessel beams, such as Mathieu-type beams, can be created by parameterizing the circumference of the phase ellipse. JPEG0007599472000331.jpg Determining the ratio of 711 radians and assigning a phase angle proportional to the circumference of the elliptical chord produces a suitable configuration for beam optimization. Image 2310 on the left shows a higher order Bessel reshaped into an ellipse near the Mathieu beam. Image 2320 on the right shows an approximate Mathieu beam with a bridged, helical structure beginning to form in the central region. Both JPEG0007599472000332.jpgIt is based on a 710 radian phase circle, which creates a helical wavefront.

[0212] Bessel family beams generated in this way may be modified to generate more exotic beams that exist in more restrictive forms in the literature. Various results may be reproduced by creating arbitrary paths in the field that are traversed by physically plausible phase functions that proceed in a manner appropriate for wave traversal, but because the optimization involves eigensystems in a more general form, the initial bandwidth cost is lower than before, since the beams need only be represented by integrals.

[0213] Taking a curved path and moving the phase forward produces a curved beam that is commonly used to steer the beam around obstacles generated by a 45° circle segment. Shown in Figure 24 is an image 2400 of a curved beam formed from a phase function that follows 45° segments of a circle connected in series along the line using a crude iterative method to equalize the sound pressure along the beam. The transducer array is positioned along the bottom edge of the image.

[0214] Further techniques involving multiple phase paths can produce a "bottle beam", i.e. a beam with a section cut out of the center, which can facilitate levitation of large objects. This is produced in FIG. 25, which shows an image 2500 of twinning two Gaussian curves along which the phase has been traced. The bottle beam produced by creating a phase function follows the Gaussian curve using a crude iterative method to homogenize the acoustic pressure along the beam. The transducer array is placed along the bottom edge of the image. The dotted nature of the lower region suggests that the phase velocities input to the optimization are too slow, and that the beam could be improved by using faster phase velocities for the initial section. If the object is sufficiently stiff, only a portion of the space needs to be a high acoustic pressure region, so a portion of the curve can also be used as in FIG. 26. FIG. 26 shows an image 2600 of an incomplete bottle beam formed from a portion of a circle defined by a progressive phase, which can be used to levitate even larger objects by optimizing the areas of high acoustic pressure in a non-uniform manner. Provided that the object or objects to be levitated have sufficient cohesive strength so as not to flow, be pushed out, or escape through the area of ​​low sound pressure, it is not necessary for the force exerted on a large object to be uniformly distributed.

[0215] Bessel-like beams can be generated by moving a path forward faster than the phase velocity. In this way, as shown in FIG. 27, the beam is optimized with a faster phase velocity than normal, resulting in a longer wavelength, creating a beam whose phase velocity appears to move faster than the speed of sound. Shown in FIG. 27 is an image 2710 of a wavefront traveling faster than the speed of sound, which can be generated by performing a phase path that moves through a corner, paradoxically generating a longer wavelength at the same frequency. The implicit optimization infers that the main beam is composed of intersecting beams, causing a specified increase in phase velocity. As can be seen in the image 2720 on the right, the wavelength is manipulated to be approximately 40% longer at the same frequency for the first part of the beam. Serpentine beams can be constructed by sampling progressive phases along a sinusoid. Shown in FIG. 28 is an example 2810, where a short portion of a sinusoid is defined, followed by a longer 2D serpentine beam of 2820, composed of progressive pressure phases that follow the sinusoid on the phase array. These show that the expansion of a meandering sinusoidal beam by simply adding further segments of progressive phase will naturally expand the beam as far as required (and allowed by the physical placement of the transducer). By making the beam meander in all three dimensions to generate a helically meandering beam, one can generate a beam that can be used to exert forces in either a normal or opposing direction, potentially against a solid object in a field much larger than the wave radiation, wavelength, in a manner similar to the Archimedes screw system. This gives a different mechanism by which the beam can be used to levitate objects. Figure 29 shows JPEG0007599472000333.jpg712 The left image 2910 shows a detailed cross section in the plane; JPEG0007599472000334.jpg712 and a right image 2920 showing the cross section of the beam in the plane. JPEG0007599472000335.jpg946 is bent, Using this technique, which follows a circular path that continues on into the 917, a true helical shape can be created, which can be used to exert forces on objects much larger than the wavelength, in a tractor beam-like fashion, to levitate them.

[0216] If the eigensystem matrix is ​​defined in a single control domain, then the matrix is ​​rank 1 and is a complex-valued vector JPEG0007599472000337.jpg68 is formed by the two-term product of itself and itself, and the matrix is In this particular case, the single eigenvector of the system is easily known. JPEG0007599472000339.jpg88, JPEG0007599472000340.jpg816 means that in this case we identify the complex-valued transducer activation coefficients ourselves by evaluating the objective function, which we define by simply evaluating the control field function for each transducer, the activation coefficients.

[0217] As before, for the quadratic form, we also seek to evaluate the validity through reweighting for each iteration. This is done by evaluating the pressure at each point or the integral of the pressure over each subregion and dividing by the desired pressure to generate new function values: JPEG0007599472000341.jpg928 small area This is achieved by updating the weights of JPEG0007599472000342.jpg76. Then The single eigenvector, also JPEG0007599472000343.jpg88, is then re-evaluated, producing an equivalent approach to repeated iterations on the eigensystem from before, but at greatly reduced computational cost.

[0218] To achieve this with the particle velocity of the medium, a test of the particle velocity can be used, which will change direction until it reaches a maximum: JPEG0007599472000344.jpg1493Pressure can also be monitored using the old technique, which produces a hybrid method, where you specify a weighting function to update the pressure as before, but with an effective control region: JPEG0007599472000345.jpg15118 Scalar but variable complex value JPEG0007599472000346.jpg825 is obtained, and the invariant weighted normal vector JPEG0007599472000347.jpg77 and particle speed Pair with JPEG0007599472000348.jpg77. This approach can be optimized to find a satisfactory pressure without changing the desired wave direction. This is shown in Figure 28 with certain modifications explained in the next section.

[0219] II. Wavefront filtering for directional control points and regions

[0220] A further technique required to create a directional output is to JPEG0007599472000349.jpg77 and normal vector or particle velocity vector If it is the opposite, discard the transducer data and set it to zero: JPEG0007599472000351.jpg2467This is about the points in the area I'm thinking about.

[0221] Figure 30 shows a left image 3010 of two velocity vector samples further defined from a transducer array, and a right image 3020 of two velocity vector samples defined close to the array, with the transducer having a wavefront in the opposite direction to the released and defined control direction. Using velocity weighting functions specified in these two disjoint domains constitutes a single domain in each case, which is optimized for a directional vector representing the particle velocity of the medium at two separate points.

[0222] Trying to find the maximum amount of velocity in 66 directions means that the resulting particle velocity direction is completely Note that this does not necessarily mean that the wavefront direction is in the JPEG0007599472000353.jpg66 direction. Figure 30 shows the result of filtering an opposing transducer providing a wavefront in the opposite direction to the desired direction, while Figure 31 shows the left image 3110 and the right image 3120, which are the results when the transducer wavefront direction is not filtered using the wavefront reference. The opposing transducers are free to interfere in the case of Figure 31, and the pathology is evident horizontally in the form of standing wave structures. This is because the wavefront direction conveys directional information that is missing from the particle velocity, which does not carry information about which direction the wave originated from, allowing for undesired interference.

[0223] This is also important not only for control regions, but also for the simple case of regular control points. The control point case is effectively a simplified control region, approximated with a single point sample of particle velocity and no weighting factors. In both cases, not excluding "useless" transducers leads to similar results regarding the component of the influence of the "standing wave" phenomenon. This is because a particle velocity vector with complex components does not allow to identify the difference between waves traveling in different directions: in this model, a wave traveling from left to right with phase 0 will have a phase JPEG0007599472000354.jpg66 is mathematically the same as a wave traveling from right to left. As a result, if this is not handled correctly in the control point or region based generation, the optimization will attempt to use both, producing solutions with standing wave components that do not contribute meaningfully, and creating regions of unnecessary, extra field-free nulls. This is an important point, because control points input as basis functions into linear optimizations, eigensystems, and even simple equations without this optimization step will exhibit these artifacts.

[0224] These described techniques are integrated into a device to create an acoustic phased array system that can create arbitrarily configured regions and beams by iterating a single eigensystem as defined herein, or in simple cases, by creating a single point or region. In the case of a single point or region, optimization is applied to generate a simple formula showing that the amplitude and phase required to drive the transducer elements are the amplitude and phase represented by evaluating the point or region. Iterative and time-dependent application of the global method is altered and modified in a time-dependent manner to produce various haptic effects. This applies to both the desired field and the derived transducer coefficients, since they are functionally related.

[0225] III. Possible claims 1. A method for monotonic acoustic fields, comprising: estimating the distribution of sound pressure to a single value by integrating the product of said sound pressure field function and a weighting function, both functions being complex-valued and derived from physical transducer properties; applying zeros or null points to be added to null positions of said matrix; using this as the objective function to be maximized in an optimization, the results of said optimization taking the form of eigenvectors of a matrix formed from cross-multiplication terms of said objective function for each transducer; The method includes summing one or more of such matrices to generate an eigensystem having one or more eigenvectors generating said particular field. 2. A method for monotonic acoustic fields, comprising the steps of: Estimating the distribution of acoustic particle velocities into a single vector value by integrating the dot product of the actual 3D weighted vector field, the particle velocity field, and a weighting function, both functions being complex-valued and derived from physical transducer properties; applying zeros or null points to be added to null positions of said matrix; using this as the objective function to be maximized in an optimization, the results of said optimization taking the form of eigenvectors of a matrix formed from cross-multiplication terms of said objective function for each transducer; The method includes summing one or more of such matrices to generate an eigensystem having one or more eigenvectors generating said particular field. 3. A method for monotonic acoustic fields, comprising: Estimating the distribution of acoustic particle velocities into a single vector value by integrating the dot product of the actual 3D weighted vector field, the particle velocity field, and a weighting function, both functions being complex-valued and derived from physical transducer properties; setting to zero the contribution of the transducer to the wavefront at points having opposite directions to the normal that are identified in each estimated subdomain of the integral; applying zeros or null points to be added to null positions of said matrix; using this as the objective function to be maximized in an optimization, the results of said optimization taking the form of eigenvectors of a matrix formed from cross-multiplication terms of said objective function for each transducer; The method includes summing one or more of such matrices to generate an eigensystem having one or more eigenvectors generating said particular field. 4. A method for monotonic acoustic fields, comprising: estimating the distribution of sound pressure to a single value by integrating the product of said sound pressure field function and a weighting function, both functions being complex-valued and derived from physical transducer properties; applying zeros or null points to be added to null positions of said matrix; using this as the objective function to be maximized in an optimization, the result of which is the eigenvectors of a matrix formed from the cross-multiplication terms of the objective function for each transducer; A method comprising: using the contribution of each transducer to describe an optimum phase and amplitude for an input signal of said transducer. 5. Estimating the distribution of acoustic particle velocities into a single vector value by integrating the dot product of the actual 3D weighted vector field, the particle velocity field, and a weighting function, both functions being complex-valued and derived from physical transducer properties; applying zeros or null points to be added to null positions of said matrix; using this as the objective function to be maximized in an optimization, the results of said optimization taking the form of eigenvectors of a matrix formed from cross-multiplication terms of said objective function for each transducer; A method comprising: using the contribution of each transducer to describe an optimum phase and amplitude for an input signal of said transducer. 6. Estimating the distribution of acoustic particle velocities into a single vector value by integrating the dot product of the actual 3D weighted vector field, the particle velocity field, and a weighting function, both functions being complex-valued and derived from physical transducer properties; setting to zero the contribution of the transducer to the wavefront at points having opposite directions to the normal that are identified in each estimated subdomain of the integral; applying zeros or null points to be added to null positions of said matrix; using this as the objective function to be maximized in an optimization, the results of said optimization taking the form of eigenvectors of a matrix formed from cross-multiplication terms of said objective function for each transducer; A method comprising: using the contribution of each transducer to describe an optimum phase and amplitude for an input signal of said transducer. 7. The method of claims 1 to 6, further comprising a device implementing the method and using the results to drive an acoustic transducer. 8. The method of claims 1-6, further comprising a device implementing the method and using the results as part of a linear basis set used to drive an acoustic transducer. 9. The method according to claims 1 to 6, wherein an acoustic beam is generated. 10. The method according to claims 1 to 6, for generating an acoustic quantity. 11. The method according to claims 1 to 6, for generating an acoustic wavefront. 12. A method for an acoustic field comprising the steps of: estimating the acoustic particle velocity as a single value at a point using a linear combination of fields obtained from physical transducer properties; setting to zero the contribution of the transducer of the wave front at the control points having an opposite direction to the normal that is specified in each estimated integral sub-domain; and using this as a basis function in an optimization problem to find the phase and amplitude for the transducer input signal. 13. A method for an acoustic field comprising the steps of: estimating the acoustic particle velocity as a single value at the point using a linear combination of fields obtained from physical transducer properties; setting to zero the contribution of the transducer of the wave front at the control points having an opposite direction to the normal that is specified in each estimated integral sub-domain; using this as the objective function to be maximized in an optimization, the result of which is the eigenvectors of a matrix formed from the cross-multiplication terms of the objective function for each transducer, and then finding the phases and amplitudes of the transducer input signals. 14. A method for an acoustic field comprising the steps of: estimating the acoustic particle velocity as a single value at the point using a linear combination of fields obtained from physical transducer properties; setting to zero the contribution of the transducer of the wave front at the control points having an opposite direction to the normal that is specified in each estimated integral sub-domain; The method includes using the contribution of each transducer to describe an optimal phase and amplitude for the transducer's input signal. 15. The method of any one of claims 1 to 14, wherein the transducer generates a haptic effect when the resulting structure is applied as an input signal and modulated. 16. The method of any one of claims 1 to 14, wherein the transducer generates a haptic effect when the obtained result is applied as an input signal and modulated over time.

[0226] 6. Evaluating complex-valued transcendental functions, multiplication and division in hardware

[0227] instead of JPEG0007599472000355.jpg1134 For the purpose of converting between JPEG0007599472000356.jpg1344, the BKM algorithm is more useful in physical hardware implementations.

[0228] I. Introduction

[0229] Many techniques exist for the computation of trigonometric and other special functions in hardware. One such method is the BKM method, which exchanges complex logarithms with complex powers and vice versa, generating these special functions as a side effect of moving from one representation to another. Of particular interest are those dealing with complex-valued functions, since one often encounters situations where such functions are known to be needed.

[0230] An efficient method of moving between logarithmic and power forms allows multiplication and division to be performed in logarithmic space as simple addition and subtraction. Due to the nature of the process, this may also be exploited to easily generate transcendental functions related to power functions. Since this can be generated from a single algorithm, implementation size is smaller as only one element needs to be implemented and only one method needs to be implemented and optimized, reducing development time.

[0231] The existing algorithm, BKM, implements this approach, but it contains three disadvantageous requirements that are largely or completely eliminated in the proposed method:

[0232] In particular, the partitioning of hardware makes it difficult to implement pre-processing and post-processing steps on data input, output, and in phase.

[0233] Non-rational values ​​and factors that cause rounding problems for efficient implementation.

[0234] Additional operations required to integrate with floating-point types.

[0235] The next section explains how the proposed method has been designed to be free of these drawbacks.

[0236] II. Motivation and narrowing the scope of comparison small

[0237] The process is Euler's equation It is concentrated in JPEG0007599472000357.jpg944, JPEG0007599472000358.jpg87 is a JPEG0007599472000359.jpg855 is a complex generalization, so it is a further coefficient. This is the complex function that the BKM algorithm embodies, and the process it performs is The algorithm is limited to the complex values ​​within the range where convergence is guaranteed. Range reduction techniques are often used to fit JPEG0007599472000361.jpg715. This is It is not a function that has a suitable implementation, insofar as it is likely easy to evaluate using the proposed method by pre- or post-processing with an irrational scaling factor to achieve JPEG0007599472000362.jpg732. In many cases, this irrational scaling factor may be removed to accommodate the intended workflow, such that inputs and outputs are multiples of rational numbers that make physical or computational sense.

[0238] When processing complex powers, the power is split into real and imaginary parts, which use the same input data but produce different results, which is in some sense similar to converting from Cartesian to polar coordinates.

[0239] If we consider only the real part of the exponentiated logarithm, this is the Cartesian-polar transformation JPEG0007599472000363.jpg contains all the coefficient information, just like the 77 components. This is traditionally an irrational number raised to the power of the real part. JPEG0007599472000364.jpg76. To feed it to the algorithm, this must first undergo range reduction, which must find the remainder, an irrational number, so that it falls within the convergence region. JPEG0007599472000365.jpg means dividing modulo 711. This involves hardware division (the 1994 paper recommends the SRT division technique) and therefore this adds a very expensive technique in terms of the logic required to perform the division along with the algorithm.

[0240] The premise of the proposed algorithm is that powers of 2 should be used to simplify range reduction. With this motivation, the real part should be interpreted as a power of 2, i.e. base 2, raised to a given power that is the real part. The integers in the real logarithm input can be referred to bit shifts on the output, and only the fractional parts of the relative binary logarithm input need to fall into the convergence region. The integer parts of the real parts of the complex logarithms are thus rounded and raised separately.

[0241] If we consider only the imaginary part of the logarithm to be exponentiated, this involves the conversion from rectangular coordinates to polar coordinates. JPEG0007599472000366.jpg87 components, all of which contain the argument or angle information. This is also traditionally JPEG0007599472000367.jpg is a power of 76, but here JPEG0007599472000368.jpg76 is not a problem. The imaginary logarithm part can be interpreted as an angle, so the problem is this time JPEG0007599472000369.jpg66, but rotational symmetry can be exploited to limit the range of convergence I need to find the remainder from the integer fraction of JPEG0007599472000370.jpg810.

[0242] Base of the logarithm Instead of JPEG0007599472000371.jpg76, you could use base 2. This would be This gives the extension JPEG0007599472000372.jpg866, which is an even less desirable irrational value, but If you use JPEG0007599472000373.jpg710 specifically to represent the number of rotations Unfortunately, the use of this basis leads to a method that does not converge. The best tradeoff between technique and convergence range is to use the base JPEG0007599472000375.jpg1110 is obtained. This means that each integer value represents a quadrant. 4 Each multiple is a rotation of the whole. Bitwise arithmetic can then be used to extract the rotation number. The quadrant number can be applied to the output as real and imaginary and a sign switch. Again, we obtain a technique in which the integer part of the imaginary component of the complex logarithm is taken by rounding and thus processed separately. Again, we obtain a technique in which the integer part of the imaginary component of the complex logarithm is taken by rounding and thus processed separately.

[0243] This is because the fractional part of the logarithm is still The powers remain within the range of JPEG0007599472000376.jpg874, which is within the convergence range of the main algorithm. JPEG0007599472000377.jpg930 Error in exponentiation plotted for input to a log-power process, where the right side is positive real numbers and the bottom of the image is positive imaginary numbers. With 50 fractional bits and a similar number of iterations, the black areas are 2 -45 The left panel 3210 shows the convergence to a precision of over 2 -45 The right panel 3220 shows in white the final required convergence region Showing JPEG0007599472000378.jpg975.

[0244] The main algorithm for the fractional part of the logarithmic input is described below.

[0245] The inverse process is also accomplished by taking the value with its exponentiated real and imaginary parts and reducing it to its original logarithm. First, find the quadrant in which the complex value resides, offset by 45 degrees. This is accomplished by determining which component has the largest absolute value; real or imaginary and its sign. This identifies which of these quadrants is offset by the angle in which the value resides. This is then used to calculate a rotation to move the value into the first quadrant, whose area is centered on the positive portion of the solid line, as shown in the center panel 3320 of Figure 33. Figure 33 shows the area of ​​the complex value within the region. JPEG0007599472000379.jpg shows the error of exponentiation plotted against the input to a power-log process of 930, where the right side is positive real numbers and the bottom of the image is positive imaginary numbers. With 50 fractional bits and a similar number of iterations, the black areas are 2 -45 The left panel 3310 shows the convergence region, where the grey indicates the convergence accuracy of 2 -45 The middle panel 3320 shows the quadrants centered on the solid line using inverted greyscale. The right panel 3330 shows, in white, the final required convergence region. Showing JPEG0007599472000380.jpg1173.

[0246] The imaginary part of the logarithmic output is preloaded with the equivalent imaginary integer part of the output complex logarithm. Further testing is required to see that the 45° offset causes a complex branch cut into the negative part of the real line, but if the two's complement system wraparound is used this can be ignored and given the required bits the imaginary part of the logarithm measured in quadrants is: This can be easily restricted to be within the range JPEG0007599472000381.jpg930.

[0247] At this point, we are sure that the value is in the first quadrant, so we know that the real part is greater than or equal to the imaginary part in absolute value, and we know that the real part is positive. To find the integer part of the real part of the logarithm, we can count the leading zeros in the fixed-point representation. Subtracting this count from the number of bits in the fixed-point representation is effectively JPEG0007599472000382.jpg1156, where k depends on the representation of the bits and is also the representation of the integer part of the binary logarithm, which can be immediately added to the real part of the logarithm register. As shown in the right panel of Figure 33, the existing data, both real and imaginary, are then scaled to the interval We can shift up by the number of bit positions needed to fit JPEG0007599472000383.jpg828, determined by the leading zeros in the initial count and k which depends on the representation of the bits. Note that this may not be the closest integer representation of the output real logarithm, but it is an integer and will get the rest of the process into the convergence region in the correct format.

[0248] At the end of this process, the integer part of the output complex logarithm is determined, and the fractional part of the output and the input range JPEG0007599472000384.jpg877 remains undetermined, which falls within the convergence range of the main algorithm. A discussion of the main algorithm for the substantially fractional part of the logarithmic output is described below.

[0249] This result is the function JPEG0007599472000385.jpg1144 is the modeled process, as can be seen in Figure 34, where the plotted error is plotted only for fractional parts, well within the power function range of the power part 3410 and logarithm process 3430, along with the corresponding parts that make up the final composite value of the power part 3420 and logarithm process 3440. The noisy grey indicates the region where correct convergence is dominated by small numerical errors in the error measurement, showing that convergence is achieved everywhere for the power function 3420 and logarithm 3440.

[0250] Changing the base of the result means that it has all the properties of logarithms and powers, i.e. you can add and subtract the real parts of logarithms, and multiply and divide by powers. Similarly, angles behave in the same way, it's just that they are represented as integer multiples of the quadrant, i.e. 4 for rotations, values ​​that are easily manipulated in hardware logic. In logarithmic form, adding and subtracting imaginary values ​​is equivalent to rotating clockwise and counterclockwise around zero in the exponentiation space of complex numbers. For complex values ​​that combine real and imaginary parts, regardless of the fact that they are composed of different bases, they can be used to implement complex multiplication and division.

[0251] These results can be trivially extended to floating-point representations, since the required preprocessing is already largely accomplished by the mantissa and power format. If the JPEG0007599472000386.jpg76 values ​​are absolutely needed, and most of the time they are not needed, then the real part of the output logarithm is JPEG0007599472000387.jpg711 and the imaginary number is JPEG0007599472000388.jpg is multiplied by 127. Similarly, the logarithmic data is To import in the format JPEG0007599472000389.jpg76, invert the value. JPEG0007599472000390.jpg1422, then we can multiply component-wise to transform, and then do an initial range reduction to exponentiate in this way.

[0252] III. General Procedure

[0253] The new algorithm proceeds in a similar manner as the BKM technique to find or use the fractional part of the logarithm. There is a set of 896 values, from which each algorithm step In JPEG0007599472000392.jpg67 JPEG0007599472000393.jpg88 is selected. The selection of whether to use JPEG0007599472000394.jpg88 is made by comparing the truncated number to the set bit number, which may be a subset of the bits from the input register. When JPEG0007599472000395.jpg88 is selected, A lookup table of pre-computed logarithms is needed to balance the bit shifts (potentially with rounding applied to the most significant unused bits) and additions / subtractions, and a multiplication is applied to a register holding the power value, with the equivalent logarithm being subtracted from the register holding the logarithm value. Depending on whether the process is to find a logarithm or to calculate a power, there may be additional shifts, additional calculations, or both, to place the system in a converged configuration. At the end of some N iterations, the desired value of either the power or the logarithm is found to within a range of precision that depends on the number of iterations N and the number of bits used for the representation.

[0254] IV. Logarithm Lookup Table

[0255] The proposed method is that these logarithms are now JPEG0007599472000397.jpg1124 powers must be reflected. This means that the real part of the logarithm table must include JPEG0007599472000398.jpg1137 and the imaginary part JPEG0007599472000399.jpg1143 or JPEG0007599472000400.jpg1047 and need to be preloaded. JPEG0007599472000401.jpg829 is the desired fixed point format.

[0256] The base of the logarithm is 2, We chose JPEG0007599472000402.jpg1113 and JPEG0007599472000402.jpg1113, but it is unlikely that the first iteration of the algorithm for the logarithmic system will be different when the choice of base is changed. These values ​​were chosen to represent a design in which both processes can be achieved using the same hardware. If the logarithmic mode uses a different implementation with different hardware than the power mode, this may not be the best choice of base.

[0257] V. Logarithm Decimal Method

[0258] The core methodology for the new algorithm largely follows that used by the BKM process, with some new key differences. New comparison thresholds have been derived to accommodate the new bases for logarithms and powers, and the change in the convergence region required for the new method. There is also a slightly modified comparison threshold for the first iteration of the logarithm process to achieve convergence within the required region.

[0259] VI. Power Mode

[0260] The power mode described in this section is also modified to provide for complex multiplication with the power value. If this is to be achieved, and the output is to be correct, this must be preloaded before the range reduction step. This replaces the power value of the output, and should not be used if this value is required. It is also preferable to wait to store and apply the answer from the integer part of the logarithm (value zinteger,output) until the end of the process. This can reduce the storage required in intermediate registers where the processing takes place, but this is to be weighed against the additional storage requirements needed for the integer part of the solution.

[0261] Alternatively, the multiplication with the final power value can be accomplished in parallel by creating additional registers and performing equivalent operations in these additional registers. In this manner, the exponentiation process can complex multiply the output power value by nearly arbitrarily many other complex values.

[0262] Assuming the input is a fractional part of the logarithm, the convergence area algorithm JPEG0007599472000403.jpg1189 becomes:

[0263] 1. Assume there are four base registers, labeled JPEG0007599472000404.jpg1184, along with two additional slave multiplication registers to demonstrate the operation of the method when used for auxiliary complex multiplication. JPEG0007599472000405.jpg1221 and JPEG0007599472000406.jpg1220. The initial values ​​of these registers are: JPEG0007599472000407.jpg4899 Now, if there is no requirement for premultiplication, JPEG0007599472000408.jpg837. The slave multiply registers can be similarly constructed with: JPEG0007599472000409.jpg23101

[0264] 2. Index Value as JPEG0007599472000410.jpg67 Repeat JPEG0007599472000411.jpg819:

[0265] 3. JPEG0007599472000412.jpg918 is truncated to 3 bits creates JPEG0007599472000413.jpg926; the sign bit of 1 and the two integer bits of two's complement are The smallest change in JPEG0007599472000414.jpg929 is 1.

[0266] 4. JPEG0007599472000415.jpg1017 is truncated, We form JPEG0007599472000416.jpg1025, which has 3 bits: 1 sign bit, 1 integer bit, and 1 fractional bit in 2's complement, with a range of JPEG0007599472000417.jpg930, the smallest change is 0.5.

[0267] 5. Test the 3 bit values ​​to determine JPEG0007599472000418.jpg88: JPEG0007599472000419.jpg49122

[0268] 6. Apply multiplication to the power register using shifts and arithmetic operations: JPEG0007599472000420.jpg1058JPEG0007599472000421.jpg49129and: JPEG0007599472000422.jpg1058JPEG0007599472000423.jpg50128To apply the multiplication to these as well, use JPEG0007599472000424.jpg1020 or Do the same for the auxiliary registers such as JPEG0007599472000425.jpg1019.

[0269] 7. Subtract the corresponding entry in the logarithm table from the register: JPEG0007599472000426.jpg20102JPEG0007599472000427.jpg22110

[0270] 8. Shift the logarithm register one place to the left, doubling its contents in value. JPEG0007599472000428.jpg2268

[0271] 9. Go back to step 2 for the next iteration until N is reached, at which point the registers contain their final values. JPEG0007599472000429.jpg26155

[0272] and JPEG0007599472000430.jpg25159

[0273] Now that we understand the format of the process, While it is easy to find other test procedures that sometimes converge in the required region, by forming JPEG0007599472000431.jpg1151 or both, we have sought to reduce complexity by specifying the required value tests in their simplest well-known form, using either different numbers of bits or different comparison values.

[0274] VII. Logarithmic Mode

[0275] The logarithmic modes described in this section may also be modified to provide complex division with the input value. If this is to be achieved, it can only occur in parallel by creating additional registers and having equivalent operations performed in these additional registers. This is in contrast to auxiliary multiplication, which may overload the original registers, and cannot be achieved here because changing the power register would prevent the algorithm from converging in this mode. However, the use of auxiliary registers allows this to be avoided, so that, if desired, the logarithmic process can perform complex division of almost any other number of complex values, with the input value to this process as the denominator.

[0276] 1. If the output is the decimal part of the logarithm, the convergence area is Here is the algorithm for JPEG0007599472000432.jpg12100:

[0277] 2. JPEG0007599472000433.jpg1053 and Assume we have four base registers, labeled JPEG0007599472000434.jpg1019, along with two additional slave divide registers to demonstrate the behavior of the method when used for auxiliary complex division. There is JPEG0007599472000435.jpg1150. The initial values ​​of these registers are: The JPEG0007599472000436.jpg5792 slave divide registers are configured similarly: JPEG0007599472000437.jpg2597 Register Note that neither the doubling nor the -1.0 is applied to either JPEG0007599472000438.jpg1251.

[0278] 2. Index Value as JPEG0007599472000439.jpg67 Iterate over JPEG0007599472000440.jpg819:

[0279] 3. JPEG0007599472000441.jpg1018 cut off JPEG0007599472000442.jpg1029 is formed as 5 bits, 1 sign bit, 1 integer bit, and 3 fractional bits in two's complement, with the range JPEG0007599472000443.jpg930, the smallest change is 0.125.

[0280] 4. JPEG0007599472000444.jpg1019 is truncated, JPEG0007599472000445.jpg10291 sign bit, 2 integer bits, 2 fractional bits in 2's complement form to 5 bits; range is JPEG0007599472000446.jpg929, the smallest change is 0.25.

[0281] 5. Test two 5-bit values ​​(fixing the comparison in the first iteration), Determine JPEG0007599472000447.jpg88: JPEG0007599472000448.jpg68139

[0282] 6. Apply multiplication to the power register using shifts and arithmetic operations: JPEG0007599472000449.jpg1158JPEG0007599472000450.jpg50129and: JPEG0007599472000451.jpg1158JPEG0007599472000452.jpg49129

[0283] Do the same for the auxiliary registers JPEG0007599472000453.jpg1050 etc. and apply the division process to these.

[0284] 7. JPEG0007599472000454.jpg88, Add to JPEG0007599472000455.jpg1047: JPEG0007599472000456.jpg2370This is, Do not add to JPEG0007599472000457.jpg1145.

[0285] 8. Shift the power register one place to the left and double its contents. JPEG0007599472000458.jpg2367This is, Do not double the value of JPEG0007599472000459.jpg1044.

[0286] 9. Subtract the corresponding entry in the logarithm table from the register: JPEG0007599472000460.jpg47112

[0287] For the next iteration, Go back to step 2 until we reach JPEG0007599472000461.jpg77, at which point the D registers will contain their final values: JPEG0007599472000462.jpg2498

[0288] and: JPEG0007599472000463.jpg24125

[0289] In some cases, the form of the process can be understood as While it is possible to find other test procedures that converge in the required region by forming JPEG0007599472000464.jpg1155 or both, or by using a different number of bits or different comparison values, we have sought to reduce complexity by specifying the required value tests in their simplest known form.

[0290] VIII. Redundant Number System

[0291] Algorithms such as the BKM method are ideal to support techniques that slow down the local flow of data, allowing for integration into high-speed logic implemented in hardware. The use of redundant number system methods such as carry-save for the arithmetic steps, coupled with the presented reduction in the number of bits for comparison at each iteration, allows algorithms such as those presented here to have high-speed low-area implementations that are ideal for incorporation into integrated circuit designs. Advantages of such implementations include the ability to achieve higher clock speeds and accommodate operands with large numbers of bits, making this design ideal for creating low-cost, modern solutions to problems involving the manipulation of complex values.

[0292] For ease of intuition, rounding has been ignored in the preceding description of the algorithm and should not be construed as meaning that correct rounding is not considered in the present invention. Utilizing the most significant bits that are not directly used in any arithmetic operation to achieve results with reduced error through the application of rounding as applied at each step may increase the accuracy of the results described herein.

[0293] IX. Higher-order bottom approach

[0294] In a similar manner to other algorithms that evaluate bitwise, steps can be merged to create a way to calculate multiple bits at once. This is called a higher base algorithm because the output is treated as having a base different from 2, i.e. one digit from a new higher base is calculated each step. The drawback of this is that the lookup tables and shift choices must also be merged, creating a larger set of choices for the algorithm to work with each iteration.

[0295] X. Real Number Algorithms

[0296] Many systems can benefit from using only the real parts of the algorithms described herein. For example, the sigmoid logistic function, often used by neural network algorithms, can be generated by first determining the real number 2x, or using a powering process, if pre-multiplied by ex. It is then possible to add 1 and divide a preset constant by this value, or, if one continues to use the E-mode followed by the L-mode described in the original paper, the logistic function becomes: JPEG0007599472000465.jpg1647

[0297] However, the numerator k and the other modifiers a and b remain configurable.

[0298] XI. Demonstration of convergence

[0299] In FIG. 35, the convergence of the power law and logarithmic law is illustrated by a graph 3500. FIG. 35 illustrates how testing each possible value in the necessarily convergent domain results in the worst case of the method graphed here. JPEG0007599472000466.jpg shows that the multiplication results in 912 errors. The graph summarizes the worst-case behavior for over one trillion complex test values ​​(each real and imaginary value in the remapped fractional domain is tested for the number of iterations shown). To make the graphs comparable, the logarithmic results were converted to power form before testing for error. This ensures that the method is convergent for this large number of significant bits. This error can be further reduced by adding rounding to the multiplication. The graph shows that as precision and iterations increase, the error decreases as well, indicating convergence.

[0300] This graph 3500 shows the worst case absolute error from all possible cases within the convergence region shrinking at a similar rate as the number of bits in the table and the number of iterations in the process. This makes clear that for all possible correct inputs for this number of bits of fixed-point data, the algorithm successfully converges linearly to the correct solution, with the exact number of bits.

[0301] XII. Possible claims 1. A method for simultaneously converting input "logarithmic" complex numbers to resulting "power" complex numbers, comprising the steps of: a. the resulting complex coefficients comprise the squares of the real input complex parts to the power of the real part; and b the resulting complex argument constitutes a square of an integer power rotation multiplied by the real value of said imaginary input complex part; A method comprising: 2. A method for simultaneously converting input "power" complex numbers into resulting "log" complex numbers, comprising the steps of: a. the real binary logarithms of the real coefficients of the input complex numbers form the real part of the resulting complex number; and b. a real number representing the square of an integer power rotation expressed in said input complex number argument forming an imaginary part of a resulting complex number; A method comprising: 3. The method of claim 1 or 2, wherein the conversion is substantially accomplished through bit shifts, comparisons and calculations on a binary representation of the complex number. 4. The method of claim 1, wherein the reduction of the input and output ranges is performed by a value manipulation that substantially corresponds to rounding components of the input complex number separately from other parts. 5. The method of claim 2, wherein the reduction of the input and output ranges is performed by a value manipulation that substantially corresponds to rounding components of the output complex number separately from other portions. 6. The method of claim 1, wherein the output complex values ​​are initialized to values ​​not equal to 1, and the application of the iterations modifies the non-1 values ​​to calculate complex multiplications of the values ​​originally intended to be output. 7. The method of claim 1, further comprising: generating complex values ​​and wherein the operations applied to the output complex values ​​are also applied to these further complex values ​​to provide multiplication of these further complex values ​​by the complex value originally intended to be output. 8. A method as claimed in claim 2, in which operations which generate further complex values ​​and which are applied to the input complex values ​​are also applied to these further values ​​to give divisions of these further complex values ​​by the complex valued input. 9. The method of claim 1, further comprising implementing the process of repeating the calculations on the remaining complex values ​​using a redundant number system. 10. The method of claim 2, further comprising implementing the process of repeating the calculations on the remaining complex values ​​using a redundant number system. 11. A device designed to carry out the method of any one of claims 1 to 10. 12. A data processor configured to carry out the method of any one of claims 1 to 10.

[0302] 7. Acoustic fields constructed from KRYLOV subspaces of a second-order problem.

[0303] This disclosure describes a method for manipulating quadratic problems when represented as eigensystems to generate a reduced-rank basis set for inclusion in a linear system. By using a reduced-rank basis, we automatically derive an optimal information-preserving set of (potentially reduced) basis functions for inclusion in a complex-valued linear system, and find optimal complex-valued reweighting of such basis functions to achieve a goal, as described by several test functions.

[0304] The derivation of the reduced set of basis functions and how they are included at appropriate stages in the implementation of the sound field linear optimization problem are described.

[0305] I. Linear Control Region

[0306] Linear control region, which is a complex weighting function, Using JPEG0007599472000467.jpg839, we define an operator that performs an integral with a linear acoustic quantity such as pressure: JPEG0007599472000468.jpg1464This is used to estimate the sound field. Note that if JPEG0007599472000469.jpg814 is paired with a similar phase distribution of the acoustic quantities, then the complex conjugate is used to maximize the function evaluation. Any approximation to this integral (such as sampling) can be used, or in fact it can be used to approximate any linear operator over the sound pressure field. If it is even more linear, then a further integral can be constructed to optimize the control region using vector particle velocities (each vector component is complex-valued due to the harmonic motion of the field), where: JPEG0007599472000470.jpg1566, where: JPEG0007599472000471.jpg77 is the particle velocity of the medium, JPEG0007599472000472.jpg66 is a three-dimensional vector field with each dimension being a complex number, which acts as a weighting function for the same purpose as above. In this way, different acoustic volume patterns may be evaluated linearly, where the phase is When aligned with a "template" generated by a scalar or vector function in JPEG0007599472000473.jpg66, the resulting particular topology maximizes the objective function.

[0307] Also, Note also that JPEG0007599472000474.jpg814 can be mostly zero (or made up of zero vectors) to create a single estimate of multiple seemingly disjoint regions.This formulation can be extended to any linear operator over a field of linear acoustic quantities, so any integer approximation (such as sampling) can be used.

[0308] II. Quadratic and eigenvalue problems

[0309] A phased array can be described by a set of complex-valued coefficients that control each transducer: To maximize the amount of audio encoded in JPEG0007599472000475.jpg78, we use the following input: Driving an array of acoustic transducers at a fixed frequency with a 66 vector can be expressed as a linear algebra problem with the following intent: JPEG0007599472000477.jpg2049

[0310] This is the eigenvector, here This is also the definition of JPEG0007599472000478.jpg66, so this is a matrix It can be solved by taking the principal eigenvector of the matrix Constructing the JPEG0007599472000480.jpg78 can be accomplished by first considering the linear acoustic quantities emitted by each transducer as basis functions in the global space of acoustic quantities of interest. This is done by JPEG0007599472000481.jpg76 describes a method for achieving maximization of pressure and other acoustic quantities at a single point. However, the use of the described control area integrals allows for dominating the acoustic quantity emitted into the field of each individual transducer element. JPEG0007599472000482.jpg812 is JPEG0007599472000483.jpg813, resulting in a column vector: JPEG0007599472000484.jpg2029

[0311] By invoking linearity, we can multiply the complex activation values ​​that drive each transducer to derive a vector: JPEG0007599472000485.jpg1939, where n is the number of transducers. Assuming that this m is the result of a linear operator acting on a linear field, it can be Hermitian transposed and multiplied by itself to produce a quadratic form. Finally, we get JPEG0007599472000486.jpg957

[0312] Weighting Function Modify JPEG0007599472000487.jpg814 and add the matrix By solving the eigenvectors of JPEG0007599472000488.jpg78, a variety of different beams, surfaces, and other acoustic phenomena can be generated and maximized to create optimally operating regions within the sound field.

[0313] The null space generation technique described above ("Eigensystem Solution for a Quadratic Problem in Aerial Haptic Systems") is also used to extend this method to describe acoustic levitation, as well as to enable null steering and mid-formation, by introducing null points into the field at positions to describe the structure, allowing null steering and the formation of aerial textures. Note that this process increases the rank of the matrix, so simplifications cannot be used.

[0314] Similarly, multiple optimizations can be performed for some control regions. This can be done by adding together 78 matrices, which may also include matrices representing control points, producing matrices with increased rank, which is not possible in the case used by the simplification method.

[0315] III. Reduced-Rank Approximation for High-Rank Matrices

[0316] The key point is that high-rank matrices can be purposefully created so that approximations of lower rank are produced. For example, a system with hundreds of independent points can be automatically reorganized to create a realization of a basis set that appears to be higher dimensional than the original one, where very few are described by independent basis functions, increasing fidelity, reducing complexity, and allowing basis functions to be created automatically as needed.

[0317] Existing approaches generate basis functions for each control point or control region. By creating an eigensystem to optimize each control point or region and summing them or adding nulls, the eigensystems gain rank. Through an iterative process of adding basis function matrices and / or nulls, a high degree of rank can be achieved by the final matrix. Linear systems, due to resource constraints, can be at most JPEG0007599472000490.jpg Using 78 test functions, up to JPEG0007599472000491.jpg may only be able to process 66 basis functions. The high rank of this optimization matrix in this case JPEG0007599472000492.jpg87, where JPEG0007599472000493.jpg77 optimization matrices were merged additionally. Using this approach, JPEG0007599472000494.jpg718 and / or The other impossible situation, that is, JPEG0007599472000495.jpg819, is solvable by this system with finite processing power. JPEG0007599472000496.jpgThis reduces the problem to a reduced basis set with 66 basis functions.

[0318] This reduced basis set is a final matrix iteration consisting of multiple optimizations of the basis functions. The Krylov subspace is found by computing the 66 most dominant eigenvectors and a Krylov subspace that occupies substantially the same matrix subspace. The Krylov subspace is found using standard techniques, such as Arnoldi iteration, in which the dominant eigenvectors are found and then removed from the remaining matrix, and then This is repeated until 66 such eigenvectors are found.

[0319] IV. Mixed approach

[0320] This approach, where basis functions are extracted from the top matrix, can be mixed and matched with existing specific basis function approaches. For example, two control regions or points can be specified as mandatory and added to the basis set input to the linear system. Then, given that they are identified by a quadratic optimization, and thus an optimization matrix of their own, the remaining basis functions can be added to a single matrix that contains the most dominant ones that make up the remaining number of basis functions in the set input to the linear system. To find the eigenvectors of JPEG0007599472000499.jpg715, the Krylov subspace method is used to decompose it. In this way, the first two basis functions specified are guaranteed to be directly imported into the linear system, while the remaining functions in the set are approximated. This system is further used, where a set of basis functions is identified from the disjoint sets of basis functions, each described by a top-ranked optimization matrix, and the set of basis functions is taken. A further technique that can be used to derive reduced basis functions is to continue taking basis functions from a matrix until a power threshold is exceeded. This threshold can be different for each matrix composed of combined basis functions, thus generating eigenvectors of different power, in this case the eigenvectors can be generated by selecting the eigenvectors, the power of each eigenvector is calculated from each matrix based on the product of the actual priority value and the magnitude of the corresponding eigenvalue.

[0321] 8. Linear system solution using test functions

[0322] It would be valuable to have a way to "test" whether a basis function set achieves certain goals, defined separately from the basis function set, and a quantitative assessment of whether these goals have been achieved could be provided to an optimization that is used to produce the desired characteristics in the sound field defined by these test functions.

[0323] Describe test functions and their role in complex-valued linear system optimization.

[0324] I. Linear Control Region

[0325] Linear control region, which is a complex-valued weighting function We define it as an operator that performs an integral with a linear acoustic quantity such as pressure, using JPEG0007599472000500.jpg838: JPEG0007599472000501.jpg1466This is used for evaluation in the acoustic field. Note that if JPEG0007599472000502.jpg814 is paired with a similar phase distribution of acoustic quantities, then the complex conjugate is used to maximize the function evaluation. Any approximation to this integral (such as sampling) can be used, or in fact it can be used to approximate any linear operator over the sound pressure field. A vector particle velocity (where each vector component is a complex value for the harmonic motion of the field) may be further integrated to optimize the control region, where: JPEG0007599472000503.jpg1466JPEG0007599472000504.jpg77 is the particle velocity of the medium, JPEG0007599472000505.jpg66 is a three-dimensional vector field with each dimension being complex-valued, and serves as a weighting function for the same purpose as above. In this method, different acoustic volume patterns may be evaluated linearly, where the phase is When aligned with a "template" generated by a scalar or vector function in JPEG0007599472000506.jpg66, the resulting particular topology maximizes the objective function.

[0326] Also, Note also that JPEG0007599472000507.jpg814 can be mostly zero (or made up of zero vectors) to create a single estimate of multiple seemingly disjoint regions.This formulation can be extended to any linear operator over a field of linear acoustic quantities, so any integer approximation (e.g., sampling) can be used.

[0327] II. Basis Function Set

[0328] As input JPEG0007599472000508.jpg Given a set of 66 basis functions, JPEG0007599472000509.jpgTo find the optimal linear combination of 66 basis functions, we use a linear system JPEG0007599472000510.jpg78 test functions were used, JPEG0007599472000511.jpg to achieve the desired objective evaluation of the 78 test functions. The 78 test functions can take the form of linear control region evaluations. Using the format JPEG0007599472000513.jpg76, we get: JPEG0007599472000514.jpg2040Or it can be evaluated through a single control point definition: JPEG0007599472000515.jpg2060Here, the overline represents the complex conjugate, JPEG0007599472000516.jpg76 represents a three-dimensional position, JPEG0007599472000517.jpg76 parameterizes all the transducer elements, JPEG0007599472000518.jpg66 is the transducer activation multiplied by the conjugate of the transducer pressure field. If the normal direction is defined, the equivalent form of the control point defined by the particle velocity can also be used as an evaluator: JPEG0007599472000519.jpg2092 but with the transducer facing the wrong direction JPEG0007599472000520.jpg76, In the case of JPEG0007599472000521.jpg1037, its contribution can be set to zero and omitted.

[0329] In both single control point definitions, the existing approach can be viewed as a special single-point integral evaluation with a unit weighting function case for the control region evaluation, where the objective function is a simple acoustic quantity without extra side effects.

[0330] III. Finding the principal eigenvector

[0331] JPEG0007599472000522.jpg Using each of the 66 test functions, JPEG0007599472000523.jpgEvaluating each of the 78 basis functions results in a rectangular matrix that is complex valued. As mentioned before, the principal eigenvector is JPEG0007599472000524.jpgBest spans the space defined by 66 test functions JPEG0007599472000525.jpg represents a combination of 78 basis functions. However, this is not necessarily a square matrix, so a singular value decomposition is more appropriate, but the only eigenvectors needed are the dominant ones. A square Hermitian matrix can be constructed whose dominant eigenvectors represent the elements of the dominant combination of the system basis functions. The eigensystem is now written as: JPEG0007599472000526.jpg39138 where, JPEG0007599472000527.jpg821 has a phase JPEG0007599472000528.jpg The coefficients that produce the coverage of the parameter space of the 78 test functions are described. Thus, the dominant eigenvector is the component that represents the best coverage of the space achievable with a single combination vector. JPEG0007599472000529.jpg821 and its phase optimizes this cover.

[0332] Since the number of test functions is equal to the number of basis functions, It should be noted that even if JPEG0007599472000530.jpg719, the existing eigensystem approach may not necessarily be applicable, since the matrix may no longer be an Hermitian, and therefore applying the Rayleigh-Ritz method directly to this matrix without completing the steps above may prevent the method from converging to the principal eigenvector.

[0333] The principal eigenvectors should be used for phase information even when phase data is not specified or is incomplete. The eigenvectors are modified at the end of each iteration to remove amplitude data, phase averaged to zero, and phase progression limited at another time step from the previous solution, or have a set of phase solutions for some basis functions with an existing stationary phase relationship, and the complements of this set are optimized in phase between themselves and the existing stationary set of phases. In this way, this rectangular matrix system is used to achieve similar effects to the existing specialized eigensystem approaches for the condition that the function generating the basis set is the same function as the function used to evaluate the basis set objective.

[0334] IV. Linear System Solutions

[0335] The complex-valued linear system has a modified role: to find a set of basis functions that best reproduces a given value of the test function. Here, the linear system takes the form of a rectangular matrix: JPEG0007599472000531.jpg1978

[0336] Because the matrix is ​​rectangular, the number of test functions may be greater than the number of basis functions, resulting in a complex-valued least-squares system, or the number of test functions may be less than the number of basis functions, resulting in a complex-valued minimum-norm system. In either case, a standard decomposition can be used to find a solution vector, where each component of the solution vector is a coefficient applied to each of the basis functions that gives rise to the acoustic field when expressing a given value of the test function.

[0337] These test functions can represent the sound pressure or directional particle velocity at a point, or a mixture of pressure and directional particle velocity over a region in space. Any combination of these functions JPEG0007599472000532.jpg78 test functions can be constructed, which are solved in the final step as a linear system.

[0338] Once the coefficients of the basis functions are known, the final values ​​for actuating the transducers are expanded as follows: JPEG0007599472000533.jpg2145

[0339] An input signal that gives the assumption of a single monotonic carrier frequency (which may be more or less valid under changing circumstances) is It is constructed by applying the components of JPEG0007599472000534.jpg66 to each transformer.

[0340] 9. conclusion

[0341] Various features of the above-described embodiments may be selected and combined to produce numerous variations of improved haptic systems.

[0342] In the foregoing specification, particular embodiments have been described. However, those skilled in the art will recognize that various modifications and changes can be made without departing from the scope of the invention as set forth in the claims. Accordingly, the specification and drawings are to be regarded in an illustrative rather than a restrictive sense, and all such modifications are intended to be included within the scope of the present teachings.

[0343] Benefits, advantages, solutions to problems, and any elements that may cause or make more pronounced any benefit, advantage, or solution are not to be construed as critical, necessary, or essential features or elements of any or all claims. The present invention is defined solely by the appended claims, including any amendments made during the pendency of this application, and all equivalents of the claims as issued.

[0344] Furthermore, relational terms such as first and second, above and below, etc. are used herein only to distinguish one entity or operation from another entity or operation without necessarily requiring or implying any actual such relation or order between such entities or operations. The terms "comprise", "comprising", "having", "having", "including", "including", "comprising", "including" or any other variation thereof are intended to cover a non-exclusive inclusion, such that a process, method, article, or apparatus that comprises, has, includes, or includes a list of elements includes not only those elements, but also other elements not expressly listed or inherent in such process, method, article, or apparatus that are expressly listed. Elements such as "comprise", "having", "including", "including" do not exclude the presence of additional identical elements in the process, method, article, or apparatus that comprises, has, includes, or includes the element. The terms "a" and "an" are defined herein as one or more, unless expressly stated otherwise. The terms "substantially," "essentially," "approximately," "about," or any other variants thereof, are defined as approximating those understood by one of ordinary skill in the art. The term "coupled," as used herein, is defined as connected, although not necessarily directly, and not necessarily mechanically. A device or structure that is "configured" in a particular way is at least so configured, but may also be configured in other ways not described.

[0345] The Abstract of the Disclosure is provided to enable the reader to quickly ascertain the nature of the technical disclosure, and is submitted with the understanding that it will not be used to interpret or limit the scope or meaning of the claims. In addition, it can be appreciated that in the foregoing Detailed Description, various features have been grouped together in various embodiments for the purpose of streamlining the disclosure. This method of disclosure should not be interpreted as reflecting an intention that the claimed embodiments require more features than are expressly recited in each claim. Instead, as the appended claims reflect, inventive subject matter lies in less than all features of a single disclosed embodiment. Accordingly, the appended claims are hereby incorporated into the Detailed Description, with each claim standing on its own as separately claimed subject matter.

Claims

1. i) generating an acoustic field from a phased transducer array having known relative positions and orientations; ii) defining a plurality of control points within the acoustic field, each of the control points having a known spatial relationship relative to the phased transducer array; iii) positioning said plurality of control points in a physical field of quantities generated by said phased transducer array by superposition of wave phases and amplitudes resulting from input drive coefficients given by an optimization process, said quantities being acoustic vector quantities; iV) optimizing to maximize the absolute magnitude of field values ​​along known directions of said control points in a physical field of quantities generated by said phased transducer array by superposition of wave phases and amplitudes resulting from input drive coefficients given by an optimization process, said quantities being acoustic vector quantities; The method includes:

2. The method of claim 1 , wherein at least one of the plurality of control points comprises a velocity of a particle at a point within a medium hosting an acoustic wave.

3. The method of claim 1 , wherein at least one of the plurality of control points comprises momentum of an acoustic wave.

4. i) generating an acoustic field from a phased transducer array having known relative positions and orientations; ii) defining a plurality of control points within the acoustic field, each of the control points having a known spatial relationship relative to the phased transducer array; iii) optimizing to maximize the absolute magnitude of field values ​​along known directions of the control points in a physical field of quantities generated by the phased transducer array by superposition of wave phases and amplitudes resulting from input drive coefficients given by an optimization process, the quantities being acoustic vector quantities; The method includes:

5. The method of claim 4 , wherein at least one of the values ​​of the plurality of control points comprises an acoustic intensity.

6. The method of claim 4 , wherein one of the plurality of control points comprises an acoustic energy flux vector.

7. The method of claim 4 , wherein at least one of the plurality of control points comprises a velocity of a particle at a point within a medium hosting an acoustic wave.

8. The method of claim 4 , wherein at least one of the plurality of control points comprises momentum in an acoustic wave.

Images