Why the NullSpace can not find solution when the Det[t]==0

Question

This question contained the problem of NullSpace, but previous ones not.

This is a problem seems like my previous one, but there is some details different. I find that perhaps one method can not solve all the cases with different parameter. And thus, this question is not a duplicate of another question.

t is a generated matrix containg kz. kz are obtained by solve Det[t]==0. I find that when s and m are set as some cases (e.g. s=5,m=2), I substitue the NSolve[Numerator[Together[N@Det[Chop[t]]]] == 0, kz] into t, of course Det[t]==0, this means there must exist NullSpace[t]. However, NullSpace[t] gives a {},Why? How to solve this problem. Many Thanks! The codes are as following:

Clear["`*"]
s = 5;
m = 2;
th = Pi/4;
fi = Pi/6;
vh = 16;
mu = 11;
delta = 8;
HBAR = SetPrecision[1.05457266*10^(-34), 100];
ME = SetPrecision[9.1093897*10^(-31), 100];
ELEC = SetPrecision[1.60217733*10^(-19), 100];
Kh = Rationalize[0.211, 0];
vKh[1] := {0, 0, 0}
vKh[2] := {Kh, 0, 0};
vKh[3] := {-Kh, 0, 0};
vKh[4] := {0, Kh, 0};
vKh[5] := {0, -Kh, 0};
vKh[0] := {0, -Kh, 0};
vKh[i_] := vKh[Mod[i, 5]];
kc = Sqrt[2*ME*ELEC/HBAR^2]*10^(-11);
ku := kc*Sqrt[mu + delta];
kd := kc*Sqrt[mu - delta];
a3 = {Pi/Kh, Pi/Kh, Sqrt[2]*Pi/Kh};
k := {-ku*Sin[th]*Cos[fi], -ku*Sin[th]*Sin[fi], kz};
vkz[i_] := 
  If[Mod[i, 5] != 0, {0, 0, (i - Mod[i, 5])/5*Kh*Sqrt[2]/(m + 1)}, {0,
     0, (i - Mod[i, 5] - 5)/5*Kh*Sqrt[2]/(m + 1)}];
f[i_, i_] := 
  Total[(k + vKh[i])^2] - ku^2 - kz^2 + (kz + Total[vkz[i]])^2;
f[i_, j_] := 
  If[i == j, f[i, i], 
   kc^2*vh*Total[
     Table[Exp[I*n*Total[(vKh[j] + vkz[j] - vKh[i] + vkz[i])*a3]], {n,
        0, m}]]];
