# Symbolic Simulations

I have 9 nonlinear equations and 10 unknowns. It is not possible to obtain a numeric solution but I do get a parametric solution. It is useful in my case because I want to observe how parameters (p1, p2, Y1, Y2, alpha, beta, gamma) effects the unknowns. I will now present you the code and explain the problems below it.

PG = 100;

p1 = 0.8;

p2 = 0.2;

(Y1 = 10 & Y2 = 15 & alpha = 0.5 & gamma = 0.4 & beta = 0.3)

Solve[C1 - 10 + PG/2 + GC1 + GL1 == 0 &&
C2 - 15 + PG/2 + GC2 + GL2 == 0 &&
GL1 - fd (a1) t 10 == 0 &&
GL2 - fd (a2) t 15 == 0 &&
GC1 - (t (1 - fd) 25 - PG) ((p1)/((p1) + (p2))) == 0 &&
GC2 - (t (1 - fd) 25 - PG) ((p2)/((p1) + (p2))) ==  0 &&
(0.3/GL1) - (0.5/C1) - (0.3/GC1) - (0.4/PG) == 0 &&
(0.3/GL2) - (0.5/C2) - (0.3/GC2) - (0.4/PG) == 0 &&
GC1 - 4 GC2 == 0,
{C1, C2, GC1, GC2, GL1, GL2, fd, a1, a2, t}]


I want to perform a simulation on the unknown variables by changing the values of parameters. In the above code parameters take only 1 value (i.e. p1 = 0.8) but keeping other values constant I want to change the values of p1 between 0 and 1 to observe the changes in the unknowns and I want to plot it on the graph. We can call it a sensitivity analysis. If I can learn to change it for p1, I will try myself to change it for the other parameters.

I tried to use "RandomReals" (i.e. p1 = RandomReals[1,5] - which creates 5 random variables between 0 and 1) but I think there was a problem of mapping each numbers to p1.

So, I need help to do the simulation and plotting the change.

I hope everything is clear in my question. Since I am new here and it is my 2nd question please tell me if I need more clarification.

• How do you propose to obtain numeric solutions? You have 10 variables and but 9 equations. – Daniel Lichtblau Apr 30 '13 at 21:59
• I would guess that the reason your question has been downvoted is because you don't identify clearly what it is you want or that you were not able to do. Personally, while I am not going to downvote, I am going to vote to close since I can't really understand what sort of an answer would satisfy you. If you can clarify this, most likely the question will not in fact be closed and an answer may even be forthcoming, – Oleksandr R. May 1 '13 at 3:17
• @DanielLichtblau I don't obtain a numeric solution. I obtain a parametric solution which works in my case. – HarveyMudd May 1 '13 at 8:28
• @OleksandrR. I will clarify the question. I believe an answer for this will be helpful for the other users too. – HarveyMudd May 1 '13 at 8:31
• @HarveyMudd A closely related approach to what you are looking for (as far as it seems) might be found here mathematica.stackexchange.com/questions/17799/…. – Artes May 1 '13 at 11:55

You might use rule replacement (instead of variable assignment) to inject specific values for the parameters.

ee = C1 - 10 + PG/2 + GC1 + GL1 == 0 &&
C2 - 15 + PG/2 + GC2 + GL2 == 0 && GL1 - fd (a1) t 10 == 0 &&
GL2 - fd (a2) t 15 == 0 &&
GC1 - (t (1 - fd) 25 - PG) ((p1)/((p1) + (p2))) == 0 &&
GC2 - (t (1 - fd) 25 - PG) ((p2)/((p1) + (p2))) ==
0 && (0.3/GL1) - (0.5/C1) - (0.3/GC1) - (0.4/PG) ==
0 && (0.3/GL2) - (0.5/C2) - (0.3/GC2) - (0.4/PG) == 0 &&
GC1 - 4 GC2 == 0;


Here I "instantiate" with the values you indicated. I also rationalize to make it easy to clear denominators. This is not absolutely necessary.

polys = Numerator[
Together[Rationalize[{40 + C1 + GC1 + GL1, 35 + C2 + GC2 + GL2,
GL1 - 10 a1 fd t, GL2 - 15 a2 fd t,
GC1 - 0.8 (-100 + 25 (1 - fd) t),
GC2 - 0.2 (-100 + 25 (1 - fd) t), -0.004 - 0.5/C1 - 0.3/
GC1 + 0.3/GL1, -0.004 - 0.5/C2 - 0.3/GC2 + 0.3/GL2,
GC1 - 4 GC2}/.{PG->100,p1->.8,p2->.2}]]]


Your system actually is underdetermined by two degrees of freedom (despite having only one less equation than variables). One way to see this is to solve it and find the solution requiring the most parameters (unsolved for variables).

solnSym = Solve[polys == 0]


The result is bigger than I want to print but does show a two dimensional component.

I'll also show how one can obtain a numerical solution set efficiently by substituting values for a couple of the variables.

soln = NSolve[polys /. {C1 -> 11/10, t -> 23/4}]

(* {{a1 -> -2.27002362868, a2 -> -1.31333144612, GL1 -> -40.4561691374,
GC1 -> -0.643830862642, fd -> 0.309946355327, C2 -> 0.270080664682,
GL2 -> -35.109122949, GC2 -> -0.16095771566}, {a1 -> -2.27002362868,
a2 -> -0.00598749869422, GL1 -> -40.4561691374,
GC1 -> -0.643830862642, fd -> 0.309946355327, C2 -> -34.6789792413,
GL2 -> -0.16006304306,
GC2 -> -0.16095771566}, {a1 -> 0.0173166796902,
a2 -> -0.887078857324, GL1 -> 0.664654466462, GC1 -> -41.7646544665,
fd -> 0.667518734491, C2 -> 26.5133900918, GL2 -> -51.0722264752,
GC2 -> -10.4411636166}, {a1 -> 0.0173166796902,
a2 -> -0.10131060352, GL1 -> 0.664654466462, GC1 -> -41.7646544665,
fd -> 0.667518734491, C2 -> -18.7260312786, GL2 -> -5.83280510478,
GC2 -> -10.4411636166}} *)


This should give some ideas on how to go about altering parameter values prior to solving.

• this was really helpful. Thank you very much. – HarveyMudd May 2 '13 at 11:17