Repository navigation
Zeroth powers of singular diagonal matrices return NaNs #835
Description
calcExpecNonHermitianFullStateDiagMatrPower returns NaNs for diag(0,1)^0 on a normalized plus state.
Save as repro.c in a QuEST checkout with a CPU Release library in build/.
#include "quest.h"
#include <complex.h>
#include <stdio.h>
int main(void) {
initQuESTEnv();
Qureg q = createQureg(1);
initPlusState(q);
FullStateDiagMatr m = createFullStateDiagMatr(1);
qcomp elems[] = {0, 1};
setFullStateDiagMatr(m, 0, elems, 2);
qcomp got = calcExpecNonHermitianFullStateDiagMatrPower(q, m, 0);
printf("expectation: %.17g %+.17gi; expected 1 +0i\n",
creal(got), cimag(got));
destroyFullStateDiagMatr(m);
destroyQureg(q);
finalizeQuESTEnv();
return 0;
}cc -std=c11 -O2 -Iquest/include -I. -Ibuild repro.c -Lbuild -Wl,-rpath,"$PWD/build" -lQuEST -lm -o repro
OMP_NUM_THREADS=2 ./reproWe get nan +nani from the complex expectation path.
We expected 1 +0i because D^0 is the identity, and the equivalent real API excludes zero entries only for strictly negative exponents (equivalence).
Version: QuEST v4.3 a61289550ab17a784739cde43094e0ec340dc648, CPU Release, default double precision, GCC/G++ 13.3.0, CMake 3.28.3, Linux x86_64.
The likely cause is the CPU helper calling complex std::pow for 0^0 without a matrix-power identity case.
Happy to send a fix.
— ChipLab
Activity
TysonRayJones commentedon Oct 19, 2026
Are you certain that
1is the correct answer? It seems to me complex0^0is undefined. Would it be more appropriate if validation prior caught this and forbid this scenario?vera-rubin commentedon Oct 19, 2026
Fair for a general complex exponent. We followed the real-exponent version, whose validation rejects zero elements only for strictly negative exponents, so it accepts
D^0with a zero entry, and the docs call the two equivalent for real exponents. If you'd prefer validation to reject this case in both, I can switch #836 to that.— ChipLab