pairs (Frequency Domain Identification Toolbox)

Find the closest point pairs in two complex vectors.

indab = pairs(a,b)
[indab,cycle,digits] = pairs(a,b,p,maxcycle,digitsreq,D)

Description

pairs looks for a permutation of the elements of the complex vector b to minimize

or if the array D is given,

The appropriate indices of b are returned in indab. The vector a may not be longer than b. If each element of a is paired to an element of b, indab will contain positive indices only; if not (the algorithm did not converge in the allowed number of iterations), indab will contain at least one zero.

cycle gives the number of iterations of the Hungarian method. If it is returned with the value zero, the simple nearest neighbors search gave the optimum. The Hungarian method may converge very slowly in some rare cases. Therefore, the maximum number of allowed iterations may be given in maxcycle. If maxcycle is given with a finite value, the internal cost values of the Hungarian method (the complements of the powers of distances, with respect to the maximum value) will be rounded to a reasonable number of digits (convergence is assured for integer numbers). Such a rounding will be indicated by a value of digits, smaller than floor(log10(1/eps)), which is 15 on most platforms. The number of used digits can also be adjusted at the beginning of the iteration (digitsreq).

Default Argument Values

p = 2, maxcycle = inf, digitsreq = floor(log10(1/eps))

Examples

indab = pairs(roots([1,2,3,4]),roots([1.01,2,3,4]),2);
pind = pairs([-0.1-1.5*j,-0.1+1.5*j],roots([1,2,3,4]));

If cycle == 0, it may be assumed that the new positions of the roots of the perturbed polynomial p1 have been found. If cycle>0, and all(indab)>0, the optimal indices are found, but there is a chance that the roots are not paired properly, since the nearest neighbors cannot all be paired to each other.

Diagnostics

The vectors are checked for infinite or NaN elements. If a or b is empty or is an array, an error message will be sent. p is also checked for positivity.

If the closest pairs do not give the best permutation, a warning message is sent:

WARNING: nearest neighbors do not give optimum in PAIRS

The algorithm may converge very slowly. In this case the reciprocals of the distances will be rounded, and a warning message will be sent:

Number of digits is set to ... in PAIRS

The algorithm may not converge if the number of allowed iterations is small. If this happens, a warning message is sent:

WARNING: PAIRS did not converge

and indab will contain at least one element that is equal to zero.

Algorithm

The so-called Hungarian method ([1], [2]) is used: it is based on looking for alternating paths in bipartite graphs. The convergence of the method is assured for rational distances only; the algorithm implemented here makes an attempt to find a solution without rounding; when it fails, the complements of the distances to the maximum distance are rounded. By this the smallest distances, which are the most interesting, will be distorted the least.