Driving techniques for phased array systems

文档序号:1590942 发布日期:2020-01-03 浏览:27次 中文

阅读说明:本技术 相控阵列系统的驱动技术 (Driving techniques for phased array systems ) 是由 本杰明·约翰·奥利弗·朗 菲利普·埃德温·琼斯 于 2017-12-13 设计创作,主要内容包括:描述了用于驱动相控阵列系统的各种技术,特别是旨在用于应用于半空触觉、参量音频、声悬浮和声成像的声学相控阵列,包括一种系统:1)能够减轻空中变化的影响以提供一致的触觉体验;2)在空中产生陷阱点;3)根据用于产生更一致的触觉效果的矢量来定义相控阵列优化;4)经由受控的声场在空间中定义一个或多个控制点或控制区域;5)使用简化的表示方法来构造声学基函数;6)针对大量的吞吐量执行复值的高效评估;7)生成矩阵的Krylov子空间;以及8)使由不同控制点和/或区域描述的目标最大化到被用于创建声学基函数的那些目标。(Various techniques are described for driving a phased array system, particularly an acoustic phased array intended for application to semi-hollow haptics, parametric audio, acoustic levitation and acoustic imaging, including a system: 1) the impact of airborne variation can be mitigated to provide a consistent haptic experience; 2) creating trap sites in the air; 3) defining a phased array optimization according to vectors for producing more consistent haptic effects; 4) defining one or more control points or control regions in space via a controlled sound field; 5) constructing acoustic basis functions using a simplified representation; 6) performing an efficient evaluation of complex values for a large amount of throughput; 7) generating a Krylov subspace of the matrix; and 8) maximizing the objectives described by the different control points and/or regions to those used to create the acoustic basis functions.)

1. A method, comprising:

creating control points in the sound field;

assigning a normal vector and direction to each control point to incorporate the amount of fluctuation into the phased array solver;

solving for the x, y, and z components of the amount of fluctuation at the point of interest.

2. The method of claim 1, wherein each is designated as ax,j0q、αy,y0qAnd alphaz,j0qIs defined as the x, y and z components of

Figure FDA0002163150250000011

Figure FDA0002163150250000013

Wherein the content of the first and second substances,

Figure FDA0002163150250000014

Technical Field

The present disclosure relates generally to improved techniques for driving phased array systems, primarily for generating designed sound fields for creating half-space haptics, parametric audio, acoustic levitation, and acoustic imaging.

Background

The continuous distribution of acoustic energy (which we refer to as the "sound field") can be used for a range of applications including mid-air haptic feedback.

In the sound field, one or more control points may be defined. These control points may be amplitude modulated with a signal and thus generate vibrotactile feedback in mid-air. An alternative way to generate feedback is to create control points that are not modulated in amplitude and move them around spatially to create a spatio-temporal modulation that can be perceived.

The minimum size and maximum power of the control points created by the device depends on the array layout, geometry and relative position of the one or more control points to be created. Thus, it may not be possible to match the behavior of a control point in one location to the behavior in another location. Additional complications may arise when dynamically moving control points through space. The resulting haptic functional change of the point experienced by the user is undesirable and often counter-intuitive, resulting in confusion and loss of experience immersion. Therefore, a system that can mitigate the effects of air changes and can reverse them to provide a consistent haptic experience is commercially valuable.

Furthermore, it has been shown that a general non-linear optimization system can be used to create trap points (trap points) in air. The mechanisms used for this purpose have a number of drawbacks, including but not limited to: potential fields (potential fields) suitable for levitation over a phased array system, use of unnecessarily expensive and time consuming optimization processes, and phase manipulation limitation are not evaluated.

Solving these problems is important because existing systems create traps that are not as strong as they could be and cannot scale properly due to these shortcomings. This allows the creation of systems that utilize acoustic levitation to provide additional visual information to facilitate alternative uses of multimodal methods and haptic systems.

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

The natural directional component of the acoustic phased array is of little concern due to the requirements of existing systems. For example, high frequency therapeutic medical systems that aim to focus energy into a small area are sufficiently approximated by scalar pressure (involving changes in the number of air molecules) to determine the best solution for the injected energy, and therefore do not require a vector solution. Since the focused phased array vector behavior (describing the direction of motion of air molecules) can be ignored in most cases, effectively achieving this has not been a priority.

However, creating haptic feedback in mid-air involves exchanging vibrations with human skin. This relies entirely on the surface in contact with the wave front and with the momentum carried by the wave, which is a result of the coupling with the vector behavior of the wave. Furthermore, the assumed linear relationship between momentum and pressure does not hold, as the system is focused such that the angle of incidence between the transducer elements is more and more different. It is not possible to ignore vector components in this system and thus scalar pressure is not sufficient to consistently find a suitable solution for creating a semi-null haptic effect. For this reason, it is of commercial interest to build this awareness of momentum effects into the device.

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

It has been shown that the computation of feature vectors maps to some very specific solutions of purely quadratic constraint problems in matrix algebra. If a matrix can be formed (solving its eigenvectors for the control coefficients of a given transducer for a particular problem in semi-null haptics), then this solution can form the basis for a simplified device implementation that presents the prior art and is therefore of value.

Furthermore, in order to describe the basis functions using a simplified representation method, it is first necessary to develop shorthand capable of simple (low entropy) description of the useful functions in terms of details of the sound field.

Ideally, the description would be something that is accessible from the perspective of the user.

After having found this available, low entropy representation, there must also be a basis for converting it into its equivalent transducer basis (a set of transducer bases that describes the basis of the transducer acoustic field for one of potentially many basis functions in later solutions). This would enable prospective users to design a variety of different structures within a sound field, enabling them to build custom systems to manipulate the energy flow in the field to the maximum, resulting in the creation of systems that can use many different effects (such as focal points, Bessel beams, local hollow beams (bowbeam), curved beams, etc.) at once. These will extend the concept of control points into control areas that encapsulate richer sound field control methods.

It has also been found that when considering momentum transfer in haptics, some problems have to be dealt with in order to use particle velocity as a substitute for wave direction.

Furthermore, performing efficient parallel complex-valued multiplication with the same second operand is necessary to achieve the large amount of throughput required for devices with performance metrics consistent with those outlined in the requirements document. In particular, for parallel computation, one element in the summation of the computations resulting in the rows of the necessary optimization matrix will suffice, but provided that the subsystem is replicated a sufficient number of times. Thus, in order for a subsystem to be useful, it must be small and efficient. This urgent need for efficiency has created a different view of the old algorithms, which has created a new approach to evaluating complex transcendental functions on top of simple multiplications.

In the wider world, this resulting algorithm is applicable anywhere the elements of the processing hardware must interact efficiently or in real-time with complex valued data whose modulus cannot be ignored. This includes, but is not limited to, communications including GPS, wireless, radio and satellite transmission and reception systems, audio processing, phased array systems including half-space haptics, time-of-flight imaging and digital interferometry, calculation of scattering effects and polarization, and the like. This may translate into a broader scope of opportunity for the present invention.

The algorithm of Volder, also known as CORDIC (Volder, j. -The CORDIC computing technique-1959), can multiply The complex unit values, but cannot include The modulus in The second operand of The computation without further computation beyond that of The basic CORDIC algorithm.

For this reason, the isolation algorithm invented in 1994, which was implemented by the authors a few years ago, is called BKM by the authors Bajard, Kla and Muller, and can also achieve the goal of complete complex-complex multiplication (Bajard, J., Kla, S. and Muller, J. -BKM: A new hardware algorithm for complex electronic functions-1994). However, this algorithm has a complex implementation and therefore does not receive as much attention if it is not complex. The algorithm may achieve this result by reducing complex-valued operands to their complex logarithms, adding them, and then exponentiating to achieve a full complex-to-complex multiplication.

In some cases, there is a more efficient way to implement such complex-valued multiplication and complex-valued division with BKM and the new algorithm, both serially and in parallel.

The 1999 complementary compilation of the BKM algorithm addresses the overly complex implementation of complex logarithms, somewhat simplifying it. (Bajard, J. and Imbert, L. -Evaluation of complex electronic functions: New version of BKM-1999). This BKM method in 1994 and the complement in 1999 form an algorithm that will be abbreviated as BKM in this document.

In attempting to optimize algorithms for incorporation into real hardware, the idea of new technology becomes apparent, which prioritizes easy integration into our or more generally any existing system, rather than evaluating the mathematical functions described. This new technology, when it is developed into proof-of-concept code, turns out to have many additional emerging properties that make it more promising for the general adoption of hardware implementations than the original algorithms. In addition, if the exact graph requested by the original technique is needed, further scaling can be applied to retrieve these results, where the algorithm as a whole is still competitive. In either case, the result is a technique with substantial benefits over the prior art.

It is known that quadratic problems with unit norm can be solved by constructing a matrix whose characteristic problem is equivalent to a quadratic problem. A quadratic problem matrix larger than rank-1 may be generated by summing many rank-1 problem matrices and/or expanding the null space of existing problem matrices. By generating the Krylov subspace of the matrix, an optimized subset of the implicit basis functions can be generated. By using the Krylov subspace as the basis set for a linear system, a correct or approximate evaluation of the basis set using a number of test functions can be achieved.

Depending on how many feature vectors are used and where they are selected, graceful degradation of a complex set of many basis functions can be achieved because the power of each subsequent feature vector is reduced.

Furthermore, the resulting complex-valued linear system for defining the shape of the sound field takes n basis functions and recombines these basis functions into a resulting field (restingfield), which describes the desired configuration of sound pressure at the carrier frequency.

Until this point, the control point basis functions were only used to actuate the control points. Described herein is how to generalize to apply different sets of control test functions to an initial control point basis function or an initial control domain basis function.

This is important because it breaks the link between the initial "palette" of behaviors in the sound field and the objective function used by the linear system. For example, this may be used to create a basis set of Bessel functions, which may then be used to actuate the control points.

Drawings

The accompanying figures, in which like reference numerals refer to identical or functionally-similar elements throughout the separate views and which together with the detailed description below are incorporated in and form part of the specification, serve to further illustrate embodiments that include the concepts of the claimed invention and to explain various principles and advantages of those embodiments.

Fig. 1 and 2 are arbitrary dynamic range haptic waveforms that can be processed by a nonlinear remapping function to create new waveforms of controlled dynamic range.

Fig. 3-12 are effective amplitude modulation waveforms generated by control points in different configurations to which different range control techniques are applied.

Fig. 13-16 are projections of the sound pressure resulting from optimization of Gor' kov correlation functions.

FIG. 17 is a diagram of the same wave with the same sound pressure level, but in a situation where the effect of the vector components may or may not produce the maximum haptic effect.

FIG. 18 is a simulation of measurements of the sound field generated by a phased array with dB SPL, dB SIL, and dB SVL.

FIG. 19 shows a control point diagram created by solving a quadratic problem.

Fig. 20 shows the intersection (cut) through the acoustic field at the edge of the near field region parallel to the flat linear transducer array.

FIG. 21 shows a graph of the Bessel function.

Fig. 22-23 show approximate cross-sections for higher order Bessel and Mathieu beams.

Fig. 24-26 show partial cross-sections of curved beams, local air core beams, and local air core beams.

Fig. 27-29 illustrate beams with increasing wavelengths, serpentine beams, and spiral beams.

Fig. 30-31 illustrate the effect of removing the transducer that contributes waves from the direction opposite to the direction of the desired waves.

Fig. 32-34 illustrate the region of convergence of the logarithmic and exponential methods without remapping.

Fig. 35 shows a convergence diagram of the logarithmic method and the exponential method.

Skilled artisans 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 improve understanding of embodiments of the present invention.

Apparatus and method components have been represented where appropriate by conventional symbols in the drawings, showing only those specific details that are pertinent to understanding the embodiments of the present invention so as not to obscure the disclosure with details that will be readily apparent to those of ordinary skill in the art having the benefit of the description herein.

Detailed Description

The foregoing techniques, devices and methods may be used in combination or individually in a semi-empty haptic system.

1.Dynamic range reduction for a semi-air haptic device

I.Psychophysical modeling

The behavior and haptic output of the system must be calibrated according to the way a person perceives touch. For example, in the case where the control point size is increased, the tactile sensation is also increased although the physical magnitude generated at the control point is not changed. Cross-modal effects (e.g., from audio) are used to enhance haptic feedback if they are correlated. This leads to differences between the physically defined control points, necessitating the use of psychophysical modeling steps in order to achieve equivalent perception.

II.Remapping using a psycho-tactile (psycho-tactile) model

Due to the physical limitations of the playback device, the precise waveform may not be reproducible as specified. The waveform can be constructed in a psychohaptically defined perception space. Or it may consist of intensity and force measurements on a substantially absolute physical scale, which may be a record of a haptic texture or arbitrarily constructed waveform potentially expressed with additional metadata. This waveform represents the behavior of the control point for a given "perfect" semi-empty haptic device. This is then perceptually evaluated using a psychophysical model that can describe and quantify the perceptual aspects of the re-creation of the control points. The control point may move while the haptic effect is recreated. The movement may be included as part of a perceptual evaluation of the haptic effect. This can then be re-interpreted as an output space available for an instance of the semi-empty haptic system constrained by user position, relative geometry, and available power limits, where the system is configured to create substantially equivalent control points.

Psycho-tactile models give the device the ability to determine how the modulation waveform of a control point can be modified in order to evoke substantially the same tactile response as a perfect device. The process of implementing this remapping may also include a psychoacoustic modeling step to prevent the resulting haptic waveform from becoming noisy, or to select a more acceptable modulation above audible that is closed (close) but less noisy in perceiving the haptic space. It is known that many waveforms can be made more haptically powerful by increasing the amount of audible noise. In some cases, this may be noisy enough to be difficult to maintain, while in other cases it may be desirable. The psychoacoustic model can then be used in conjunction to provide an audio to haptic compromise setting that biases the output to be quieter or stronger.

