Find complex conjugate pairs
findConjugates.Rd
Find all complex conjugate pairs in a vector of complex numbers and return one number from each pair.
Arguments
- x
a vector of complex numbers
- complex.eps
a small positive number used to identify complex conjugates:
x[i]
andx[j]
are considered conjugates if(abs(x-Conj(x)) / max(abs(x[i], x[j]))) < complex.eps
and
(abs(x[i] - x[j]) > complex.eps
. The latter condition excludes repeated roots.
Details
1. Compute normalization m2 = outer(abs(x), abs(x), max)
.
2. Compute complex differences
c2 = abs(outer(x, Conj(x), "-")) / m2
.
3. If any abs(c2) < complex.eps
, make sure the numbers are not
duplicate reals via
(d2 = abs(outer(x, x, "-"))) > complex.eps
.
Examples
# none
findConjugates(NULL)
#> complex(0)
findConjugates(numeric(0))
#> complex(0)
findConjugates(0:4)
#> complex(0)
findConjugates(rep(0:1,each=3))
#> complex(0)
# some
findConjugates(c(1+1i, 0, 1-1i, 2-2i, 3, 2+2i, NA))
#> [1] 1-1i 2+2i
# Testing with polyroot and solve(polynomial(...))
set.seed(1234)
if(require(polynom)){
p <- polynomial(sample(1:10, size=45, rep=TRUE)) # degree 44
z <- solve(p)
findConjugates(z, complex.eps=.Machine$double.eps)
# this identifies all 21 conjugate pairs, R 2.6.0 for Windows
z1 <- polyroot(p)
findConjugates(z1, complex.eps=.Machine$double.eps)
# this only identifies only 3 conjugate pairs, R 2.6.0 for Windows
}
#> Loading required package: polynom
#> complex(0)