Solve for parameter value (complex) for zero determinant

I am new to mathematica (less then a month), with the discussion available on this platform and with some tutorial available I have manage to write a code (not a smart one) to solve for the parameter values (values can be real or complex) for which the determinant is zero using Reduce function.

Objective: Search for value of 'sb' for which the det(DCoeffMat) is zero for different values of 'Pb'

(*%%% Parameters %%%%*)
ClearAll[sb, Pb, b, gama, L, w, h, Y, MI, EI, mu, rho, A, CN]
L = 10;
w = 0.02;
h = 0.02;
Y = 209*10^9;
MI = w*(h^3)/12;
EI = Y*MI;
mu = 0.001;
rho = 7800;
A = w*h;
CN = 4*Pi*mu/(Log[2*L/w] + 0.5) // N;
gama = rho*A*EI/(CN^2*L^4)
Do[B = b /. NSolve[b^4 + Pb*b^2 + sb + gama*sb^2 == 0, b] // N;
DCoeffMat = {{1, 1, 1, 1}, {B[]^3 + Pb*B[],
B[]^3 + Pb*B[], B[]^3 + Pb*B[],
B[]^3 + Pb*B[]}, {B[]*Exp[B[]], B[]*Exp[B[]],
B[]*Exp[B[]],
B[]*Exp[B[]]}, {B[]^2*Exp[B[]] + Pb*Exp[B[]],
B[]^2*Exp[B[]] + Pb*Exp[B[]],
B[]^2*Exp[B[]] + Pb*Exp[B[]],
B[]^2*Exp[B[]] + Pb*Exp[B[]]}};
DA = Det[DCoeffMat];
Print[Reduce[
DA == 0 && -1000 < Re[sb] < 1000 && -1000 < Im[sb] < 1000, sb] //
Quiet // N]
, {Pb, 0.1, 4., 0.1}]

Here is the matrix whose determinant is equated to zero On running the program its showing an error.

Am I missing something???

As per the literature from where the problem is pick up it says that for pb<2.4674 the roots of sb has small real part (with non zero imaginary part) and after that its real value increases.

• You can't use NSolve with sb as a parameter. You're missing a ; after the // N and before DCoeffMat. And it's unlikely that Reduce will be able to solve your equation. You will most likely need to use FindRoot instead. – Carl Woll Jan 26 '18 at 5:49
• Thanks for the correction. I found that NSolve can give the four roots of the polynomial equation. Can Findroot function give complex roots??? – T S Singh Jan 26 '18 at 5:58
• Sure FindRoot can give complex roots, e.g., FindRoot[x^2 - 3 I, {x, 1}]. – Carl Woll Jan 26 '18 at 6:25
• @Carl, I have replaced the solver with the code FindRoot[DA, {sb, I}] , its working. But I want vary the initial point (in the above case I have given 0+0*I). How can I do that within the Do loop?? With the motive to find different solutions. – T S Singh Jan 26 '18 at 8:00

Not a good idea to start variable names with capitals, hence the renaming below.

ll = 10;
w = 0.02;
h = 0.02;
yy = 209*10^9;
mi = w*(h^3)/12;
ei = yy*mi;
mu = 0.001;
rho = 7800;
aa = w*h;
cn = 4*Pi*mu/(Log[2*ll/w] + 0.5) // N;
gama = rho*aa*ei/(cn^2*ll^4);

Now define the functions of interest so they only operate on explicit numeric input.

bB[pb_?NumericQ, sb_?NumericQ] :=
b /. NSolve[b^4 + pb*b^2 + sb + gama*sb^2 == 0, b];
dCoeffMat[pb_?NumericQ, sb_?NumericQ] :=
With[{b = bB[pb, sb]}, {{1, 1, 1, 1}, {b[]^3 + pb*b[],
b[]^3 + pb*b[], b[]^3 + pb*b[],
b[]^3 + pb*b[]}, {b[]*Exp[b[]], b[]*Exp[b[]],
b[]*Exp[b[]],
b[]*Exp[b[]]}, {b[]^2*Exp[b[]] + pb*Exp[b[]],
b[]^2*Exp[b[]] + pb*Exp[b[]],
b[]^2*Exp[b[]] + pb*Exp[b[]],
b[]^2*Exp[b[]] + pb*Exp[b[]]}}]
dA[pb_?NumericQ, sb_?NumericQ] := Det[dCoeffMat[pb, sb]]

From here one can use FindRoot to get a value for sb that makes the determinat vanish given a value for pb. I use the real part of the determinant to avoid imaginary fuzz.

sbvals =
sb /. Table[FindRoot[Re[DA[pb, sb]] == 0, {sb, 1}], {pb, .1, 4., .1}]

(* Out= {0.0000893249467205, 0.000180282234353, 0.000271244538972, \
0.00036220809792, 0.000453172158726, 0.000544136470395, \
0.000635100925702, 0.000726065470127, 0.000817030074505, \
0.000907994720715, 0.000998959397339, 0.00108992409677, \
0.00118088881389, 0.00127185354457, 0.00136281828633, \
0.00145378303708, 0.0015447477952, 0.00163571255949, \
0.00172667732895, 0.00181764210281, 0.00190860688046, \
0.00199957166138, 0.00209053644524, 0.00218150123126, \
0.00227246601965, 0.00236343081001, 0.002454395602, 0.00254536039557, \
0.00281298528292, 0.00312767097318, 0.00341574060368, \
0.00368348997729, 0.00393511426993, 0.00417357535177, \
0.00440105694776, 0.00461922434613, 0.00482938236888, \
0.00503257643349, 0.00522965991584, 0.00542134060187} *)
• @ Daniel Thanks. The results shows real values of sb. I am interested in both the Re[sb] and Im[sb]. How and where do I have to modify in your code ? – T S Singh Jan 26 '18 at 18:22
• Try removing the Re[...] wrapper in that FindRoot and maybe give initial guess {sb, 1+I} (I don't have time to test this at the moment). – Daniel Lichtblau Jan 26 '18 at 18:31
• @ Daniel I have implemented your suggestion, now the imaginary parts of sb are visible. One last doubt (hope not troubling you much :-)) suppose in implementing FindRoot if I want to vary the initial guess (as complex number) with the purpose to search roots in wide range of complex plane. How can it be done?? – T S Singh Jan 27 '18 at 3:22
• You could use a Table that accounts for all the starting values. Something like Table[FindRoot[Re[DA[pb, sb]] == 0, {sb, sbval}], {pb, .1, 4., .1}], {sbval, sbinits}] where sbinits is the list of initial values you want to try. – Daniel Lichtblau Jan 27 '18 at 15:01
• @ Daniel Thanks – T S Singh Jan 29 '18 at 16:40