III.Psychological tactile unit

If a waveform is defined in the psycho-tactile space, where the amplitude is intended to be perceptually equivalent, a unit must be designed that takes into account the psycho-tactile effect of the generated control point in mid-air. Most notable of these is the integration of tactile sensations across the hand. The focal region involved by each control point changes size and shape again depending on the nature of the physical device and this must be counteracted in the construction of the psycho-tactile unit. Such a psycho-tactile unit changes what the magnitude of the physical control point it evaluates to depending on the size and shape of the focal region, thus counteracting the effect of integration of the tactile sensation over a large area. In this way, the method allows the same number of psycho-tactile units to encode equivalent sensations, keeping the haptic equivalence on the device and physical instance.

IV.Dynamic waveform remapping

To achieve an optimal mapping to the target device, the complete waveform and its trajectory may be considered. However, due to the interactivity and thus the requirement of low latency, only a portion of the waveform or even individual samples may be available. While the entire waveform can be remapped into the capabilities of the target device, doing so requires that the waveform must be known in advance, and it can therefore be considered a pre-processing step for the known waveform. An alternative is to remap the physics of each waveform sample in amplitude, since remapping amplitude samples can be done in real time with an interactive waveform, making it more suitable for providing haptic feedback in real time.

Various functions may be considered as good candidates for providing remapping for the control point amplitude samples. These include log-linear curves and sigmoidal functions (sigmoidal functions), which can convincingly reproduce the perceived sensory input on devices with low dynamic range in other modalities. Although these curves generally do not require specification of an upper limit in the mapping function, as shown in fig. 1 and 2, the dynamic range can be used more efficiently, provided that there is a limited amount of available power. FIG. 1 illustrates a haptic waveform of arbitrary dynamic range 100 that can be processed by a nonlinear remapping function 110 to create a new waveform of controlled dynamic range 120. This can be used to acquire and convert any dynamic range of haptic waveforms to a smaller range. Since the waveform can be infinite in amplitude, this does not fill the entire dynamic range of the output, resulting in an unused output range. In fig. 2, similar to fig. 1, a non-linear function 210 remaps the amplitude waveform 200 to an amplitude available to a target device 220. Awareness of the range of the wave allows the remapping function to be accurately adjusted to fit the output range specific to the device. This allows the full range of device outputs to be utilized.

To illustrate, a minimum linear rescaling (minimalist linear rescaling) is used in the following example to describe how mapping prevents undesirable behavior in a half-air haptic device.

A graph 300 is shown in fig. 3, where the first step towards achieving remapping of the waveform into a reproducible perception is to feed data into a solver function describing the waveform and its limits 300. In FIG. 3, the solid line is A'1And the dotted line is A'1,range. This is the input to the simple waveform remapping function. The y-axis here is expressed in units inside the device.

This is achieved by defining the maximum amplitude that needs to be reproduced to transmit the waveform without clipping-this is denoted as A 'at each sample time'1,range. Then it is taken as input and the input waveform amplitude A'1Fed together into the solver.

A graph 400 is shown in fig. 4, which shows the result of attempting to reproduce the waveform without applying any changes to the input data. The result is that the waveform is clipped due to the lack of available power to recreate the necessary amplitude at the control point. In FIG. 4, the solid line is A'1,simulatedIt is amplitude (A'1,χ1T) and attempts to directly output the results of the user-specified magnitudes in fig. 1 regardless of the dynamic range of the device.

Defining the physical configuration of the active transducers in the system as T and position χ1Desired control point amplitude A'1Function magnitude (A ') can be defined using the assistance capability of solver function of semi-air haptic system'1,χ1T), if a given device is physically capable of providing this at a given time, the function magnitude (A'1,χ1T) produces the desired carrier amplitude at the control point, or alternatively the closest amplitude that can be reproduced at the control point location.

Evaluation A1,limit:=amplitude(A′1,range,χ1T), the maximum value of the range to which the input amplitude can be mapped at that point in time is generated. This need not be accurate, since the control point does not move much between samples, so a from a previous iteration of the control point can be used1,limitThe value of (c). The value may also be filtered to remove rounding errors and other noise sources. In this case, the remapping function takes a single sample of the waveform and drives it from the interval 0, A'1,range) Remapping to interval [0, A1,limit) The above. Shown in FIG. 5 is a graph 500 that results from A 'through a remapping function 500'1∈[0,A′1,range) To A1∈[0,A1,limit) Is directly mapped. By varying the amplitude from A'1Remap as by A1,limit:=amplitude(A′1,range,χ1T) to a value of a1, clipping of the waveform is prevented and a'1,simulated=A1

V.Remapping multiple simultaneous control points

To use the system for multiple control points at a time, the ranges must be evaluated simultaneously so that the solver can guarantee that: the waveform can be represented correctly, as in the two sine wave case example of fig. 6. A graph 600 is shown with different amplitudes A 'at different points in time'1And A'2Two control points are created. With respect to the transducer array, one control point moves toward it and the other control point moves across it. They are arranged to oscillate in amplitude with a sine wave of two different frequencies.

Unlike a single control point (where violating the power limits of the device results in power starvation), clipping waveforms with multiple control points that exceed the limits of the device promote crosstalk between points, confounding their haptic effects. Shown in fig. 7 is a graph 700 illustrating an example of two movement control points. Fig. 7 illustrates an attempt to recreate two control points using the transducer column, which results in complex and unpredictable behavior because different resonances between the two points are created and replaced, but eventually fail to reach the desired level.

Shown in FIG. 8 is a graph 800 that uses an auxiliary solver evaluation to determine the power limit present at two points, which can then be written as { A }1,limit,A2,limit}:=amplitude(A′1,range,A′2,range,χ1T), resulting in the waveform shown. Here, one control point is moving towards the device with the potential to become stronger, and another control point is moving across the face of the device and away from the device with the potential to become weaker. As can be seen from fig. 8, their waveforms can then be represented as expected and without crosstalk due to the dynamic range limitations generated by the solver solution. Thus, fig. 8 shows two points whose amplitude envelope is reduced to a level that is guaranteed by the solver to be faithfully recreated without error or with a controlled known level of error. They are the same power because this is the configuration that minimizes the squared magnitude error at the idealized magnitude level.

VI.Power redistribution across multiple simultaneous control points

The amplitude range may be set to different values and may even be set dynamically over time in order to manage the power level. In this way, the steering range can determine the proportion of the array output that is concentrated at a given control point. This is illustrated in fig. 9 and 10 for simple manipulation of power from one control point to another, but is also a more complex out-of-phase relationship that enables behavior dynamics between control points such as those illustrated in fig. 11 and 12.

A graph 900 of input waveforms and ranges with two input control points is shown in fig. 9. These have different waveform ranges that are also doubled as hints to allow each point to draw more or less power available from the device. Shown in fig. 10 is a graph 1000 having an output derived from the inputs in fig. 9. In this configuration, the power generated by the array moves from one control point to another as the control points move away from and towards the array, respectively. Thus, more power may be used than if the control points were expected to have equal magnitudes, as expected for fig. 6, 7 and 8.

Shown in fig. 11 is a graph 1100 with an alternating function (alternating function) in which one control point or group of control points ramps up, where another control point or group of control points that gradually decreases in amplitude can be effectively modeled by matching the range to a sinusoid. The result is to reproduce each waveform according to the power requirements of each side of the modulation system while ensuring that each side of the system gets full array power. The result is time slicing (time slicing) without the potential for arbitrary output clipping. Shown in fig. 12 is a graph 1200 with the output resulting from the input in fig. 11 but without the described grouping. This results in a reduction of the output dynamic range, since the power budget is unnecessarily shared between all control points.

2.Efficient creation of acoustically suspended trap sites

I.Optimization target of suspended spherical particles

The optimization problem required for suspended particles is centered on the use of Gor' kov energy potential (energy potential). This is a scalar function that evaluates the energy at each possible location of the particle in the acoustic field. Finding a minimum or low energy state at a spatial location corresponds to finding a potentially stable trap point. However, there are acoustic effects that are not considered by the Gor 'kov potential (Gor' kov potential), such as acoustic streaming, where sound waves cause bulk motion (bulk motion) in a fluid. For this reason it is also preferable to ensure that traps do not appear in the region of fluctuating pressures, which are usually determined to be the smallest possible among the energy potentials, since these fluctuations are disruptive since other factors are not included in the definition of potential.

Since the direction and magnitude of external forces such as gravity and acoustic flow effects cannot be derived from Gor 'kov potential alone, it is generally preferred to compute the Laplacian of Gor' kov scalar potential (Laplacian). This gives a more stable arrangement of a degree of control over which axes are prioritized over others in order to generate internal thrust.

II.Evaluation of Gor' kov potential

The Gor' kov potential field is an energy field that enables spherical particles to float due to the presence of a minimum of energy, which is described as (in Bruus 2012):

Figure BDA0002163150260000121

Figure BDA0002163150260000122

Figure BDA0002163150260000123

wherein c is0Is the speed of sound through the acoustic medium, cpIs the speed of sound through the spherical particles, p0Is the density of the acoustic medium, ppIs the density of spherical particles, related to the time-averaged pressure

Figure BDA0002163150260000124

The left term of (a) is a representation of potential energy and relates to the particle velocity vector

Figure BDA0002163150260000125

The right hand terms of (a) are inversely related to potential energy. Since the Gor' kov potential makes the assumption of being present in the active or far field true, these pressure and velocity terms can be replaced with time harmonic terms that are in phase. This indicates that the time harmonic Gor' kov potential field can be rewritten as:

extending this to many sources yields:

Figure BDA0002163150260000127

Figure BDA0002163150260000131

the problem with this is that the cross term actually means that this becomes:

Figure BDA0002163150260000132

resulting in a term that is quadratic and grows in the number of squares of the transducer.

III.Optimizing potentials by feature problems

Ignoring time harmonic effects and considering pressure and particle velocity inputs in each transducer, the potential definition can be rewritten as:

Figure BDA0002163150260000134

U=xHMx,

where j is the index of the desired trap, and additionally

Figure BDA0002163150260000141

Representing the phase and amplitude changes as one moves from transducers q and r to trap point j. Using constraint xHx is 1 to ensure that the transducer power is at a non-zero level, the optimum complex transducer activation coefficient to achieve a potential well is given by the eigenvector of M with the largest eigenvalue.

The selection of x is such that the objective function xHMx maximization due to Mq,rThe negative sign in the definition yields the minimum in Gor' kov potential field. This is because xHx is 1 is the definition of the feature vector, and xHMx=xHX, x when x (eigenvalue of matrix) is maximized by definitionHThis must yield a maximum value if it exists when x is 1. The standard statement of the problem then becomes:

maximizing xHMx

Is subject to xHx=1,

And is

Figure BDA0002163150260000142

However, as previously discussed in the literature, minimizing the energy of the Gor' kov field creates a focal point rather than an actual trap point. To create trap sites, pressure fluctuations at the traps are minimized. This is important because the definition of pressure fluctuations is complex-valued and represents an amount of fluctuation that cannot be meaningfully minimized but instead made zero.

IV.Constructing nulls within a feature problem to generate pressure nulls

The zero pressure condition may be constructed as a series of constraints that evaluate to zero pressure at a given point in the sound field. If there is a portion of the solution vector that is parallel to the known best pressure focus solution in the high dimensional phase space, the pressure can only be evaluated as a non-zero value at a given location. For trap point j, there is a vector

Figure BDA0002163150260000143

It is as follows:

Figure BDA0002163150260000144

thus, the constraint of zero pressure at trap point j can be described as the constraint equation:

Figure BDA0002163150260000145

since what amplitude is not important for this and the constraint needs to be unit length, therefore:

Figure BDA0002163150260000151

thus, the standard statements of this problem are:

maximizing xHMx

Is subject to xHx=1,

Figure BDA0002163150260000152

And is

Figure BDA0002163150260000153

The eigensystem results were obtained from Golub (1973) and re-expressed as complex values, giving the following parameters: complex value

Given the above system with m 1, the characteristic equation for this system can be written as:

Figure BDA0002163150260000154

where λ and μ are Lagrange multipliers (Lagrange multipliers). This can be distinguished by the conjugate of the x vector, resulting in a system:

Mx-λx+μc=0,

by multiplying by cHGenerating:

cHMx-cHλx+cHμc=0,

cHMx-cHxλ+cHcμ=0,

assuming c is a unit length, c is generatedHc 1, earlier definition of the problem yields cHx is 0, and replacing both yields:

μ=-cHMx。

replacing this back into the earlier derivative yields:

Mx-cHMxc=λx,

it is as follows:

(I-ccH)Mx=λx,

since this is a statement of the feature system, however, (I-cc)H) M may not be symmetric, making feature vector computation difficult. However, I-ccHCan be shown as idempotent ((I-cc)H)2=I-ccH) And the eigenvalues are invariant through the ordering, yielding:

(I-ccH)M(I-ccH)x=λx,

extending this to more constraint vectors, further multiplications on the left side do not present a problem:

which after multiple multiplications ultimately yields:

Figure BDA0002163150260000162

C=[c1 … cj … cm],

wherein

Figure BDA0002163150260000163

Represents the Moore-Penrose pseudo-inverse.

In summary, by determining to have

Figure BDA0002163150260000164

