StackOverflow2013

Note that there are some explanatory texts on larger screens.

plurals

PO
text
Body
copied!<p><strong>Updated</strong>: see below.</p> <p>I'd approach this first by visualizing the imaginary parts of the roots:</p> <p><img src="https://i.stack.imgur.com/6wY1v.png" alt="plot of the imaginary parts of the roots"></p> <p>This tells you three things immediately: 1) the first root is always real, 2) the second two are the conjugate pairs, and 3) there is a small region near zero in which all three are real. Additionally, note that the exclusions only got rid of the singular point at <code>x=0</code>, and we can see why when we zoom in:</p> <p><img src="https://i.stack.imgur.com/aU06u.png" alt="zoom in of above photo with x between 0 and 1.5"></p> <p>We can then use the <code>EvalutionMonitor</code> to generate the list of roots directly:</p> <pre><code>Map[Module[{f, fcn = #1}, f[x_] := Im[fcn]; Reap[Plot[f[x], {x, 0, 1.5}, Exclusions -> {True, f[x] == 1, f[x] == -1}, EvaluationMonitor :> Sow[{x, f[x]}][[2, 1]] // SortBy[#, First] &];] ]&, geyvals] </code></pre> <p>(Note, the <code>Part</code> specification is a little odd, <code>Reap</code> returns a <code>List</code> of what is sown as the second item in a <code>List</code>, so this results in a nested list. Also, <code>Plot</code> doesn't sample the points in a straightforward manner, so <code>SortBy</code> is needed.) There may be a more elegant route to determine where the last two roots become complex, but since their imaginary parts are piecewise continuous, it just seemed easier to brute force it. </p> <p><strong>Edit</strong>: Since you've mentioned that you want an automatic method for generating where some of the roots become complex, I've been exploring what happens when you substitute in <code>y -> p + I q</code>. Now this assumes that <code>x</code> is real, but you've already done that in your solution. Specifically, I do the following</p> <pre><code>In[1] := poly = g.RotationMatrix[Pi/2].h.g /. {y -> p + I q} // ComplexExpand; In[2] := {pr,pi} = poly /. Complex[a_, b_] :> a + z b & // CoefficientList[#, z] & // Simplify[#, {x, p, q} \[Element] Reals]&; </code></pre> <p>where the second step allows me to isolate the real and imaginary parts of the equation and simplify them independent of each other. Doing this same thing with the generic 2D polynomial, <code>f + d x + a x^2 + e y + 2 c x y + b y^2</code>, but making both <code>x</code> and <code>y</code> complex; I noted that <code>Im[poly] = Im[x] D[poly, Im[x]] + Im[y] D[poly,[y]]</code>, and this may hold for your equation, also. By making <code>x</code> real, the imaginary part of <code>poly</code> becomes <code>q</code> times some function of <code>x</code>, <code>p</code>, and <code>q</code>. So, setting <code>q=0</code> always gives <code>Im[poly] == 0</code>. But, that does not tell us anything new. However, if we </p> <pre><code>In[3] := qvals = Cases[List@ToRules@RReduce[ pi == 0 && q != 0, {x,p,q}], {q -> a_}:> a]; </code></pre> <p>we get several formulas for <code>q</code> involving <code>x</code> and <code>p</code>. For some values of <code>x</code> and <code>p</code>, those formulas may be imaginary, and we can use <code>Reduce</code> to determine where <code>Re[qvals] == 0</code>. In other words, we want the "imaginary" part of <code>y</code> to be real and this can be accomplished by allowing <code>q</code> to be zero or purely imaginary. Plotting the region where <code>Re[q]==0</code> and overlaying the gradient extremal lines via</p> <pre><code>With[{rngs = Sequence[{x,-2,2},{y,-10,10}]}, Show@{ RegionPlot[Evaluate[Thread[Re[qvals]==0]/.p-> y], rngs], ContourPlot[g.RotationMatrix[Pi/2].h.g==0,rngs ContourStyle -> {Darker@Red,Dashed}]}] </code></pre> <p>gives</p> <p><img src="https://i.stack.imgur.com/N5eo7.png" alt="x-y plot showing gradient extremals and region where there are 3 real roots"></p> <p>which confirms the regions in the first two plots showing the 3 real roots.</p>

Querying!

Guidance

An individual column

Larger individual text columns get their own page to allow for proper reading.

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload