I believe this is a question of more list manipulation than plotting, but I hope you could help me. I solve numerically a system of equations of 3 variables, take only real solutions where one of the parameters is within the unit circle and evaluate its stability for different cost parameter C as

costValues = {};
solutionValues = {};
stabilityValues = {};
(*real steady state solutions where m is within the unit circle*)
ssOut = Select[NSolve[{x == xEq, p == pEq, m == mEq} /.{C -> costVal}, {x, p, m}], 
Element[{x, p, m} /. #, Reals] && 0 < ({x, p, m} /. #)[[3]] < 1 &];
(*calculate eigenvalues for the Jacobian for each solution*)
eigenValues = Table[N[Eigenvalues[jacobian /. {C -> costVal} /.solution]], {solution, ssOut}];
(*check how many eigenvalues are within the unit circle*)
stability = Table[Length[Cases[eigenValue, x_ /; -1 < Re[x] < 1]], {eigenValue,  eigenValues}];
(* append the results to the lists *)
AppendTo[costValues, costVal];
AppendTo[solutionValues, ssOut];
AppendTo[stabValues, stability];
, {costVal, 0, 1}]

My result lists look the following

costVals= {0., 1};
solutionValues= {{}, {{x -> 28.8952, p -> 0.059159, m -> 0.010203}, {x -> -28.8946, 
p -> -0.0577157, m -> 0.010203}, {x -> 28.8886, p -> 0.0432996, 
m -> 0.010203}, {x -> -28.888, p -> -0.0418563, m -> 0.010203}}};
stabValues= {{}, {3, 2, 1, 1}};

Of course, there could be many more costVar values. The point is, that for each costVar value I could have a number of solutions with different stability values.

I would like to plot on separate diagrams different solution values x, p and m against the costValues and with a condition that if a stabValue is equal 2 the color should be green (stable) and if not then it should be red (unstable). So for example above, on horizontal axis I would like to have values 0 and 1 and on vertical axix for instance for variable x I would like to 0 points for costVal=0 and 4 points for costVal=1 where the second solution should be marked green and the other red. I spent a couple of days on that without a success, do you have any idea how to do that efficiently?


To get the x values for example so you can plot them, you can do the following:

listx = {};

AppendTo[listx, x /. solutionValues[[2]][[i]]], {i, 1, 

After that, you use Transpose to get whatever you want to plot against listx. You use If[] to get all your conditional statements on the length of costVals.

I hope it helps, L.


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.