The solution of the initial problem can be found by using the feature vector of the maximum feature value. Multiplication by

Figure BDA0002163150260000165

Can be viewed as adding a rank-1 variation, adding a base vector to the null space of the solution matrix.

Finally, it is useful to add support for the laplacian to compute Gor' kov potential in a simple manner. The Gor' kov potential solution found herein creates trap sites, but there is no control over which directions force is applied from, and there may even be no direction of the force. To achieve this, it is necessary to optimize for the point of maximum divergence in the gradient of the potential field and then re-weight the gradient of the field. This is known as the laplacian of Gor' kov field.

V.Modification of the Laplacian for generating Gor' kov fields

The modification of the laplacian to minimize Gor' kov field is only necessarily applicable to the M matrix. If the effects of wavefront curvature from each transducer are discarded and they are treated as plane waves only, the singular point (singularity) that complicates the formula is removed and the particle velocity becomes constant. The plane wave obtained by each transducer can be rewritten as follows:

Figure BDA0002163150260000171

then, without Gor' kov for the laplacian, the matrix elements can be represented as an equivalent system:

Figure BDA0002163150260000172

wherein M isq,rAll constants in are combined as Aj,q,r. Taking the Laplace operator of the formula, then generating a modified matrix formula:

Figure BDA0002163150260000173

using this value instead of the corresponding component of the M matrix yields a solution optimized for the laplacian of Gor' kov potential. Due to the item

Figure BDA0002163150260000174

Can be rewritten as:

Figure BDA0002163150260000175

it can therefore also be divided into contributions from the x, y and z directions:

in this way, weighting may be applied to generate the constraint axis and bias the forces generated from each direction.

VI.Implementation in a device

The system in this document represents an apparatus that can be used to compute one or more floating traps based on a simple feature system solver. This can potentially be implemented in device firmware to create embodiments that can use acoustic levitation for confinement, trapping, and transporting items.

This is illustrated in fig. 13-16.

Fig. 13 shows the laplacian of equivalent direction weighted Gor' kov. The left side 1310 is the x-y plane parallel to the transducer array. The right side 1320 is the x-z plane perpendicular to and cut through (slicing through) the flat transducer array.

Fig. 14 shows the laplacian of z-direction heavily weighted Gor' kov. The left side 1410 is the x-y plane parallel to the transducer array. The right side 1420 is the x-z plane that is perpendicular to and cuts through the planar transducer array.

Fig. 15 shows Gor' kov maximization under zero pressure conditions. The left side 1510 is an x-y plane parallel to the transducer array, depicting the spiral trap as shown by the non-uniform shading around the ring of uniform strength. The right side 1520 is the x-z plane that is perpendicular to and cuts through the underlying planar transducer array.

In fig. 16, the additional feature models 1610, 1620, 1630 correspond to trapping solutions that also solve the core problem. These are less effective because they do not maximize the objective function. This is because the value of the objective function (the modified Gor' kov potential) is small for these modes orthogonal to the dominant feature vector. The value of the implemented objective function is reflected in the feature value associated with each feature vector. The feature vectors associated with smaller feature values represent these alternative less efficient solutions.

3.Phased array device focusing with vector effects applied to semi-null haptics

I.Active sound intensity

The active acoustic intensity I can be defined using the pressure p and the particle velocity vector u of the acoustic mediuma. For acoustic waves, this can be described as:

Figure BDA0002163150260000181

where the active intensity is a purely real vector quantity. In the near field of the wave source, the reactive field (reactive field) instead dominates

Figure BDA0002163150260000182

However, in the far field, and therefore more than one or two wavelengths, Ir→ 0. For this reason, only I will be consideredaAnd assume that it is guaranteed to be real.

Consider a monochromatic plane wave source, where:

Figure BDA0002163150260000183

Figure BDA0002163150260000187

resulting in:

Figure BDA0002163150260000184

wherein

Figure BDA0002163150260000185

Is the pressure amplitude, and

Figure BDA0002163150260000186

is the particle velocity amplitude in each direction.

The acoustic intensity and hence energy flux of a wave propagating along a given axis is of interest, so the dot product used to measure the scalar contribution is:

Figure BDA0002163150260000191

due to the form of the acoustic wave equation:

Figure BDA0002163150260000192

Figure BDA0002163150260000193

the principle of superposition must be applied to both pressure and particle velocity vectors. Assuming n sound sources and approximating each sound source as a plane wave at the point of interest, the equation results:

for efficiency, it is desirable to solve for the appropriate variables using a linear system. This means that it cannot be solved for the sound intensity in a given direction, since it consists of two independent variables. Prior art solution to complex-valued scalar pressures

Figure BDA0002163150260000195

And then assume that the direction of each of the individual waves is insignificant or approximately parallel. Instead, the particle velocity vector is solved

Figure BDA0002163150260000196

Since this directly affects momentum (and kinetic energy), it is important for mid-air haptic feedback applications.

This results in extending the current solution mechanism to take advantage of complex-valued particle velocity vectors

Figure BDA0002163150260000197

Two potential techniques to create a solution that respects the directionality of the wave and its resulting energy transfer.

II.Solving method of complex value pressure

The solution method for complex-valued pressures involves the use of an acoustic model of the pressure. Given the transducer position and the control point position, the model can be queried and the resulting phasor at that point in the field determined. This is a complex number of pressures that,

Figure BDA0002163150260000198

wherein j corresponds to the index of a sample point in the sound field, anq corresponds to the index number of the actuated transducer. This enables the creation of an underdetermined linearity system:

Figure BDA0002163150260000201

where b is the phasor in air at position 1, …, j, …, m and x is the transducer activation coefficient from transducer 1, …, q, …, n, the initial phase and amplitude required to create the air behavior specified by b.

However, since b is complex and contains phase and amplitude specifications, the phase offset can be freely modified. This can be used to find the phase value that opens the largest possible amplitude subspace. This is achieved by finding the dominant eigenvector y of the system, Ly ═ λ y:

Figure BDA0002163150260000202

where a represents the desired amplitude at each point in the sound field and k represents the per-column normalization applied, which cancels out the dot product of 1 in each diagonal of the matrix for each control point in the sound field 1, … …, j, …, m. When multiplied by the initial amplitude a, the normalized component of the eigenvector y of the matrix will then correspond to the phasor vector b, which promotes enhancement between each pressure point in the field.

The minimum norm solution for the underdetermined problem can be solved using the following form:

xleastnorm=AH(AAH)-1b。

the key to this solution is the inverse of the matrix, C ═ AAHIt is the first step to solve the system via Cholesky decomposition and is written as:

Figure BDA0002163150260000203

meaning that it may also help to derive the eigensystem matrix L. However, it also works with the following form of system:

wfocus=C-1b,

resulting in a very small linear system and subsequent solution size that can be solved using any technique (QR decomposition, Cholesky factorization) and held in this minuscule state until transmitted. The C matrix and the wfocusThe vectors are not indexed in the transducer until the following is calculated:

xleastnorm=AHwfocus

usually worth passing through multiplication by a at the last possible momentHTo convert back to xleastnormThe transducers are activated, which multiplies the waves emitted from each transducer, since the transducers can be numerous.

III.Solving for one or more anisotropic control points

Creating control points in the sound field such that they each have a normal vector or direction specified with them is one way to incorporate the momentum requirement into the phased array solver. Giving each control point a normal vector in effect means solving for the x, y and z particle velocities and hence the components of the amount of fluctuation at the point of interest, rather than the pressure at each control point. To achieve this, we extend the linear system equation from the previous to include three components of each particle velocity (instead of the previously defined components of pressure). These three components (designated as alpha) are thenx,j0q、αy,j0qAnd alphaz,j0q) Is defined as:

Figure BDA0002163150260000212

Figure BDA0002163150260000213

in each of the cases, the number of the cases,

Figure BDA0002163150260000214

is a vector representing the direction of the wave propagating from transducer q at control point j, andis the normal vector direction of control point j, the momentum is measured across its position. The core linear system to be solved then becomes:

which specifies three components of the velocity of the particles of the acoustic medium at each control point. However, when the solution to the system uses a simplified representation approach, i.e. the specification of the C matrix to reduce the required workload, the number of elements (entries) of the matrix is increased by a factor of nine, making the problem significantly more computationally intensive. For each of the elements of the previous matrix C,

Figure BDA0002163150260000222

becoming a sub-matrix:

Figure BDA0002163150260000223

this still performs well when using the existing method in the feature system matrix L, although the user input amplitude becomes

Figure BDA0002163150260000224

And

Figure BDA0002163150260000225

or

Figure BDA0002163150260000226

But it can still be constructed using the same method. When considering the effect of adding waves on another axis on another control point, the cross terms of the eigensystem matrix can be considered to be performed on one axis at the position of one control pointAnd (6) measuring the rows. Albeit wfocusThe vector now contains three times more elements, but multiplied by AHThe system can be returned to the transducer activation coefficient required to reproduce the desired momentum vector in the air above the device.

This solution is particularly interesting for the application of semi-empty haptic feedback, as it opens up the possibility of creating waves approaching the hand from different angles. This may enable the acoustic waves to induce shear waves in the skin on the boundary, which brings a completely new set of sensations that can be created using the half-space haptic technique, at the cost of the extra required calculations.

IV.Solving for particle velocity at one or more control points in an acoustic medium

Another technique to achieve a conversion from a pressure value-targeted solver to a momentum-capable solver is to create a linear system solver that works entirely with acoustic particle velocity, which (similar to the original pressure method) is also complex-valued due to monochromatic harmonic times. Further, it is a scalar. This can be achieved by: calculating the direction of maximum particle velocity at the control point and thus the composite wave direction

Figure BDA0002163150260000231

And then multiplied to make each element a single complex-valued scalar. Since the dot product affects all sampling directions equally across x, y and z by the scaling factor, this can be applied as post-processing to the initial matrix construction. However, the direction of the composite wave must be known before the solver can function, and the direction of the composite wave cannot therefore be affected by the solver. When applied to the construction of a semi-empty haptic device, this approach would be most suitable for small devices where semi-empty but relatively close-range feedback is optimal, and is an extension of the current solution mechanism.

To achieve this, the acoustic model of each transducer is re-set to the problem of determining the velocity of particles in the acoustic medium. Therefore, a system that needs to solve for a given control point has

Figure BDA0002163150260000232

In the form of alpha vectors, to produce complex-valued particle velocities measured along the direction of maximum motion,

Figure BDA0002163150260000233

unfortunately, however, it is not known a priori

Figure BDA0002163150260000234

And is required to construct a linear system matrix:

Figure BDA0002163150260000235

thus, the acoustic model may be modified to yield the following values:

Figure BDA0002163150260000236

Figure BDA0002163150260000237

Figure BDA0002163150260000241

wherein

Figure BDA0002163150260000242

These values must then be used to construct the dot product of the C matrix:

Figure BDA0002163150260000243

there are two potential ways to achieve this. The required dot product is defined by:

Figure BDA0002163150260000244

this can be written as:

wherein

Figure BDA0002163150260000246

Is a dyadic product. This requires that nine separate summations be maintained for each matrix element (three diagonal) in order to subsequently calculate the normal vector. This is beneficial because it results in the initial matrix configuration of C having many similarities to the initial matrix configuration of the anisotropic control point method, and thus may allow for a common implementation infrastructure that supports both technologies. However, in many embodiments, an alternative to having a two-stage modeling system may instead be more desirable than maintaining the nine separate summing variables. In such a system, the first stage would run a simplified model that calculates only velocity amplitude and wave direction, yielding three velocity directions summed per control point

Figure BDA0002163150260000251

And

Figure BDA0002163150260000252

since this is separate from the phase calculation, this can be a simplified model. These sums each give a component of the final strongest amplitude and can therefore be used to calculate the direction in which the sums are made in the second stage.

To sum up:

Figure BDA0002163150260000253

wherein in a second phase, the model is capable of calculating:

Figure BDA0002163150260000254

inside the transducer model before starting the other stages of the process. The remaining phases proceed in a very similar manner to before, except that the direction of the velocity effect is different.

This approach allows the application of mid-air haptic feedback to more directly handle momentum transfer to the hand, rather than pressure. This leads to an improvement in the repeatability of the sensation and thus to an improvement in the measurement of the perception behavior, which is crucial for the extension of the half-space haptic as a commercial technology.

V.Drawings

Fig. 17 shows that the same wave with the same Sound Pressure Level (SPL) produces different effects at the surface depending on the angle of incidence. At the left side 1710, waves (left) traveling parallel to the surface (right) strike with some energy E, reflecting and exiting with slightly less energy due to momentum exchange. On the right side 1720, waves traveling perpendicular to the same surface do not reflect and exchange energy, and therefore exit with the same amount of energy. Therefore, it is impossible to impart a haptic effect.

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 dB SPL, dB SIL, and dB SVL created by the time-of-flight method, respectively, with a light gray of +10dB, a white of 0dB, and a dark gray of-10 dB. The bottom row 1830(dB SPL-dBSIL), (dB SPL-dB SVL), and (dB SIL-dB SVL), respectively, with the same color scheme, reveals the difference between each measurement method. Because the SVL field measurement models the velocity of the medium-based particles, the SVL field measurement simulation may allow the waves to cancel further due to additional vector component cancellation. Thus, the dB SVL is always less than the dB SIL, which is less than the dB SPL, but they are separated by similar distances in the log dB space due to the multiplication of their respective components with the effective geometric mean of the SIL measurements.

4.Eigensystem solution for quadratic problem in semi-air haptic system

I.Rayleigh-Ritz theorem for quadratic optimization

