When I try to solve my question according to @DanielLichtblau's comments, I encountered another issue, with which I have been struggling for a whole night.
Taking a matrix for example,
mat[a_, b_, c_, k_] = {{Cos[a*k], f^3 + 2, -I*f + k/4}, {0, b*Sin[a*k], f}, {1/f^2, (k^2 - 1)/c, k*Cos[a*k]}};
where a
, b
, c
and k
are all real numbers, while f
is a complex variable
Defining its determinant with a given parameter c=1
det[a_, b_, k_] := Det[mat[a, b, 1, k]]
Because I need the determinant to be zero, taking its numerator should be enough. Also, I assign values to the remaining parameters (a
, b
, k
).
fexpr = Numerator[det[1, 1, 0.1]]
In view of the variable f
in the matrix being complex-valued, expanding the numerator and making explicit replacement of its real and imaginary parts
polys = ComplexExpand[ReIm[fexpr], f] /. Re[f] -> ref /. Im[f] -> imf
Solving for the real and imaginary parts of f
that make the numerator of the determinant vanishing
NSolve[polys, {imf, ref}, Reals]
Here the constraint of Reals
was used since both imf
and ref
are real-valued, we only need the real-valued solutions. But, it reports an error which appears to relate to the Reals
constraint.
NSolve::nddc: The system (-0.00249584+0. I)-0.0998334 imf-0.00988384 imf^2+1. imf^4+2. ref-2.95516 imf^2 ref+0.00988384 ref^2-6. imf^2 ref^2+0.985054 ref^3+1. ref^4==0&&(0. +0. I)+<<7>>+2.95516 imf ref^2+4. imf ref^3==0 contains a nonreal constant -0.00249584+0. I. With the domain [DoubleStruckCapitalR] specified, all constants should be real.
Strangely, if I try a slightly modified polys1
polys1 = {-0.1 - 0.1*imf - 0.01*imf^2 + 1.1*imf^4 + 2.01*ref -
3.9*imf^2 ref + 0.01*ref^2 - 6.1*imf^2 ref^2 + 0.98*ref^3 +
1.01*ref^4, 0.01 + 2.01*imf - 0.985*imf^3 + 0.02*imf ref - 4.01*imf^3 ref + 2.955*imf ref^2 + 4*imf ref^3};
It turns to work well with NSolve
NSolve[polys1, {imf, ref}, Reals]
(*{{imf -> 0.00113843, ref -> -1.68558}, {imf -> 1.05537, ref -> 0.278329}, {imf -> -1.03992, ref -> 0.310137}, {imf -> -0.00495377, ref -> 0.0494334}}*)
Can anyone please help with this. Thank you very much!
Update (Dec22)
I am trying to combine @Bob's answer into a loop but the modified code always gave many errors. I did slight modifications, which can give correct results without looping. Could you give any suggestion?
det[a_, b_, k_] = With[{a = a, b = b, k = k}, Det[mat[a, b, 1, k]]];
CurvePoints = {};
fr[k_?NumericQ, a_?NumericQ, b_?NumericQ] := Module[{sol, solrAbs, fexpr},
fexpr = Numerator[ComplexExpand[det[a, b, k], {f}]] /. f -> fre + I*fim;
sol = {f -> fre + I*fim} /. NSolve[Thread[(ReIm[fexpr] // ComplexExpand) == 0], {fre, fim}, Reals, WorkingPrecision -> 20];
solrAbs = Abs[Re@sol[[All, 1, 2]]];
Table[If[solrAbs[[n]] < 0.35,
AppendTo[CurvePoints, {k, a, b} -> sol[[n, 1, 2]]];
fr[k, a, b] = {k, a, b, sol[[n, 1, 2]]}, Sequence[]], {n, 1, Length[sol]}]]
This is the loop with a temporary monitor
pts = Monitor[Table[fr[k, a, b], {a, 0, 1, 0.5}, {b, 0, 1, 0.5}, {k, 0, 0.1, 0.1}], {k, a, b, s}];
a,b,k
? $\endgroup$Solve[fexpr == 0, f]
? $\endgroup$f
separately in my real problem. Please see my previous question and discussion therein. $\endgroup$