t := Array[f, {5*s, 5*s}];
slu := Select[
   kz /. NSolve[Numerator[Together[N@Det[Chop[t]]]] == 0, kz], 
   Re[#] >= 0 && Im[#] >= 0 &];
td[i_] := t /. kz -> slu[[i]];
nu[i_] := NullSpace[Chop[td[i]]];

Det[Chop[td[1]]]
nu[1]

Answer 1

First, I made some changes, including making everything exact:

HBAR = Rationalize[1.05457266*10^(-34), 0];
ME = Rationalize[9.1093897*10^(-31), 0];
ELEC = Rationalize[1.60217733*10^(-19), 0];

t = ExpToTrig@Array[f, {5*s, 5*s}];  (* ExpToTrig eliminates N::meprec warning *)

A warning arises when (1 + E^(-((2 I π)/3)) + E^((2 I π)/3)) is numerically evaluated at arbitrary precision. Symbolically it is zero, and ExpToTrig is one way to convert it to zero (and much faster than Simplify[t]).

Next, the following finds the nullspace:

wp = 100;
nt = N[t, wp];                             (* numericize matrix  t  *)
nd = Det[nt]; // AbsoluteTiming            (* symbolic det. with numeric coeffs. *)
With[{obj = Numerator@Together@nd},
   sol = NSolve[obj == 0 && 0 <= Re@kz && 0 <= Im@kz, kz, 
     WorkingPrecision -> MachinePrecision] (* or WorkingPrecision->Precision@obj *)
   ]; // AbsoluteTiming
sol = If[Precision@sol === MachinePrecision,
  Select[sol, Abs@First@Ratios@Eigenvalues[nt /. #][[{1, -1}]] < 10^-13 &],
  Select[sol, Precision[#] > 10 &]
  ]
NullSpace[nt /. #] & /@ sol
(*
  {0.409474, Null}              <-- timing of Det[nt]

  {0.132046, Null}              <-- timing of NSolve[]

  {{kz -> 0.105997 + 0.00517193 I}, {kz -> 0.414762 + 0.178031 I}} <-- sol

  {{{-0.0525609 + 0.0371921 I,..., -0.00461674 - 0.00426698 I}},   <-- null spaces
   {{-0.631073  + 0.195673 I,...,  -0.056275   + 0.0518164 I}}}    <--  "     "
*)

NSolve might (and does) produce spurious solutions. I used Select to eliminate them. If using MachinePrecision, the spurious ones won't have a (nearly) zero eigenvalue. If using arbitrary precision, the only spurious ones I found had an extremely low precision compared to Precision@obj, which is around 90.

Answer 2

This is basically due to @MichaelE2, but I am making a few changes which might remove potential trouble spots.

One change is to set precision high, then rationalize. It might give the same results as rationalizing with second arg set to 0, I didn't check. If not, it might be a hair more faithful to the decimals. Probably makes no difference in the long run either way.

s = 5;
m = 2;
th = Pi/4;
fi = Pi/6;
vh = 16;
mu = 11;
delta = 8;
HBAR = Rationalize[SetPrecision[1.05457266*10^(-34), 100]];
ME = Rationalize[SetPrecision[9.1093897*10^(-31), 100]];
ELEC = Rationalize[SetPrecision[1.60217733*10^(-19), 100]];
Kh = Rationalize[0.211, 0];
vKh[1] = {0, 0, 0};
vKh[2] = {Kh, 0, 0};
vKh[3] = {-Kh, 0, 0};
vKh[4] = {0, Kh, 0};
vKh[5] = {0, -Kh, 0};
vKh[0] = {0, -Kh, 0};
vKh[i_] := vKh[Mod[i, 5]];

kc = Sqrt[2*ME*ELEC/HBAR^2]*10^(-11);
ku = kc*Sqrt[mu + delta];
kd = kc*Sqrt[mu - delta];
a3 = {Pi/Kh, Pi/Kh, Sqrt[2]*Pi/Kh};
k = {-ku*Sin[th]*Cos[fi], -ku*Sin[th]*Sin[fi], kz};
vkz[i_] := 
  If[Mod[i, 5] != 0, {0, 0, (i - Mod[i, 5])/5*Kh*Sqrt[2]/(m + 1)}, {0,
     0, (i - Mod[i, 5] - 5)/5*Kh*Sqrt[2]/(m + 1)}];
f[i_, i_] := 
  Total[(k + vKh[i])^2] - ku^2 - kz^2 + (kz + Total[vkz[i]])^2;
f[i_, j_] := 
  kc^2*vh*Total[
    Table[Exp[I*n*Total[(vKh[j] + vkz[j] - vKh[i] + vkz[i])*a3]], {n, 
      0, m}]];
t = Array[f, {5*s, 5*s}];

Instead of computing Det with the matrix numericized, we'll interpolated with the variable kz numericized. The advantage is this will avoid creation of denominators which in turn allows us to bypass doing polynomial algebra (e.g. Together) on an approximate rational function. I use Quiet due to the hidden zeros issue already noted. One could instead use the ExpToTrig tactic. In a production code setting I would interpolate at points on the unit circle. The below quick-and-dirty interpolation at integers is perhaps easier to follow.

Quiet[dets = Table[Det[N[t, 200]], {kz, -25, 25}]];
det = InterpolatingPolynomial[Transpose[{Range[-25, 25], dets}], kz];

We now can extract the roots and select the ones of interest. Probably better to do at high precision but I wanted to show that machine arithmetic suffices for purposes of next computing null vectors.

ns = 
 Select[NSolve[det, kz, WorkingPrecision -> MachinePrecision], 
  With[{val = kz /. #}, Re[val] >= 0 && Im[val] >= 0] &]

(* Out[390]= {{kz -> 0.105997 + 0.00517193 I}, {kz -> 
   0.414762 + 0.178031 I}} *)

td[i_] := t /. ns[[i]];

NullSpace[td[1]]

(* Out[300]= {{-0.0491714 - 0.0415701 I, -0.00914742 + 
   0.000528511 I, -0.0373394 + 0.00133131 I, -0.0108843 + 
   0.000614003 I, -0.022615 + 0.00106742 I, -0.391521 + 0.0582159 I, 
  0.0118632 - 0.420233 I, 0.00762543 - 0.0142246 I, 
  0.028503 - 0.0612858 I, 
  0.00949882 - 0.017915 I, -0.0176812 + 0.0501965 I, -0.0550383 + 
   0.204854 I, -0.020311 + 0.0653341 I, -0.0406274 + 
   0.141211 I, -0.0234697 + 0.0763275 I, 
  0.631418 + 0.247233 I, -0.0681042 - 0.324485 I, 
  0.0121311 - 0.00123308 I, 0.0563468 - 0.0137243 I, 
  0.0153808 - 0.00171599 I, -0.0300175 + 0.00143636 I, 
  0.00417262 - 0.00790418 I, 0.002762 - 0.00512705 I, 
  0.0037664 - 0.00709248 I, 0.00297463 - 0.00553832 I}} *)