Suppose that the matrix M is a square-and-half positive definite Hermitian (Hermitian) d × d matrix. Since the matrix is a square matrix and a hermitian matrix, each eigenvalue must be a real number. Since it is semi-positive, these eigenvalues must be positive, e.g., xHMx≥0。

The feature values may then be organized as:

λmin≤λ1≤λ2≤…≤λd≤λmax

lambda of the matrixminAnd λmaxIt can then be shown as a solution to the optimization problem:

Figure BDA0002163150260000261

and is

Figure BDA0002163150260000262

This is the Rayleigh-Ritz theorem. It follows that the x-vector from each of these can be found as the corresponding feature vector for the found feature values.

II.Quadratic optimization of phased array systems

When considering a phased array system, x is defined as followsHIt is useful that:

Figure BDA0002163150260000271

where n is the number of monochromatic sound sources so that the vector solution can be used to drive the array of source transducers. Matrix M is then defined as an n N Hermite matrix such that when x is a unity norm (| x21) feature vector:

Figure BDA0002163150260000272

this then represents the optimized form. What is left to do is to fill the matrix M with a particular pure quadratic (since the matrix M needs to be a symmetric matrix or a hermitian matrix) objective function.

There are five potentially useful pure quadratic forms for the objective function that produces the hermitian matrix M, and thus it can be maximized using existing feature vector solvers to obtain a solution field that can be used to drive the half-space ultrasonic haptics:

1. maximum pressure optimization (total SPL): a. maximization

Figure BDA0002163150260000273

It is proportional to the square of the pressure; b. this is an objective function commonly used to drive phased array systems.

2. Maximum particle velocity optimization (total SVL): a. maximization

Figure BDA0002163150260000274

Proportional to the square of the particle velocity; b. this takes into account the directionality of the wave, as this is necessary to obtain a consistent tactile perception.

3. Maximum particle velocity is optimized, but in a given direction (SVL available): a. maximization

Figure BDA0002163150260000275

Which is proportional to the square of the particle velocity along a given direction; b. this takes into account the directionality of the wave and the (skin) surface

Figure BDA0002163150260000276

The impact is orthogonal.

4. Maximum sound intensity (energy flux density) optimization (total SIL): a. maximizationIt is proportional to the sound intensity; b. this takes into account the directionality of the wave as it is necessary to obtain a consistent tactile perception.

5. Maximum sound intensity (energy flux density) optimization along a given direction (usable SIL): a. maximization

Figure BDA0002163150260000282

Proportional to the sound intensity along a given direction; b. this takes into account the directionality of the wave and the (skin) surface

Figure BDA0002163150260000283

The impact is orthogonal.

III.Optimized to maximize pressure

A quadratic form optimizer that can maximize the pressure at each control point will model the following quantities:

Figure BDA0002163150260000284

at each control point j, a summation of the objective functions is generated

Figure BDA0002163150260000285

It makes the principal eigenvector of each element of the n × n matrix solve for:

Figure BDA0002163150260000286

wherein

Figure BDA0002163150260000287

Is a positive and real amplitude of the wave pressure, and

Figure BDA0002163150260000288

is the complex phase offset on the wave traveling from the transducer q at the control point j of the initial transmission. Weight wjIs used for the quadratic form

Figure BDA0002163150260000289

The terms are re-weighted to provide relative amplitude control between the control points. Global amplitude scaling is provided by a unit norm condition with respect to the input power, xHx is 1. This (in combination with the quadratic reweighting factor and the control loop) can be used to accurately process the individual pressure amplitudes.

IV.Optimized to maximize particle velocity

A quadratic form optimizer that can maximize particle velocity will model the following quantities:

Figure BDA0002163150260000291

for each control point j. This will result in a summation of the objective function:

Figure BDA0002163150260000292

it makes the principal eigenvector of each element of the n × n matrix solve for:

Figure BDA0002163150260000293

wherein u isxjqIs the x-direction component of the induced particle velocity amplitude (which is real and positive), and

Figure BDA0002163150260000294

is the complex phase offset on the wave traveling from the transducer q at the control point j of the initial transmission. Weight wjIs used for the quadratic form

Figure BDA0002163150260000295

The terms are re-weighted to provide relative particle velocity magnitude control between control points. Global amplitude scaling is provided by a unit norm condition with respect to the input power, xHx is 1. This (in combination with the quadratic re-weighting factor and the control loop) can be used to accurately process the individual velocity amplitudes.

V.Optimized to maximize particle velocity in a given set of directions

A quadratic form optimizer that can maximize particle velocity in a given set of directions will model the following quantities:

Figure BDA0002163150260000301

for each control point j. This will result in a summation of the objective function:

Figure BDA0002163150260000302

it makes the principal eigenvector of each element of the n × n matrix solve for:

Figure BDA0002163150260000303

wherein

Figure BDA0002163150260000304

Is the x-direction component of the induced particle velocity amplitude (which is real and positive), and

Figure BDA0002163150260000305

is the complex phase offset on the wave traveling from the transducer q at the control point j of the initial transmission. Weight wjIs used for the quadratic form

Figure BDA0002163150260000306

The terms are re-weighted to provide relative amplitude control between the control points. Global amplitude scaling is provided by a unit norm condition with respect to the input power, xHx is 1. This (in combination with the quadratic re-weighting factor and the control loop) can be used to accurately process the respective amplitudes.

VI.Active sound intensity

The active acoustic intensity I can be defined using the pressure p and the particle velocity vector u of the acoustic mediuma. For acoustic waves, this can be described as:

where the active intensity is a purely real vector quantity. In the near field of the wave source, the reaction field is instead dominant

Figure BDA0002163150260000312

However, in the far field, and therefore more than one or two wavelengths, Ir→ 0. NeedleFor this reason, only I is consideredaAnd assume that it is guaranteed to be real. Consider a monochromatic plane wave source, where:

Figure BDA0002163150260000313

we have:

wherein

Figure BDA0002163150260000316

Is the pressure amplitude, and

Figure BDA0002163150260000317

is the particle velocity amplitude in each direction.

The acoustic intensity and hence energy flux of a wave propagating along a given axis is also of interest, so the dot product for measuring scalar contributions is formed:

Figure BDA0002163150260000318

due to the form of the acoustic wave equation:

Figure BDA00021631502600003110

the principle of superposition must be applied to both pressure and particle velocity vectors. Assuming n sound sources and approximating each sound source as a plane wave at the point of interest, the equation results:

Figure BDA00021631502600003112

it is known that, in the modality in question,

Figure BDA0002163150260000321

and p and

Figure BDA0002163150260000322

in phase, it is known that the placement of the conjugate does not matter, and the order of the conjugate can be reversed to conjugate the pressure rather than the particle velocity (this can be verified by extending the above equation and using commutativity to reverse the order of the complex exponentials in each term). According to a corollary, this can be written in such a way as to generate a hermitian matrix M:

Figure BDA0002163150260000323

Figure BDA0002163150260000324

VII.optimized to maximize sound intensity

The first quantity to generate a scalar value is:

Figure BDA0002163150260000325

but assuming it is in the far field, this can also be written as:

this is the fourth power of the pressure sum, but due to the structure of the complex exponentials, two of these are identical and can therefore be removed by the square root.

Thus, a quadratic form optimizer that can maximize the sound intensity will model the following quantities:

Figure BDA0002163150260000331

for each control point j. This will result in a summation of the objective function:

Figure BDA0002163150260000332

it makes the principal eigenvector of each element of the n × n matrix solve for:

Figure BDA0002163150260000333

wherein

Figure BDA0002163150260000334

Is the x-direction component of the normal vector of the wave front,

Figure BDA0002163150260000335

is a positive and real amplitude of the wave pressure, and

Figure BDA0002163150260000336

is the complex phase offset on the wave traveling from the transducer q at the control point j of the initial transmission. Weight wjIs used to re-weight the quadratic term, which is still effectively the square root of a function involving the fourth power of pressure,

Figure BDA0002163150260000337

to provide relative particle velocity amplitude control between control points. Global amplitude scaling is provided by a unit norm condition with respect to the input power, xHx is 1. This (in combination with the quadratic re-weighting factor and the control loop) can be used to accurately process the individual velocity amplitudes.

VIII.Optimized to maximize sound intensity in a given set of directions

A quadratic form optimizer that can maximize the sound intensity along a given direction will model the following quantities:

Figure BDA0002163150260000341

but then assuming it is in the far field, this can also be written as:

Figure BDA0002163150260000342

for each control point j. This will result in a summation of the objective function:

Figure BDA0002163150260000343

it makes the principal eigenvector of each element of the n × n matrix solve for:

Figure BDA0002163150260000344

wherein

Figure BDA0002163150260000345

Is the x-direction component of the induced particle velocity amplitude (which is real and positive),is a positive and real amplitude of the wave pressure, and

Figure BDA0002163150260000347

is the complex phase offset on the wave traveling from the transducer q at the control point j of the initial transmission. Weight wjIs used for the quadratic form

Figure BDA0002163150260000348

The terms are re-weighted to provide relative particle velocity magnitude control between control points. Global amplitude scaling is provided by a unit norm condition with respect to the input power, xHx is 1. This (in combination with the quadratic re-weighting factor and the control loop) can be used to accurately process the individual velocity amplitudes.

IX.Finding weighting factors

Rayleigh-Ritz theorem shows that any of the above matrices can be maximized by finding the dominant eigenvector and allowing the power to increase up to the physical limits of the transducer. This can be achieved by power iteration. By constructing the matrix M while individually defining each of the sums of the control points, respectively, the evaluation of the objective function can be performed as an additional side effect when multiplying by the vector x. This can be used to find the weight wjSuch that they allow maximization to find the magnitude of the solution that fits the system when scaled.

As described in the power iteration:

Figure BDA0002163150260000351

wherein:

Figure BDA0002163150260000352

we can use power iteration to determine the objective function at each control point.

Defining:

Figure BDA0002163150260000353

and:

Figure BDA0002163150260000354

at each iteration, it may instead calculate:

xa,j=Mjxa-1

thus, the objective function at each control point can be evaluated as:

Figure BDA0002163150260000355

finally, x can be calculated as beforea

Figure BDA0002163150260000361

With the final feature vector solution as above. However, calculating it in this way allows oj(xa-1) An intermediate result for control point j is generated. To determine the level at which this represents the desired control point output, the average objective function value is taken to allow the set of values oj(xa-1) A direct comparison is made with the target value given by the user. In this case, any method for optimizing the derivative-free of the loss function (such as, for example, the Nelder-Mad algorithm) may be used to update the weighting factors for the next iteration.

Intermediate results x after several power iteration methods and coupled derivative-free optimizationaShould be used to drive the transducer. This will therefore allow the system to track the expected value of the objective function over time, while allowing some error. The advantage of this system is that it is a very simple way of driving the transducer array. The biggest disadvantage is that many calculations need to be performed, but the calculation is actually only a complex multiplication, which is implemented directly with enough hardware.

X.Drawings

Shown in fig. 19 are: top left 1910: the control points created by solving the quadratic particle velocity problem with non-perpendicular vectors are shown as intersections through the x-y plane. Right top portion 1920: shown as the same control point passing through the intersection of the x-z planes. Left bottom 1930: for comparison, control points were created by solving a quadratic pressure problem, shown as intersection points through the x-y plane. Right bottom 1940: by solving the intersection of the pressure problem on the x-z plane.

5.Control region of phased array specified by quadratic problem

The present disclosure describes a method for automatically generating basis functions using user-provided design parameters. These parameters are input in an online fashion to an efficient optimization of a quadratic function in space through a given transducer basis set. This gives the user an accessible design space through online optimization, which can generate a variety of different effects. Using a linear operator, called the 'quadratic control region', since it is the control region of the quadratic problem, a matrix can be constructed by evaluating the sum-product of the transducer effects. The principal eigenvector of the matrix may be selected as the best solution in a set of orthonormal solutions that are equivalent to the eigenvector of the matrix. This optimization yields as output the transducer activation coefficients — the original complex values needed to drive the individual transducer elements to recreate an acoustic field or phenomenon with the desired properties.

These original complex values can then be used to further include the basis functions of later linear systems, the solution of which enables the user to produce all desired effects simultaneously, while re-weighting the output power available for each.

I.Linear control region

A linear control region defined as an operator that integrates a linear acoustic quantity, such as pressure, with a complex-valued weighting function w (x, y, z):

Figure BDA0002163150260000371

may be used to evaluate the sound field. It should be noted that when w (x, y, z) is paired with a similar phase distribution in the acoustic quantities, the complex conjugate is used for the maximization function evaluation. Any approximation of this integral (such as sampling) may be used, or indeed this may be used to approximate any linear operator on the acoustic pressure field. Additional integrals can be constructed to optimize the control region using vector particle velocity, as long as it is also linear, where:

Figure BDA0002163150260000372

where u is the particle velocity of the medium and w is a three-dimensional vector field, each dimension being a complex value, which is used as a weighting function with the same purpose as above. In this way, different acoustic magnitude patterns can be evaluated in a linear fashion, maximizing the particular phase of the result when the phase is aligned with a 'template' created by a scalar or vector function in w.

It should also be noted that w (x, y, z) may be mostly zero (or consist of a zero vector) to create a single evaluation of multiple seemingly disjoint regions. Any integral approximation (such as sampling) can be used, as the formula can be extended to any linear operator on the field of linear acoustic quantities.

A phased array can be described by controlling a series of complex valued coefficients for each transducer. Driving an acoustic transducer array at a fixed frequency with such an input x-vector to maximize the acoustic quantity encoded in the matrix M can be represented as a linear algebraic problem, so the intention is:

maximizing xHMx

