1
$\begingroup$

First off, I am new to using Mathematica. What I am trying to do is the following:

  1. I have a matrix that is a function of two variables. I want to plot the region in the space of these 2 variables where the matrix is positive definite.
  2. How I would do this in any other programming language is checking for positive definiteness of the matrix while iterating over the 2 variables in a nested loop and then try to plot the regions where the condition comes out to be true.
  3. I am not sure how I would go about doing this in Mathematica. I am able to write a similar code as to what I described above using Do and I get an output. But I am not able to figure out how to plot that data.

Also, is there a better way to do it in Mathematica than the method I described above?

This is the code I have so far.

K = 7;
Array[x, 2 K + 1, 0];
x[0] = 1; 
x[1] = 0;
Do[
 x[2] = xa;
 x[n] = 0;
 x[n + 3] = (n/(n + 2)) (Ener) x[
     n - 1] + (n (n - 1) (n - 2)/(4 n + 8)) x[
     n - 3] - ((n + 1)/(n + 2)) x[n + 1];
 (M = Table[x[i + j - 2], {i, K + 1}, {j, K + 1}]) // MatrixForm;
 pred = PositiveSemidefiniteMatrixQ[M];
 Print[pred, xa, Ener], {xa, 0.295, 0.310, 0.001}, {Ener, 1.34, 1.44, 
  0.01}, {n, 1, 2 K - 2, 2}
 ] 

Any help is highly appreciated. Thanks!

$\endgroup$
6
  • $\begingroup$ Please show the code you have so far from point 3 by editing your question and adding the code as text, appropriately formatted. $\endgroup$
    – MarcoB
    May 16, 2020 at 21:48
  • $\begingroup$ I have added it now $\endgroup$
    – adithya
    May 16, 2020 at 21:58
  • 1
    $\begingroup$ there was not a single point where pred was positive semidefinite when running your code. $\endgroup$
    – Nasser
    May 16, 2020 at 23:30
  • $\begingroup$ Do you know you could use HankelMatrix[] to construct your M? $\endgroup$ May 16, 2020 at 23:54
  • $\begingroup$ In your code, every single matrix M is the same -- which is why Nasser sees pred as always False. Also your M is only a function of one variable which seems to be called x[13]. $\endgroup$
    – bill s
    May 17, 2020 at 0:32

1 Answer 1

2
$\begingroup$

You have an off-by-one error in your final value of n. I think you actually should iterate $n$ until it reaches $2K-1$, not $2K-2$. You also do not need to define x[1] = 0 at the start; your loop code defines that at every step automatically anyway. Then Array[x, 2 K + 1, 0]; does nothing, since you do not store its output. Additionally, MatrixForm does nothing as well, since its result is never displayed. Finally, much better to collect the results you want in a list for further work using Sow and Reap, rather than Print them to the console.

You also should not use K or any other single-letter capital variable, since most of those are reserved for internal use by Mathematica. In this case, K is used for summation indices in symbolic sums (see ?K).

Here is an amended version of your code:

k = 7;
x[0] = 1;

results = 
  Reap[
    Do[
      x[2] = xa;
      x[n] = 0;
      x[n + 3] = (n/(n + 2)) (Ener) x[n - 1] + (n (n - 1) (n - 2)/(4 n + 8)) x[n - 3] - ((n + 1)/(n + 2)) x[n + 1];
      M = Table[x[i + j - 2], {i, k + 1}, {j, k + 1}];
      If[PositiveSemidefiniteMatrixQ[M], Sow[{xa, Ener, M}], Nothing],
      {xa, 0.295, 0.310, 0.001},
      {Ener, 1.34, 1.44, 0.01},
      {n, 1, 2 k - 1, 2}
    ]
  ][[2, 1]];

You can check that results indeed contains values of xa, Ener for which M is positive semidefinite, as well as the value of M itself:

Dimensions@results (* Out {249, 3} *)

PositiveSemidefiniteMatrixQ[#[[3]]]& /@ results (*Out: {True, True, .... True} *)
$\endgroup$

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.