2
$\begingroup$

I am trying to find multiple complex and real roots of the determinant of symbolic matrix. Someone advised me to use Muller's method. I am wondering if there is a package already written to do this method?

The problem is: The data is given and the function to be solved for is NsymDisp[w_, k_], we solve NsymDisp for k if some values of w is given.

Cnorm = 129.5;
hnorm = 5*10^-6;
data = {a -> 50*10^-6/hnorm,
 h -> 5*10^-6/hnorm,
  \[Mu] -> (84.7)/Cnorm,
  \[Lambda] -> (129.5)/Cnorm,
  Subscript[C, 11] -> (268)/Cnorm,
  Subscript[C, 13] -> (95)/Cnorm,
   Subscript[C, 33] -> (186)/Cnorm,
   Subscript[C, 55] -> (49)/Cnorm,
   Subscript[C, 44] -> (37)/Cnorm,
   Subscript[C, 66] -> (95)/Cnorm,
   Subscript[\[Rho], c] -> 1};
  NsymDisp[w_, k_] := 
  Module[{\[Omega] = w, Cp, Cs, p, 
   q, \[Alpha]11, \[Beta]11, \[Gamma]11, \[Alpha]12, \[Beta]12, \
 \[Gamma]12, \[Alpha]21, \[Beta]21, \[Gamma]21, \[Alpha]22, \[Beta]22, \
 \[Gamma]22},
  Cp = Sqrt[(\[Lambda] + 2 \[Mu])/Subscript[\[Rho], c]] /. data;
   Cs = Sqrt[\[Mu]/Subscript[\[Rho], c]] /. data;
   kp = \[Omega]/Cp;
   ks = \[Omega]/Cs;
   p = Sqrt[kp^2 - k^2];
   q = Sqrt[ks^2 - k^2];   
   \[Alpha]11 = 2 \[Mu] k p;
   \[Beta]11 = (Subscript[\[Rho], 
   c] \[Omega]^2 - (Subscript[C, 11] + 
     Subscript[C, 
      13] (\[Lambda] - Subscript[C, 13])/Subscript[C, 33]) k^2 - 
      Subscript[C, 13] (\[Lambda] + 2 \[Mu]) p^2/Subscript[C, 33]) k;
    \[Gamma]11 = (Subscript[C, 13]^2 (Subscript[C, 55] - 2 \[Mu]) k^2 - 
      Subscript[C, 
       33] (Subscript[C, 55] - 2 \[Mu]) (Subscript[C, 11] k^2 - 
         Subscript[\[Rho], c] \[Omega]^2) - 
      Subscript[C, 13] Subscript[C, 
       55] (2 \[Mu] k^2 - 
         Subscript[\[Rho], c] \[Omega]^2)) k p/(2 Subscript[C, 33]
        Subscript[C, 55]);
    NAs11 = ( \[Alpha]11 + h^2 \[Gamma]11)*Sin[p a] + 
     h \[Beta]11 Cos[p a];
    \[Alpha]12 = \[Mu] (q^2 - k^2);
    \[Beta]12 = (Subscript[\[Rho], 
        c] \[Omega]^2 - (Subscript[C, 11] - 
          Subscript[C, 
           13] (Subscript[C, 13] + 2 \[Mu])/Subscript[C, 33]) k^2) q;
     \[Gamma]12 = (Subscript[C, 
        13]^2 ((\[Mu] - Subscript[C, 55]) k^2 - \[Mu] q^2) k^2 + 
       Subscript[C, 
        33] ((Subscript[C, 55] - \[Mu]) k^2 + \[Mu] q^2) (Subscript[C, 
           11] k^2 - Subscript[\[Rho], c] \[Omega]^2) - 
        Subscript[C, 13] Subscript[C, 
        55] (Subscript[\[Rho], 
           c] \[Omega]^2 + \[Mu] (q^2 - k^2)) k^2)/(2 Subscript[C, 33]
         Subscript[C, 55]);
      NAs12 = (\[Alpha]12 + h^2 \[Gamma]12)*Sin[q a] + 
      h \[Beta]12 Cos[q a];          
     \[Alpha]21 = -\[Lambda] k^2 - (\[Lambda] + 2 \[Mu]) p^2;
     \[Beta]21 = (Subscript[\[Rho], c] \[Omega]^2 - 2 \[Mu] k^2) p;
       \[Gamma]21 = ((\[Lambda] + 
         2 \[Mu]) (Subscript[\[Rho], c] \[Omega]^2 + 
            Subscript[C, 13]
             k^2) p^2 - ((Subscript[C, 13] + Subscript[C, 
             33] - \[Lambda]) Subscript[\[Rho], 
            c] \[Omega]^2 - (Subscript[C, 11] Subscript[C, 33] + 
              Subscript[C, 
              13] (\[Lambda] - Subscript[C, 
                13])) k^2) k^2)/(2 Subscript[C, 33]);
      NAs21 = (\[Alpha]21 + h^2 \[Gamma]21) Cos[p a] + 
      h \[Beta]21 Sin[p a];
    \[Alpha]22 = 2 \[Mu] k q;
    \[Beta]22 = (\[Mu] (k^2 - q^2) - 
         Subscript[\[Rho], c] \[Omega]^2) k;
      \[Gamma]22 = ((Subscript[C, 11] Subscript[C, 33] - 
           Subscript[C, 
           13] (Subscript[C, 13] + 2 \[Mu])) k^2 - (Subscript[C, 13] + 
            Subscript[C, 33] + 2 \[Mu]) Subscript[\[Rho], 
        c] \[Omega]^2) k q/(2 Subscript[C, 33]);
   NAs22 = (\[Alpha]22 + h^2 \[Gamma]22 ) Cos[q a] + 
       h \[Beta]22 Sin[q a];             
     symMat = {{NAs11, NAs12}, {NAs21, NAs22}} /. data;
    Det[symMat // Chop]
      ]
$\endgroup$
  • 2
    $\begingroup$ may be this muller-method-in-mathematica $\endgroup$ – Nasser Nov 12 '17 at 4:30
  • 1
    $\begingroup$ If the "symbolic matrix" is not too complicated, can you include it here? $\endgroup$ – J. M.'s discontentment Nov 12 '17 at 4:34
  • $\begingroup$ If you are satisfied to find the roots by any means, there are several routines described on this site. For instance, if the determinant is real, then question [16439] (mathematica.stackexchange.com/q/16439/1063) and its answers provide multiple approaches. $\endgroup$ – bbgodfrey Nov 12 '17 at 14:26
  • $\begingroup$ Absent a specific example, there is not really much of a question here. No clear indication of why Roots or NRoots would not suffice, for example. $\endgroup$ – Daniel Lichtblau Nov 12 '17 at 15:30
  • 1
    $\begingroup$ Apparently i had to be validated as "Not a robot". Now I can add comments. my hat's off to MSE for such an informative indication (read: none whatsoever) of what was the problem. $\endgroup$ – Daniel Lichtblau Nov 12 '17 at 17:23
2
$\begingroup$

[Not a full answer but I m having browser problems and cannot add a comment.]

If a range is specified then NSolve can find roots in that range.

In[22]:= Chop[NSolve[NsymDisp[w, 5] == 0 && Abs[w] < 6, w]]

(* Out[22]= {{w -> -5.83531}, {w -> -5.65654}, {w -> -5.31554}, {w -> \
-5.1549}, {w -> -5.00189}, {w -> -4.85735}, {w -> -4.72215}, {w -> \
-4.59719}, {w -> -4.48334}, {w -> -4.38144}, {w -> -4.29225}, {w -> \
-4.21635}, {w -> -4.15414}, {w -> -4.10577}, {w -> -4.07127}, {w -> \
-4.05058}, {w -> -4.04368}, {w -> 4.05058}, {w -> 4.07127}, {w -> 
4.10577}, {w -> 4.15414}, {w -> 4.21635}, {w -> 4.29225}, {w -> 
4.38144}, {w -> 4.59719}, {w -> 4.72215}, {w -> 5.00189}, {w -> 
5.1549}, {w -> 5.31554}, {w -> 5.65654}, {w -> 5.83531}} *)`
| improve this answer | |
$\endgroup$
  • $\begingroup$ Thanks Daniel, actually the function is NsymDisp[5,w]. Reduce and NSolve don't give the same results. $\endgroup$ – qahtah Nov 13 '17 at 18:07
  • $\begingroup$ Are you sure they are different? NSolve[NsymDisp[5, w] == 0 && Abs[w] < 3, w] and Reduce[NsymDisp[5, w] == 0 && Abs[w] < 3, w] look likew the same results to my eye. This is using Mathematica 11.2, in case that matters. $\endgroup$ – Daniel Lichtblau Nov 13 '17 at 19:50
  • $\begingroup$ Thanks Daniel, I am not talking about a single point, I varied the value of w from 0 to 10 with an increment of 0.1. $\endgroup$ – qahtah Nov 15 '17 at 3:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.