Is subject to xHx=1,

And is

Figure BDA0002163150260000381

This can be solved by taking the principal eigenvector of the matrix M, since the statement of the problem is also the definition of the eigenvector (here x). The construction of the matrix M can be achieved by first considering the linear acoustic quantities emitted by each transducer as basis functions in the global space of the acoustic quantities in question. This is the way in which maximization of pressure and other acoustic quantities is achieved at a single point described by the absolute position in x, y and z. However, using the described control region integration, q (x, y, z) which controls the amount of acoustic radiation in the field transmitted to each individual transducer element can be replaced by f (x, y, z), resulting in a column vector:

Figure BDA0002163150260000382

by invoking linearity, the complex activation values that drive each transducer can be multiplied, resulting in a vector:

Figure BDA0002163150260000383

where n is the number of sensors. Given that m is the result of a linear operator operating on a linear field, it can be hermitian transposed and multiplied by itself to generate a quadratic form. Finally, the process is carried out in a batch,

mHm=xHMx=xHffHx。

by modifying the weighting functions w (x, y, z) and solving the eigenvectors of the matrix M, various beams, surfaces and other acoustic phenomena can be created and maximized, resulting in regions of optimal actuation in the acoustic field.

The previously discussed null-space generation technique ("eigensystem solution to the quadratic problem in the semihollow haptic system") can also be used to extend the method to describe acoustic levitation and allow null steering (null steering) and formation of the semihollow texture by introducing a null into the field at the location of the described structure. It should be noted that this procedure increases the rank of the matrix, so simplification is not possible if this is used.

Likewise, for many control regions, multiple optimizations may be performed by adding the M matrices together. This may also include a matrix representing the control points. This generates a matrix with an increased rank, so some simplified methods are not possible using this.

By taking the result into account and re-weighting each point in an iterative manner, the system can be modified so as to converge to the desired geometrical configuration of the acoustic quantity. This may typically be required when the amplitudes are not uniform in the output.

This can be achieved by summing the generated acoustic quantities and re-weighting the system closer to the result in the next iteration. For example, the iterative re-weighting of any sub-divided volume or single pressure point for integration becomes:

Figure BDA0002163150260000391

where t +1 refers to the next coarse iteration of the technique, t is the current iteration, prIs the desired pressure (which may be complex) and the denominator is the total pressure obtained by the current iteration. The fraction can be presented only as a real magnitude, so as to affect only the amplitude of the resulting region and not the phase.

Bessel beams may be created by: the weighting function w (x, y, z) is expressed in the beam traveling direction, and the unit pressure is created while modifying the phase angle forward as the beam moves forward. This can be sampled with a sample book as long as the angle traversed by each subsequent spot is small enough to sample the phase along the beam axis.

The sampling of the phase must have significantly more than two phase samples for every 2 pi angle of phase shift. This is because if the increment along the path is pi, the direction of wave travel is indeterminate. This also means that in most cases the phase must be incremented forward so that it completes a 2 pi radian rotation for each wavelength of the monochromatic frequency distance it covers. This is an approximation to the path integral of the form described in the control region definition.

As a result of maximizing the integration written in this manner, the transducer will generate an optimized Bessel-like beam in the case of a sample beam progressively defined in a straight line moving away from the array, due to its theoretical ability to generate optimized pressure in this configuration. The intersection of the passing fields at the edge of the near field region parallel to the planar linear transducer array is shown in fig. 20, revealing the structure of the complex valued transducer activation coefficients. A simple progressive phase beam produces a Bessel function approximation (left 2000) and an undersampled helical beam produces a Bessel-like function with dipole behavior (right 2010), where the dipoles create two structures that are activated in opposite phase to each other. Due to limited transducer resolution, limited array size, and the effects of other physical considerations that optimization will need to face, the beam may not be an accurate Bessel beam, but may work more efficiently than a simple evaluation of the Bessel beam function due to optimization that takes into account the physical limitations of the real transducer array. On large arrays, which closely follow the Bessel function, fig. 21 shows a plot 2100 of the Bessel function plotted along pressure when considered radially from the center of the array in fig. 20.

An additional property of the system is that transducers of any position and orientation with any generated monochromatic field function can be used. This is a significant advantage of generating Bessel beams by the sampling method, which by definition will only work in the case of planar arrays.

By creating multiple out-of-phase paths, phase shifts and phase singularities can be implicitly created by wrapping a helicon wave around the paths. The path may also be curved, resulting in a large degree of configurability, while using simple quadratic optimization embedded in the eigensystem approach. This is achieved by defining separate phase windings in many parallel beam lines that are closer to each other than the beam width but are out of phase with each other to cause rotation in the beam. Then, by setting the nearby phases to cancel each other in the center line of the beam, a high-order Bessel beam can be formed, as shown in fig. 22. In fig. 22, more than two samples of progressive phase around a circle (hence, the phase in a spiral arrangement) results in the formation of a Bessel function with angular momentum, as seen at the left side 2210. By using more than 2 pi radians to specify the phase around the helix, higher order Bessel beams can be generated, such as the 4 pi beam at the right side 2220.

These phase lines may be sampled using points or any other integral approximation as before, as these are examples of how disjoint regions may be defined as a single region optimization, where each may be added to the optimized region by setting the weighting function to non-zero only in the volume to be optimized. Higher order Mathieu-type elliptical spiral beams may be created by incrementally phased beams that are incrementally phased and spaced around the ellipse. As shown in fig. 23, an elliptical helical beam can be created by parameterizing the circumference of the ellipse in phase, which has similar properties to higher order Bessel beams (such as Mathieu-type beams). Determining the ratio of 2n pi radians and assigning a phase angle proportional to the elliptical chordal circumference generates an appropriate configuration for optimization of the beam. The left image 2310 shows a higher order Bessel that is reshaped into an elliptical, near Matheiu beam. Right image 2320 shows an approximate Matheiu beam with a cross-linked helical structure that begins to form in the central region. Both based on a phase circle of 12 pi radians, creating a helical wavefront.

The family of Bessel beams created in this way can be modified to produce more exotic beams that exist in the literature in a more constrained form. By creating an arbitrary path in the field traversed by a physically feasible phase function, which moves forward in a manner suitable for the traversal of the wave, various results can be reproduced, but with more general and lower formation on initial bandwidth cost than before due to the optimizations implied by the eigensystem, since the beam only needs to be represented by an integral.

Taking a curved path and moving the phase forward can produce a curved beam, typically used to steer the beam around obstacles created by a 45 ° arc followed by a path that pushes the formed beam out to infinity. Shown in fig. 24 is an image 2400 of a curved beam formed by a phase function following a 45 degree segment connected to a circle continuing along a line, using a coarse iterative approach to equalize the acoustic pressure along the beam. The transducer array is disposed along a lower edge of the image.

Additional methods involving multiple phase paths may create "local hollow beams," beams with a section cut from the center, which may facilitate levitation of large objects. This has been generated in fig. 25, which shows an image 2500 that pairs two gaussian curves along which the phase has been made to follow. The local hollow beams generated by creating the phase function follow a gaussian curve using a coarse iterative approach to equalize the sound pressure along the beam. The transducer array is disposed along a lower edge of the image. The dashed nature of the lower region indicates that the phase speed input to the optimization is too slow and that the beam can be improved by using a faster phase speed for the initial portion. Portions of curves may also be used, such as in fig. 26, because if the object is sufficiently rigid, only some portions of the space need be made to be high acoustic pressure regions. Fig. 26 shows an image 2600 of an incomplete local hollow beam formed by portions of a circle defined with progressive phase, which can be used to levitate even larger objects by optimizing the high acoustic pressure region in a non-uniform manner. If it is assumed that the object or substance to be suspended has sufficient cohesion without flowing, extruding or leaking through the low acoustic pressure region, it is not necessary to evenly distribute the force exerted on the large object.

A Bessel-like beam can be created by moving the path forward faster than the phase velocity. In this way, the beam will be optimized at a faster phase speed than usual, creating a beam that appears to have a longer wavelength and that travels faster than the speed of sound due to the phase speed, as shown in fig. 27. Shown in fig. 27 is an image 2710 of a wavefront that travels faster than the speed of sound, which can be created by implementing a phase path that is angularly shifted in a manner that appears paradoxically generates longer wavelengths with the same frequency. Implicit optimization infers that the main beam is composed of cross beams, resulting in a specified increase in phase velocity. This can be seen from the right image 2720, where the wavelength is steered approximately 40% longer at the same frequency for the initial portion of the beam. A serpentine beam can be constructed by sampling the progressive phase along a sinusoid. Shown in fig. 28 is an example 2810 where a short portion of a sinusoid is defined, followed by 2820 of a longer two-dimensional serpentine beam constructed from progressive phases of pressure over a phased array following the sinusoid. These show that extending a serpentine sinusoidal beam by simply adding an additional portion of progressive phase causes the beam to be naturally extended as long as needed (and allowed by the physical arrangement of the transducers). By creating a serpentine beam in all three dimensions to give a helical serpentine beam, a beam can be created that can be used in a manner similar to an archimedean spiral system to potentially apply forces to solid objects in the field that are much larger than a wavelength, in a direction perpendicular to or opposite to the emission of the wave. This gives a different mechanism whereby the beam can be used to levitate an object. FIG. 29 shows a left image 2910 detailing a portion in the x-z plane and a right image 2920 showing a portion of a beam in the x-y plane. Using this technique, which follows a path in a circle wound in the x-direction and the y-direction (but which is also progressive in the z-direction), a true spiral formation can be created. This can be used to levitate and apply forces to objects much larger than a wavelength in a manner similar to a tow beam.

If the eigensystem matrix is defined by a single control region, the matrix has rank 1 and is formed by the sum-vector product of the complex-valued vector m and itself, defining the matrix as fHf. In this particular case, the single feature vector of the system is known and is f, meaning that f ═ x. By evaluating the objective function in this case, the objective function itself specifies the coefficients of the desired complex-valued transducer activation. This is defined by simply evaluating the control region function of each transducer, which becomes the activation coefficient.

As before, due to the quadratic form, evaluation of validity by re-weighting each iteration may also be performed. This is done by evaluating the pressure at each point or the pressure integral over each sub-area and dividing by the required pressure, thereby generating a new function value to update wt,i(x, y, z), the weight of sub-region i. The single feature vector, also f, can then be reevaluated, generating an equivalent method from repeated iterations on the previous feature system, but at a greatly reduced computational cost.

To achieve this with the particle velocity of the medium, we can use a test of the particle velocity, which is equivalent to changing direction until a maximum is reached:

Figure BDA0002163150260000431

although prior art techniques may also be used to monitor the pressure to create the mixing process. This hybrid technique involves specifying a weighting function and updating the pressure as before, but with a valid control region:

Figure BDA0002163150260000432

a scalar but variable and complex-valued w (x, y, z) is generated, which is paired with a constant weighted normal vector v and a particle velocity u. The method can be optimized to find a satisfactory pressure without changing the direction of the waves sought. This is also shown in fig. 28, in addition to the specific corrections described in the next section.

II.Wavefront filtering of directional control points and regions

An additional technique necessary to create a directional output is to wavefront the transducer output

Figure BDA0002163150260000433

Transducer data is discarded and set to zero as opposed to the normal vector or particle velocity vector u, thus satisfying:

Figure BDA0002163150260000434

for points in the area under consideration.

Fig. 30 shows a left image 3010 of two velocity vector samples further defined from the transducer array and a right image 3020 of two velocity vector samples defined proximate the array with the transducers having their wave fronts removed in a direction opposite the defined control direction. The velocity weighting functions specified in the two disjoint regions are used to form a single region in each case, wherein this region is optimized for the direction vector representing the particle velocity of the medium at two separate points.

Note that when trying to find the maximum amount of velocity in the x-direction, it does not necessarily mean that the resulting particle velocity direction will be entirely in x. While fig. 30 shows the results of filtering out conflicting transducers that provide a wavefront opposite to the desired direction, fig. 31 shows a left image 3110 and a right image 3120 that are the results of not filtering the transducer wavefront direction using the wavefront standard. In the case of fig. 31, the conflicting transducers are free to interfere, with pathology (pathology) evident in the form of standing wave structures on the horizontal. This is because the wavefront direction conveys the missing directional information in the particle velocity, which does not yield information about which direction the wave originates from, thus producing undesirable interference.

This is important not only for the simple case of control areas but also for normal control points. The case of a control point is actually a simplified control region where it is approximated with a single point sample of particle velocity without a weighting factor. In either case, not filtering out 'useless' transducers produces results similar to the effect of the 'standing wave' phenomenon. This is because the particle velocity vector with complex components cannot determine the difference between waves traveling in different directions: in this model, a wave traveling from left to right with a phase 0 is used for the same point sample mathematically as a wave traveling from right to left with a phase pi. Thus, if this is not handled correctly on the basis of generating control points or control regions, the optimization will attempt to use both, generating a solution with standing wave components that do not make meaningful contributions, creating additional null, field-free regions, which are not needed. This is an important point, since control points input as basis functions into linear optimization, eigensystems or even simple formulas without the filtering step will exhibit these artifacts.

The techniques described may be integrated into a device to create an acoustic phased array system that is capable of producing arbitrarily configured regions and beams by iterating through a single feature system as defined herein or in the simple case of creating a single point or region. In the case of a single point or region, optimization can be 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. The iteration of the global method and the repeated application of the transit time may be changed and modified in a time-dependent manner to produce various haptic effects. This can be applied to the desired field and derived transducer coefficients as they are functionally related.

III.The latent claims

1. A method on a monochromatic sound field:

-evaluating the distribution of the acoustic pressure as a single value by integrating the product of the acoustic pressure field function and the weighting function, wherein any of said functions can be complex-valued and derived from the physical transducer properties.

-applying zero or more zeros to be added to the null space of the matrix.

Using this as the objective function to be maximized in the optimization, wherein the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer.

-summing one or more such matrices to generate a system of features, wherein one or more feature vectors produce a specified field.

2. A method on a monochromatic sound field:

-evaluating the distribution of acoustic particle velocities as a single vector value by integrating the dot product of a real three-dimensional weighted vector field and a particle velocity field with a weighting function, any function of which can have a complex-valued component and is derived from the physical transducer properties.

-applying zero or more zeros to be added to the null space of the matrix.

Using this as the objective function to be maximized in the optimization, wherein the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer.

-summing one or more such matrices to generate a system of features, wherein one or more feature vectors produce a specified field.

3. A method on a monochromatic sound field:

-evaluating the distribution of acoustic particle velocities as a single vector value by integrating the dot product of a real three-dimensional weighted vector field and a particle velocity field with a weighting function, any function of which can have a complex-valued component and is derived from the physical transducer properties.

-setting to zero the contribution of a transducer whose wave front at that point has a direction opposite to the normal established at each evaluated subregion of the integration.

-applying zero or more zeros to be added to the null space of the matrix.

Using this as the objective function to be maximized in the optimization, wherein the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer.

-summing one or more such matrices to generate a system of features, wherein one or more feature vectors produce a specified field.

4. A method on a monochromatic sound field:

-evaluating the distribution of the acoustic pressure as a single value by integrating the product of the acoustic pressure field function and the weighting function, wherein any of said functions can be complex-valued and derived from the physical transducer properties.

-applying zero or more zeros to add to the null space of the matrix.

Using this as the objective function to be maximized in the optimization, wherein the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer.

-using the contribution of each transducer to describe the optimal phase and amplitude of the transducer input signal.

5. A method on a monochromatic sound field:

-evaluating the distribution of acoustic particle velocities as a single vector value by integrating the dot product of a real three-dimensional weighted vector field and a particle velocity field with a weighting function, any function of which can have a complex-valued component and is derived from the physical transducer properties.

-applying zero or more zeros to be added to the null space of the matrix.

Using this as the objective function to be maximized in the optimization, wherein the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer.

-using the contribution of each transducer to describe the optimal phase and amplitude of the transducer input signal.

6. A method on a monochromatic sound field:

-evaluating the distribution of acoustic particle velocities as a single vector value by integrating the dot product of a real three-dimensional weighted vector field and a particle velocity field with a weighting function, any function of which can have a complex-valued component and is derived from the physical transducer properties.

-setting to zero the contribution of a transducer whose wave front at that point has a direction opposite to the normal established at each evaluated subregion of the integration.

-applying zero or more zeros to be added to the null space of the matrix.

Using this as the objective function to be maximized in the optimization, wherein the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer.

-using the contribution of each transducer to describe the optimal phase and amplitude of the transducer input signal.

7. A method as claimed in claims 1 to 6, wherein a further device implements the method and uses the result to drive an acoustic transducer.

8. The method of claims 1 to 6, wherein a further device implements the method and uses the result as part of a linear basis set for driving an acoustic transducer.

9. The method of claims 1 to 6, wherein an acoustic beam is generated.

10. The method of claims 1 to 6, wherein volume is generated.

11. The method of claims 1 to 6, wherein a sound surface is generated.

12. A method on an acoustic field:

-using a linear combination of fields derived from the physical transducer properties, evaluating the acoustic particle velocity as a single value at a point.

-setting the contribution of a transducer whose wave front has a direction at said point opposite to the normal at said control point to zero.

Using this as a basis function in an optimization problem for finding the phase and amplitude of the transducer input signal.

13. A method on an acoustic field:

-using a linear combination of fields derived from the physical transducer properties, evaluating the acoustic particle velocity as a single value at a point.

-setting the contribution of a transducer whose wave front has a direction at said point opposite to the normal at said control point to zero.

Using this as the objective function to be maximized in the optimization, where the optimization result takes the form of the eigenvector of the matrix formed by the cross-multiplication terms of the objective function of each transducer, which then finds the phase and amplitude of the transducer input signal.

14. A method on an acoustic field:

-using a linear combination of fields derived from the physical transducer properties, evaluating the acoustic particle velocity as a single value at a point.

-setting the contribution of a transducer whose wave front has a direction at said point opposite to the normal at said control point to zero.

-using the contribution of each transducer to describe the optimal phase and amplitude of the transducer input signal.

15. The method of claims 1-14, wherein the transducer produces a haptic effect when the resulting structure is applied as an input signal and modulated.

16. The method of claims 1-14, wherein the transducer produces a haptic effect when the resulting structure is applied as an input signal and modified by time.

6.Evaluation of transcendental functions, multiplications, and divisions of complex values in hardware

By aiming at

Figure BDA0002163150260000481

Rather than in

Figure BDA0002163150260000482

To convert between, the BKM algorithm can be redesigned to yield a greatly simplified technique, which is more useful in physical hardware implementations.

I.Introduction to

There are many techniques for the computation of trigonometric and other special functions in hardware. One such method is the BKM method, which replaces the complex logarithm with a complex exponential and vice versa, generating these special functions as side effects of moving from one representation to another. Of particular interest are those cases that deal with complex valued functions, as such functions are found to be frequently encountered where necessary.

An efficient method of moving between logarithmic and exponential forms enables multiplication and division to be performed as simple additions and subtractions in logarithmic space. This can also be exploited to easily generate an exponential-related transcendental function due to the nature of the process. Since this can be generated from a single algorithm, this results in reduced implementation size (since only one element needs to be implemented) and reduced development time (since only one method needs to be implemented and optimized).

The existing algorithm BKM implements this method, but includes three disadvantageous requirements that are almost or completely eliminated in the proposed method:

among other things, hardware division has difficulty implementing functions to pre-and post-processing steps with respect to data input and output.

Irrational values and factors create rounding problems for efficient implementation.

Additional operations are required to integrate with the floating point type.

The next section will describe how the proposed method is designed without these drawbacks.

II.Motivation and comparison range reduction

The process is implemented by EulerFormula as center, eCos θ + i sin θ, but wherein erIs another coefficient because the formula has a complex generalized form that yields er+iθ=er(cos θ + i sin θ). The complex function embodied by the BKM algorithm is exactly such that the process of implementation is

Figure BDA0002163150260000491

The algorithm heavily uses range reduction techniques to fit the complex value r + i θ into a range where convergence can be guaranteed. This is not a function with clean implementation, and to some extent is achieved by pre-or post-processing with irrational scaling factorsIt may be easier to evaluate using the proposed method. In many cases, this irrational scaling factor can be eliminated to provide a workflow whose inputs and outputs are intended to be multiples of physically or computationally meaningful rational numbers.

When dealing with complex exponentials, the exponent may be split into real and imaginary parts, which use the same input data but produce different results. This is similar in some respects to a cartesian to polar transformation.

This contains all the modulus information, similar to the r component of the cartesian to polar transform, as long as the real part of the logarithm of the number to be taken is considered. This is conventionally the real part of the irrational constant e to the power. To feed it into the algorithm, this first has to be range reduced, meaning we have to find the remainder, divide the modulus by the irrational value ln2, so that it fits inside the convergence domain. This involves hardware division (the 1994 paper recommends SRT division techniques) and therefore adds a very expensive technique in terms of the logic required to implement the division with the algorithm.

The premise of the proposed algorithm is that a power of 2 should be used to make the range reduction straightforward. With this motivation, the real part should be interpreted as a given power of 2, which is the real part and therefore the base 2 exponent. Then, an integer of a real logarithmic input may refer to a bit shift on the output, and only the fractional part of the relative binary logarithmic input needs to fit into the convergence domain. The integer part of the real component of the complex logarithm is then obtained by rounding to the nearest value and exponentiating, respectively, in this manner.

As long as the imaginary part of the logarithm of the number to be indexed is considered, this contains all argument or angle information, similar to the theta component of a cartesian-polar transformation. This is traditionally the power e, but e is not a problem here. The problem is pi because the imaginary logarithmic part can be interpreted as an angle, but to constrain the convergence domain we have to find the remainder again from the integer fraction of 2 pi so that rotational symmetry can be exploited. This in turn involves more division processes, which are expensive in terms of logical size.

If not a logarithm of the base e, a base 2 may be used. This is extended by 2Cos (θ ln 2) + i sin (θ ln 2). This is again a more undesirable irrational value, however, if we use the base e exclusivelyThen we can use θ to represent the number of rotations. Unfortunately, using this base results in a non-converging approach. The best compromise between technique and convergence range yields a base

Figure BDA0002163150260000501

This means that each integer value represents a quadrant. Every fourth multiple is a full rotation. The number of rotations can then be extracted using a bitwise operation. The number of quadrants may be applied to the output as a switch of real and imaginary numbers and signs. Again, this results in a technique in which the integer part of the imaginary component of the complex logarithm is taken by rounding to the nearest value and processing separately in this manner.

This allows the fractional part of the logarithm to still be indexed in the range of R [ -0.5, +0.5) + i [ -0.5, +0.5), which fits within the convergence range of the main algorithm, as shown in the right panel 3220 of fig. 32. In fig. 32, the exponential error plotted for the input of the log-to-exponential process in the region [ -2.0, +2.0), where the right side is positive real numbers and the lower part of the image is positive imaginary numbers. With 50 decimal places and a similar number of iterations, the black area indicates that convergence to greater than 2 has occurred-45The value of accuracy of (c). The left panel 3210 shows a convergence domain, where gray is mapped to greater than 2-45The error of (2). The right panel 3220 shows the final desired convergence region R [ -0.5, +0.5) + i [ -0.5, +0.5) in white.

A discussion of the main algorithm for the fractional part of the logarithmic input is described later.

The inverse process is also implemented, taking the values with their real and imaginary parts indexed and reducing them to their original logarithms. First, we find the quadrant to which the complex value belongs, offset by 45 °. This is achieved by determining which component has the largest absolute value; whether it is real or imaginary and its sign. This will locate which of these quadrants is offset by the angle to which the value belongs. 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 real number line, as shown by the center panel 3320 of fig. 33. FIG. 33 shows the exponential error plotted for the input of the exponential-to-logarithmic process in the region [ -2.0, +2.0), where the right side is positive real and the lower part of the image is positive imaginary. With 50 decimal places and a similar number of iterations, the black area indicates that convergence to greater than 2 has occurred-45The value of accuracy of (c). The left panel 3310 shows the convergence domain, where gray is mapped to greater than 2-45The error of (2). The middle panel 3320 shows quadrants centered on the real number line using inverted gray scale. The right panel 3330 shows the final desired convergence area for white

Figure BDA0002163150260000511

The imaginary part of the logarithmic output is preloaded with the equivalent imaginary integer part of the output complex logarithm. Due to the 45 ° offset, further testing will be necessary for complex branch intersections to fall on the negative part of the real line, but this can be neglected if 2's complement system wrap-around is used, since the imaginary part of the logarithm measured in the quadrant can easily be constrained to lie in the interval [ -2.0, 2.0) in view of the necessary bits.

At this time, when the determination value is located in the first quadrant, the absolute value of the known real part is greater than or equal to the imaginary partThe absolute value of the part, and the real part is known to be positive. To find the integer part of the real part of the logarithm, we can count the leading zeros in the fixed-point representation. The value of subtracting this count from the number of bits in the fixed-point representation is actually

Figure BDA0002163150260000512

Where k depends on the representation of the bit and is thus a representation of the integer part of the binary logarithm, which can be added directly to the real part of the logarithm register. Then, the existing data (both real and imaginary) can be shifted up by the number of bit positions required to cause the real part to lie in the interval [ +0.5, +1.0), as shown in the right panel of fig. 33, which is determined by the initial count leading zero and k depending on the representation of the bit. It should be noted that while this may not be the closest integer representation to output the real number logarithm, it is an integer and brings the remainder of the process into the converged region in the correct format.

At the end of the process, the integer part of the output complex logarithm has been determined, with the fractional part of the output and the input range not determined

Figure BDA0002163150260000513

Which fits within the convergence range of the principal logarithm. A discussion of the main algorithm for the substantially fractional part of the logarithmic output will be described below.

The result is a functionIs now the process modeled, as can be seen from fig. 34. In fig. 34, the exponential 3410 and logarithmic process 3430 domains of exponential error are plotted for only the sufficiently fractional portion, and the corresponding values that constitute the reduction in the range of final composite values of the exponential 3420 and logarithmic 3440 processes. The noisy gray represents the region where proper convergence results in a small numerical error dominating the error measurement, showing that convergence has been achieved anywhere in the exponent 3420 and the logarithm 3440.

Changing the base of the result still means that it has all the properties of logarithm and exponent, meaning that the real parts of logarithms can be added and subtracted to provide multiplication and division of exponents. Likewise, the angles appear the same, except that they are represented as integer multiples of a quadrant, so 4 times each rotation, a value that is easy to manipulate in hardware logic. In their logarithmic form, the addition and subtraction of imaginary values is the same when rotated clockwise and counterclockwise around zero in complex-valued exponential space. For complex values, taking the real and imaginary parts of their combination regardless of the fact that they consist of different base numbers, adding and subtracting the real and imaginary parts can be used to implement complex multiplication and division.

These results can be extended mathematically simply to floating point representations, since in these formats the required pre-processing has been achieved mainly by means of mantissa and exponent forms. If the base e values are absolutely required and most of the time they are not, then the real part of the output logarithm may be multiplied by ln2 and the imaginary part by

Figure BDA0002163150260000521

Likewise, to employ a base e logarithmic data format, the value may be reciprocal before the initial range reduction to take the exponent using the method

Figure BDA0002163150260000522

And converted by component multiplication.

III.General procedure

The new algorithm proceeds in a similar manner to the BKM technique to find or use the fractional part of the logarithm. There is a set of values dnE.g. {0,1, i, -1, -i, 1+ i, 1-i, -1+ i, -1-i }, from which d is selected at each algorithm step nn. Which d to use by comparing the truncated number with the set bit countnAlternatively, the set bit count may be a subset of bits from the input register. Once d is selectednA look-up table of pre-computed logarithms is required to balance 1+ dn2-nMultiplication (bit-shifting (potentially applying rounding to the most significant unused bits) and addition/subtraction operations) applied to a holding exponentA register of values, wherein the equivalent logarithm is to be subtracted from the register holding the logarithm value. Depending on whether the process finds a logarithm or computes an exponent, there may be additional shifts, additional operations, or both, in order to place the system into a converged configuration. At the end of some N iterations, the desired value (exponential or logarithmic) is found to be within a precision that depends on the number of iterations N and the number of bits used for the representation.

IV.Lookup table of logarithms

The proposed method differs from the prior art in that the logarithms must now be reflected with an exponent of 2r(eπ/2). Thus, for N-1, …, N, in the desired fixed-point format, the real part of the log table must be preloaded with log2||1+dn2-n| | and the imaginary part must be preloaded with 2/π arg (1+ d)n2-n) Or

Figure BDA0002163150260000531

Although we have selected 2 and eπ/2As a logarithmic base, but it is likely that different base selections may provide an algorithm for a logarithmic system that does not have a different first iteration. These values have been chosen because they represent a design where both processes can be implemented using the same hardware. This may not be the best base choice if separate implementations using different hardware are to be used for logarithmic mode, such as from exponential mode.

V.Method for logarithmic fractional part

The core methods of the new algorithm mainly follow those used by the BKM process, with some new key differences. The comparison thresholds have been re-derived to accommodate the new bases of logarithms and exponents and the varying convergence fields required by the new technique. For the first iteration of the logarithmic process, there is also a slightly modified comparison threshold to achieve convergence within the desired domain.

VI.Exponential mode

The exponent mode described in this section may also be modified to provide a complex with the exponent valueAnd (4) number multiplication. If this is to be done, this must be preloaded before the range reduction step in case the output is correct. It should also be noted that this will replace the output exponent value and therefore should not be used if this value is required. Store and wait to separate the solution from the integer part of the logarithm (value Z)integerOutput) is also applicable to the end of the process. This may reduce the storage required for intermediate registers where processing occurs, although this should be balanced against the additional storage requirements required for the integer part of the solution.

Alternatively, multiplication with the final exponent value may be accomplished in parallel by creating additional registers and ensuring that equivalent operations occur in these additional registers. In this manner, the exponential process can complex multiply the output exponent value by almost any number of other complex values.

Assuming that the fractional part of the input logarithm is input, for the convergence field ZinputE.r [ -0.5, +0.5) + i [ -0.5, +0.5) is:

1. assume that there are four basic registers, labeledAnd

Figure BDA0002163150260000542

in addition to this, there are two additional dependent multiply registers

Figure BDA0002163150260000543

And

Figure BDA0002163150260000544

to demonstrate how the method operates when used to assist in complex multiplication. The initial values of these registers are:

Figure BDA0002163150260000545

Figure BDA0002163150260000546

Figure BDA0002163150260000547

Figure BDA0002163150260000548

wherein Zpremultiply: 1.0 if no pre-multiplication requirement exists. The slave multiply register may also be similarly configured to:

Figure BDA0002163150260000549

Figure BDA00021631502600005410

2. iteration is performed with the value 1, …, N as index N:

3. cut off

Figure BDA00021631502600005411

To form

Figure BDA00021631502600005412

So that it has three bits; one sign bit and two integer bits in the 2's complement, such that the range is [ -4.0, +4.0), with a minimum variance of 1.

4. Cut off

Figure BDA00021631502600005413

To form

Figure BDA00021631502600005414

So that it has three bits; one sign bit, one integer bit and one decimal place in the complement of 2, so that the range is [ -2.0, +2.0), with a minimum variation of 0.5.

5. Testing 3-bit values to determine dn

Figure BDA0002163150260000551

Figure BDA0002163150260000552

6. The multiplication is applied to the exponent register using a shift sum operation to compute:

and:

for any auxiliary register (such as

Figure BDA0002163150260000555

And

Figure BDA0002163150260000556

) The same operation is performed to apply the multiplication process to these registers as well.

7. Subtract the corresponding element in the log table from the register:

Figure BDA0002163150260000561

Figure BDA0002163150260000562

8. the logarithmic register is shifted to the left by one position, doubling the value of its contents.

Figure BDA0002163150260000563

9. Return to step 2 of the next iteration until N is reached, at which point the registers will contain their final values:

Figure BDA0002163150260000565

Figure BDA0002163150260000566

and:

it is understood that the form of the process is formed by using different numbers of bits or different comparison values

Figure BDA0002163150260000569

Or both, it is easy to find other test procedures that are even sometimes convergent in the desired domain, although we have tried to reduce complexity by specifying the desired value test in the simplest known form.

VII.Logarithmic mode

The logarithmic pattern described in this section may also be modified to provide complex division with the input value. If this is to be done, it can only be done in parallel by creating additional registers and ensuring that equivalent operations occur in these additional registers. This is in contrast to auxiliary multiplication, where the original register may be overloaded, which cannot be achieved here, because modifications to the exponent register in this mode would prevent convergence of the algorithm. However, this can be circumvented using an auxiliary register, allowing the logarithmic process to produce complex-valued divisions of almost any number of other complex values with the value input to the process as the denominator, if desired.

1. Assuming that the fractional part of the output logarithm is output for the convergence domain

Figure BDA0002163150260000571

The algorithm of (1) is as follows:

2. assume that there are four basic registers, labeled

Figure BDA0002163150260000572

And

Figure BDA0002163150260000573

in addition to this, there are two additional dependent divide registers

Figure BDA0002163150260000574

And

Figure BDA0002163150260000575

to demonstrate how the method operates when used to assist complex division. The initial values of these registers are:

Figure BDA0002163150260000577

Figure BDA0002163150260000578

Figure BDA0002163150260000579

the slave divide register can also be similarly constructed as:

Figure BDA00021631502600005710

Figure BDA00021631502600005711

it is noted that addNeither multiple nor-1.0 is applied to registers

Figure BDA00021631502600005712

And

Figure BDA00021631502600005713

2. iteration is performed with the value 1, …, N as index N:

3. cut off

Figure BDA00021631502600005714

To form

Figure BDA00021631502600005715

So that it has five bits; one sign bit, one integer bit and three decimal bits in the 2's complement, such that the range is [ -2.0, +2.0), with a minimum variance of 0.125.

4. Cut off

Figure BDA00021631502600005716

To form

Figure BDA00021631502600005717

So that it has five bits; one sign bit, two integer bits and two decimal bits in the 2's complement, such that the range is [ -4.0, +4.0), with a minimum variance of 0.25.

5. Testing two 5-bit values (where modified to compare for the first iteration) to determine dn

Figure BDA0002163150260000581

Figure BDA0002163150260000582

6. The multiplication is applied to the exponent register using a shift sum operation to compute:

Figure BDA0002163150260000583

and:

Figure BDA0002163150260000584

for any auxiliary register (such as

Figure BDA0002163150260000585

And

Figure BDA0002163150260000586

) The same operation is performed to apply a division process to these registers.

7. Will dnIs added to

Figure BDA0002163150260000587

And

Figure BDA0002163150260000588

Figure BDA0002163150260000589

Figure BDA0002163150260000591

do not add this to

Figure BDA0002163150260000592

And

Figure BDA0002163150260000593

8. the exponent register is shifted one position to the left, doubling the value of its contents.

Figure BDA0002163150260000594

Figure BDA0002163150260000595

Need not to change

Figure BDA0002163150260000596

And

Figure BDA0002163150260000597

the value of (a) is doubled.

9. Subtract the corresponding element in the log table from the register:

Figure BDA0002163150260000598

Figure BDA0002163150260000599

10. return to step 2 of the next iteration until N is reached, at which point the registers will contain their final values:

Figure BDA00021631502600005911

and:

Figure BDA00021631502600005912

Figure BDA00021631502600005913

it is understood that the form of the process is formed by using different numbers of bits or different comparison valuesOr both, it is easy to find other test procedures that are even sometimes convergent in the desired domain, although we have tested by specifying the desired values in the simplest known formEfforts are made to reduce complexity.

VIII.Redundancy number system

Like the BKM approach, this algorithm is ideal for supporting techniques that slow down local data flow, allowing integration into high-speed logic implemented in hardware. The use of a redundant number system scheme (such as carry preservation of arithmetic steps coupled with a reduced number of bits for comparison at each iteration) allows for an implementation of an algorithm such as that shown herein with a high speed, low footprint, which is ideal for inclusion in an integrated circuit design. Benefits of such an implementation would include the ability to achieve higher clock speeds while scaling operands with larger numbers of bits well, making the design ideal for producing low-cost, modern solutions to problems involving manipulation of complex values.

For ease of visualization, rounding has been omitted from the foregoing algorithmic description and should not be construed as meaning that proper rounding is not considered in the present invention. It should be clear that the precision of the results described herein can be increased when applied at each step by using the most significant bits that are not directly used in any arithmetic operation to achieve a result with reduced error through the application of rounding.

IX.Higher radix method

In a similar manner to other algorithms that evaluate on a per-bit basis, the steps can be combined to create a method of computing multiple bits at once. This is referred to as a higher radix algorithm because the output is seen as having a base different from 2, meaning that one digit from the new higher base is calculated per step. The disadvantage of this is that the look-up table and shift choices will also have to be combined, creating a larger set of choices for the algorithm to complete each iteration.

X.Real number algorithm

Many systems may benefit from using only the real part of the algorithm described herein. For example, a sigmoid logic function typically used by neural network algorithms may be implemented by first determining a real number 2xOr if pre-multiplying by exIs generated using an exponential process. Then, add 1 and preset constantDividing by this value is possible, or continuing to use the E-mode and then the L-mode described by the original paper can produce the following logical function:

Figure BDA0002163150260000601

but the molecule k and the further modifiers a and b are still configurable.

XI.Display convergence

In fig. 35, the convergence of the exponential method and the logarithmic method is illustrated by the graph 3500. FIG. 35 shows that testing each possible value in the necessary convergence domain yields the worst case log2 error plotted here. The graph summarizes worst-case behavior from over one trillion complex test values (each real and imaginary value in the remapped decimal domain has been tested for the number of iterations shown). To make the graphs comparable, the logarithmic result is converted to its exponential form before testing for errors. This ensures that the method is convergent for so many significant bits. This error can be further reduced by adding full rounding in the multiplication. The graph shows that as the accuracy and number of iterations increase, the error similarly decreases, showing convergence.

The graph 3500 shows the worst case absolute error from all possible cases in the convergence domain that contract at a similar rate as the number of bits in the table and iteration in the process. This reveals that the algorithm successfully converges linearly to the correct solution with bit precision for all possible correct inputs of that number of bits of the setpoint data.

XII.The latent claims

1. A method for simultaneously converting an input 'logarithmic' complex number into a resulting 'exponential' complex number, wherein:

a. the modulus of the resulting complex number consists of a real power of 2, which is the real power of the real input complex number part, and;

b. the argument of the resulting complex number consists of an integer power rotation of 2 multiplied by the real value of the imaginary input complex number part.

2. A method for simultaneously converting an input 'exponential' complex number into a resultant 'logarithmic' complex number, wherein:

a. inputting a real binary logarithm of a real modulus of the complex number forming a real part of the resulting complex number, and;

b. the integer power rotation of a real number of 2 in terms of the argument of the input complex number forms the imaginary part of the resulting complex number.

3. The method of claim 1 or 2, wherein the conversion is substantially achieved by bit shifting, comparing and operating on a binary representation of a complex number.

4. The method of claim 1, wherein the reduction of the input and output ranges is performed by processing rounded values corresponding substantially to the components of the input complex number separately into further sections.

5. The method of claim 2, wherein the reduction of the input and output ranges is performed by processing rounded values of components substantially corresponding to the output complex number separately into further sections.

6. The method of claim 1, wherein the output complex value can be initialized with a value not equal to 1, such that application of an iteration will modify the non-uniform value to compute a complex multiplication from the value originally intended to be output.

7. A method according to claim 1, wherein further complex values are instantiated such that operations applied to output complex values can also be applied to these further values, the multiplication of these further complex values resulting from the complex values originally intended to be output.

8. A method according to claim 2, wherein further complex values are instantiated such that operations applied to input complex values can also be applied to these further values, the division of these further complex values resulting from the input complex values.

9. The method of claim 1, wherein a redundant number system is used to implement a process of repeated operations applied to existing complex values.

10. The method of claim 2, wherein a redundant number system is used to implement the process of repeated operations applied to existing complex values.

11. An apparatus designed to perform any of the methods of claims 1-10.

12. A data processor configured to perform any of the methods of claims 1 to 10.

7.Sound field constructed from Krylov subspace of quadratic form problem

The present disclosure describes a method of manipulating quadratic problems when represented as a eigensystem to produce reduced rank basis sets for inclusion in a linear system. By using a reduced rank basis, an optimal information-preserving set (potentially reduced) of basis functions can be automatically derived for inclusion in a complex-valued linear system to find the optimal complex values of such basis functions to re-weight to achieve a set of goals described by a number of test functions.

The derivation of the reduced set of basis functions and how they are included in the appropriate phase of the implementation of the field linearity optimization problem will be described.

I.Linear control region

A linear control region, which is defined as an operator that performs integration with a linear acoustic quantity such as pressure, has a complex-valued weighting function, w (χ), χ ═ x, y, z }:

Figure BDA0002163150260000631

may be used to evaluate the sound field. It should be noted that when w (χ) is paired with a similar phase distribution in the acoustic quantity, the complex conjugate is used for the maximization function evaluation. Any approximation of this integral (such as sampling) may be used, or indeed this may be used to approximate any linear operator on the acoustic pressure field. Additional integrations can be constructed to optimize the control region using vector particle velocities (where each vector component is complex-valued due to harmonic motion of the field), as long as it is also linear, where:

Figure BDA0002163150260000632

where u is the particle velocity of the medium and w is a three-dimensional vector field, each dimension being a complex value, which is used as a weighting function with the same purpose as above. In this way, different acoustic magnitude patterns can be evaluated in a linear fashion, where a particular phase of the result will maximize the objective function when the phase is aligned with a 'template' created by a scalar or vector function in w.

It should also be noted that w (χ) may be mostly zero (or consist of a zero vector) to create a single evaluation of multiple seemingly disjoint regions. Any integral approximation (such as sampling) can be used, as the formula can be extended to any linear operator on the field of linear acoustic quantities.

II.Quadratic form problem and characteristic problem

A phased array can be described by controlling a series of complex valued coefficients for each transducer. Driving an acoustic transducer array at a fixed frequency with such an input x-vector to maximize the acoustic quantity encoded in the matrix M can be represented as a linear algebraic problem, so the intention is:

maximizing xHMx

Is subject to xHx=1,

And is

Figure BDA0002163150260000633

This can be solved by taking the principal eigenvector of the matrix M, since the statement of the problem is also the definition of the eigenvector (here x). The construction of the matrix M can be achieved by first considering the linear acoustic quantities emitted by each transducer as basis functions in the global space of the acoustic quantities in question. This is the way in which the maximization of pressure and other acoustic quantities is achieved at a single point described by the absolute position χ. However, using the described control region integration, q (χ) controlling the acoustic mass in the field launched into each individual transducer element can be replaced by f (χ), resulting in a column vector:

Figure BDA0002163150260000641

by invoking linearity, the complex activation values that drive each transducer can be multiplied, resulting in a vector:

Figure BDA0002163150260000642

where n is the number of sensors. Given that m is the result of a linear operator operating on a linear field, it can be hermitian transposed and multiplied by itself to generate a quadratic form. Finally, the process is carried out in a batch,

mHm=xHMx=xHffHx。

by modifying the weighting function w (χ) and solving the eigenvectors of the matrix M, various beams, surfaces, and other acoustic phenomena can be created and maximized, resulting in regions of optimal actuation in the acoustic field.

The previously discussed null-space generation technique ("eigensystem solution to the quadratic problem in the semi-empty haptic system") can also be used to extend the method to describe acoustic levitation and allow null steering and formation of semi-empty textures by introducing a null into the field at the location of the described structure. It should be noted that this procedure increases the rank of the matrix, so simplification is not possible if this is used.

Likewise, for many control regions, multiple optimizations may be performed by adding the M matrices together. This may also include a matrix representing the control points. This generates a matrix with an increased rank, so some simplified methods are not possible using this.

III.Reduced rank approximation to high rank matrix

The key point is that a high rank matrix can be intentionally created so that a reduced rank approximation can be generated. For example, a system with hundreds of individual points may be automatically re-grouped so that very few points are described by independent basis functions, creating a real number of basis groups that appear to be more dimensional than they originally were, to increase fidelity, reduce complexity and allow basis functions to be automatically created as needed.

Existing methods generate a new basis function for each control point or control region. The feature system rank is increased by creating a feature system that optimizes each control point or control region and sums or adds nulls to them. By repeatedly adding basis function matrices and/or nulls, the final matrix may obtain a high rank. Due to resource constraints, a linear system may only be able to process up to n basis functions using up to m test functions. The high rank of the optimization matrix may in this case be n', where the n "optimization matrices used to generate the individual basis functions have been merged together by addition. Using this approach, other unlikely cases where n' > n and/or n "> n can be reduced to a problem represented in a simplified basis set with n basis functions that can be solved with the system with limited processing power.

This simplified basis set is found by computing a Krylov subspace that occupies substantially the same matrix subspace as the n most dominant eigenvectors of the final matrix composed of multiple optimizations of the basis functions. The Krylov subspace can be found using standard techniques such as the Arnoldi iteration, where the dominant feature vector is found and then removed from the remaining matrix and then applied in an iterative fashion until n such feature vectors are found.

IV.Mixing method

This method of extracting basis functions from a high rank matrix can thus be mixed and matched with existing specific basis function methods. For example, two control regions or points may be designated as having to have (must-have) and added to the base set input to the linear system. The remaining basis functions (assuming they can be specified by quadratic optimization and thus their own optimization matrix) can then be added to the single matrix, which is then deconstructed using this Krylov subspace approach to find the most dominant n-2 eigenvectors to form the remaining number of sets of basis functions input to the linear system. In this way, the pre-specified two basis functions are guaranteed, since they are directly imported into the linear system, but the remaining functions of the set are approximated. The system may even further be employed wherein a set number of basis functions may be specified and retrieved from a set of disjoint basis functions, each described by a high rank optimization matrix and then the set number of basis sets retrieved. Additional techniques that may be used to derive the reduced basis functions continue to take the basis functions from the matrix until the power threshold is crossed. The threshold may be different for each matrix of intersecting basis sets and thus produce a different power of feature vectors, where the power of each feature vector is described by its corresponding feature value magnitude, and in this way allows a feature vector to be produced by selecting a feature vector calculated from each matrix based on the product of the real number priority value and the corresponding feature value magnitude.

8.Linear system solution using test functions

A method of 'testing' the set of basis functions to achieve some goal defined separately from the set of basis functions would be valuable. A quantitative assessment of whether these objectives are achieved is provided to the optimization and this is used to generate the desired properties defined by these test functions in the sound field.

The test function and its role in complex-valued linear system optimization will be described.

I.Linear control region

A linear control region, which is defined as an operator that performs integration with a linear acoustic quantity such as pressure, has a complex-valued weighting function, w (χ), χ ═ x, y, z }:

Figure BDA0002163150260000661

may be used to evaluate the sound field. It should be noted that when w (χ) is paired with a similar phase distribution in the acoustic quantity, the complex conjugate is used for the maximization function evaluation. Any approximation of this integral (such as sampling) may be used, or indeed this may be used to approximate any linear operator on the acoustic pressure field. Additional integrations can be constructed to optimize the control region using vector particle velocities (where each vector component is complex-valued due to harmonic motion of the field), as long as it is also linear, where:

where u is the particle velocity of the medium and w is a three-dimensional vector field, each dimension being a complex value, which is used as a weighting function with the same purpose as above. In this way, different acoustic magnitude patterns can be evaluated in a linear fashion, where a particular phase of the result will maximize the objective function when the phase is aligned with a 'template' created by a scalar or vector function in w.

It should also be noted that w (χ) may be mostly zero (or consist of a zero vector) to create a single evaluation of multiple seemingly disjoint regions. Any integral approximation (such as sampling) can be used, as the formula can be extended to any linear operator on the field of linear acoustic quantities.

II.Set of basis functions

Given a set of n basis functions as inputs, an optimal linear combination of n basis functions can be found by a linear system using m test functions to achieve a desired target evaluation of the m test functions. These m test functions may take the form of linear control region evaluations. Using the form of target f as given above, these would be:

Figure BDA0002163150260000671

or they may be evaluated by a single control point definition as follows:

where the upper dashed line represents the complex conjugate, χ is the three-dimensional position, q parameterizes all transducer elements, and x is the transducer activation, which is multiplied by the conjugate of the transducer's pressure field. An equivalent form of the control point defined by the particle velocity may also be used as an estimator, although with a defined normal direction:

Figure BDA0002163150260000673

in which any transducer q pointing in the wrong direction can also be omitted (when

Figure BDA0002163150260000674

Its contribution is set to zero).

In the two single control point definitions, the existing method can be seen as a special single point integral evaluation with a unit weighting function case of the control area evaluation. The objective function here is then a simple acoustic quantity with no additional side effects.

III.Finding dominant feature vectors

Each of the m basis functions is evaluated by using each of the n test functions, resulting in a rectangular matrix that is complex-valued. As before, the dominant feature vector represents a combination of m basis functions that best span the space defined by the n test functions. However, since this is not necessarily a square matrix, singular value decomposition is more appropriate, but the only required feature vector is the dominant one. A square hermitian matrix may be constructed whose principal eigenvectors represent the elements of the principal combination of system basis functions. The characterization system is then described by:

Figure BDA0002163150260000681

wherein b is1,…,bnCoefficients describing the coverage of the parameter space whose phase generates the m test functions. Thus, the dominant feature vector has a component b1,…,bnIt describes the optimal coverage of the space achievable with a single combined vector, and its phase thus optimizes this coverage.

It should be noted that even in the case where the number of test functions is equal to the number of basis functions (so m ═ n), the existing signature system method may not necessarily be applicable. This is due to the fact that: the matrix may no longer be a hermitian matrix and thus applying the Rayleigh-Ritz method directly to the matrix without completing the above steps may result in the method failing to converge to the dominant eigenvector.

If no phase data is specified or the phase data is incomplete, the dominant eigenvector should be used for phase information. At the end of each iteration, the feature vector may be modified to strip the amplitude data (which averages zero in phase, is finite in phase progression) from the previous solution at a separate time step, or to have a set of phase solutions for some basis functions that hold an existing and constant phase relationship, with the complement of the set being optimized with the phase between itself and a fixed set of existing phases. In this way, the rectangular matrix system can be used to achieve an effect similar to existing eigensystem methods, which are dedicated to the case where the function generating the basis set is the same as the function used to evaluate the objective of the basis set.

IV.Linear system solution

Complex-valued linear systems now have the modified effect of finding the set of basis functions that best replicate the given values of the test function. The linear system that must be used here takes the form of a rectangular matrix:

since the matrix is rectangular, it can be a least squares system with a number of test functions greater than the number of basis functions, resulting in complex values, or a minimum norm system with a number of test functions less than the number of basis functions, resulting in complex values. In either case, a solution vector may be found using a standard decomposition, the interpretation of which is that each component of the solution vector is a coefficient to be applied to each basis function that causes the sound field, when tested, to express a given value of the test function.

These test functions may represent acoustic 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 can then constitute m test functions, which are solved for the final step of the linear system.

Once the coefficients of the basis functions are known, the final value of the actuation transducer can be extended to:

Figure BDA0002163150260000692

an input signal that gives the assumption of a single monochromatic carrier frequency (which may be more or less effective due to environmental changes) may be constructed by applying this x component to each transducer.

9.Conclusion

The various features of the foregoing embodiments may be selected and combined to produce many variations of the improved haptic system.

In the foregoing specification, specific embodiments have been described. However, one of ordinary skill in the art appreciates that various modifications and changes can be made without departing from the scope of the present invention as set forth in the claims below. Accordingly, the specification and figures 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 present teachings.

The benefits, advantages, solutions to problems, and any element(s) that may cause any benefit, advantage, or solution to occur or become more pronounced are not to be construed as a critical, required, or essential feature or element of any or all the claims. The invention is defined solely by the appended claims including any amendments made during the pendency of this application and all equivalents of those claims as issued.

Moreover, in this document, relational terms such as first and second, top and bottom, and the like may be used solely to distinguish one entity or action from another entity or action without necessarily requiring or implying any actual such relationship or order between such entities or actions. The terms "comprises," "comprising," "includes," "having," "has," "having," "includes," "including," "contains," "containing," "contains" or any other variation thereof, are intended to cover a non-exclusive inclusion, such that a process, method, article, or apparatus that comprises, contains a list of elements does not include only those elements but may include other elements not expressly listed or inherent to such process, method, article, or apparatus. An element followed by "comprising … … a," "having … … a," "including … … a," "including … … a," does not exclude, without limitation, the presence of additional identical elements in the process, method, article, or apparatus that comprises, has, contains, or contains the element. The terms "a" and "an" are defined as one or more unless the context clearly dictates otherwise. The terms "substantially", "essentially", "approximately", "about" or any other version thereof are defined as being proximate as 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 certain way is configured in at least that way, but may also be configured in ways that are not listed.

The Abstract of the disclosure is provided to allow the reader to quickly ascertain the nature of the technical disclosure. The abstract is submitted with the understanding that it will not be used to interpret or limit the scope or meaning of the claims. In addition, in the foregoing detailed description, it can be seen that the various features are grouped together in various embodiments for the purpose of streamlining the disclosure. This method of disclosure is not to be interpreted as reflecting an intention that the claimed embodiments require more features than are expressly recited in each claim. Rather, as the following claims reflect, inventive subject matter lies in less than all features of a single disclosed embodiment. Thus the following claims are hereby incorporated into the detailed description, with each claim standing on its own as a separately claimed subject matter.

71页详细技术资料下载
上一篇:一种医用注射器针头装配设备
下一篇:检测和抑制话音查询

网友询问留言

已有0条留言

还没有人留言评论。精彩留言会获得点赞!

精彩留言,会给你点赞!