# Plot “Findroot” with very complex variables

I ran into an extremely complex "zero problem" and don't know how to solve and plot it. Can anyone help with this, please? [Pi]s and [Pi]r are objective functions, pd, po, qd, qo, and [Phi] (0<[Phi]<1) are variables. pd, po, qd, qo are determined to maximize [Pi]s and [Phi] is determined to maximize [Pi]r. This problem needs to be addressed by backward induction. First, I need to solve pd and po. Then I intend to solve qd and qo (with some constraints), but results are so much complicated and thus I use "Reduce" afterward to select solutions. Last, [Phi] can be obtained by using the above solutions pd, po, qd, and qo. I want to Plot [Phi] when the parameter c is within a certain range, say 0<c<0.1. I tried several times but get no feasible solutions, which is weird. Then how should I do? I am a beginner so really hope to get some answers. Thanks a lot!

dD := (pd - po)/(qd - qo) - pd/qd

do := 1 - (pd - po)/(qd - qo)

\[Pi]s :=
dD (pd - 1/2 (qd*qd) - c) + do (1 - \[Phi]) (po - 1/2 (qo*qo))

\[Pi]r := do (\[Phi]) (po - 1/2 (qo*qo))

s.t. qo > qd > 0 && 0 < c < pd - 1/2 (qd*qd) &&
1 - (pd - po)/(qd - qo) > 0 && pd/qd - po/qo < 0 &&
0 < \[Phi] < 1, \[Phi]\[Element]Reals

Step 1     FullSimplify[Solve[{D[\[Pi]s, po] == 0, D[\[Pi]s, pd] == 0}, {po, pd}]]
which gives po -> (qo (2 qo (2 + qo) (-1 + \[Phi]) + qd (-4 + qo (-2 + \[Phi])) (-1 + \[Phi]) + 2 c \[Phi] + qd^2 \[Phi]))/(2 qd (-2 + \[Phi])^2 + 8 qo (-1 + \[Phi])), pd -> (-2 c qd (-2 + \[Phi]) - qd^3 (-2 + \[Phi]) + 4 c qo (-1 + \[Phi]) + 2 qd^2 (-1 + \[Phi]) (-2 + qo + \[Phi]) + qd qo (-1 + \[Phi]) (4 + (-2 + qo) \[Phi]))/(2 qd (-2 + \[Phi])^2 + 8 qo (-1 + \[Phi]))

Step 2   Solve[{D[\[Pi]s, qo] == 0, D[\[Pi]s, qd] == 0}, {qo, qd}]   (*very complex*)

Step 3   Reduce[qo > qd > 0 && 0 < c < pd - 1/2 (qd*qd) && 1 - (pd - po)/(qd - qo) > 0 && pd/qd - po/qo < 0 && 0 < \[Phi] < 1, \[Phi], Reals]         (*also complex*)

Step 4   g[c_?NumericQ] = \[Phi] /. FindRoot[D[\[Pi]r, \[Phi]] == 0, {\[Phi], 0, 0, 1}];
Assuming[\[Phi] \[Element] Reals, Plot[{g[c]}, {c, 0, 0.1}]]

• The line beginning s.t. qo>qd>0 ... I guess s t are not further variables. Is this a statement of conditions? Your step 2 does not solve for me in a reasonable time. Did you actually get a solution? – Hugh Jun 6 at 16:44
• The letters "s.t." proceed with constraints that the functions must follow, namely qo>qd>0 ... . And I can get 16 pairs of solutions in Step 2, however, only 8 of them may be feasible. And then, it is extremely hard to proceed. So, I am wondering if there is any other way to solve this kind of problem (Plot [phi] and see how it changes with c). – Qiaoke Zhang Jun 6 at 16:56
• Crossposted here. – Rohit Namjoshi Jun 6 at 19:14

I have had to make several assumptions about what you are doing but this seems to be a way forward. First I am going to ignore your constraints to begin with and just try and find some solutions.

I start by defining your equations and then solving the equations in your Step 1. I then look at the equations you wish to solve in your Step 2 and substitute in the solutions from Step1. This is what I get

ClearAll[dD, do, πs, πr]
dD = (pd - po)/(qd - qo) - pd/qd;
do = 1 - (pd - po)/(qd - qo);
πs = dD (pd - 1/2 (qd*qd) - c) +
do (1 - ϕ) (po - 1/2 (qo*qo));
πr = do (ϕ) (po - 1/2 (qo*qo));
sol1 = FullSimplify[
Solve[{D[πs, po] == 0, D[πs, pd] == 0}, {po, pd}]];
eqn1 = {D[πs, qo] == 0, D[πs, qd] == 0} /. First@sol1 //
FullSimplify


The equation I get out is not too bad

{((2 c (-2 + ϕ) - 2 qd (-2 + ϕ) + qd^2 (-2 + ϕ) -
2 (-2 + qo) qo (-1 + ϕ)) (-1 + ϕ) (2 c (-2 + ϕ) +
qd^2 (-2 + ϕ) + 2 qo (-2 + 3 qo) (-1 + ϕ) +
2 qd (-2 + ϕ) (1 - 2 qo + (-1 + qo) ϕ)))/(
qd (-2 + ϕ)^2 + 4 qo (-1 + ϕ)) == 0, (
qo (-1 + ϕ) (4 c +
qd (2 qd - 2 qo + (-2 + qo) ϕ)) (-2 c (qd (-2 + ϕ)^2 +
2 qo (-1 + ϕ)) +
qd (qd^2 (-2 + ϕ)^2 + 6 qd qo (-1 + ϕ) +
qo (-1 + ϕ) (-2 qo + (-2 + qo) ϕ))))/(
qd (qd (-2 + ϕ)^2 + 4 qo (-1 + ϕ))) == 0}


I now try and solve this equation and determine the number of solutions

sol2 = Solve[eqn1, {qo, qd}];
Dimensions@sol2


(* {16, 2} *)

So lets look at some of these solutions. Solutions 11 to 14 seem very complicated but the rest are reasonable. Here are the reasonable ones.

sol2[[{1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 15, 16}]] // TableForm


Now let's look at the equation you wish to find roots for. Start by substituting for the first set of solutions.

eqn2 = D[πr, ϕ] /. First@sol1 // Simplify

(* (qo (2 c (-2 + ϕ) - 2 qd (-2 + ϕ) + qd^2 (-2 + ϕ) -
2 (-2 + qo) qo (-1 + ϕ)) (qd (4 + qo (-2 + ϕ) -
4 ϕ) + 4 qo (-1 + ϕ) - 2 qo^2 (-1 + ϕ) +
2 c ϕ + qd^2 ϕ))/(4 (qd (-2 + ϕ)^2 +
4 qo (-1 + ϕ))^2)  *)


Now we want to substitute in the solutions for the second set of equations. Start by looking at the simple solutions in the 16

 eqns3 = Table[{n,
eqn2 /. sol2[[n]]}, {n, {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 15,
16}}] // Simplify;
eqns3 // TableForm


So these are trivial. Always zero or never zero.

So now you have to look at the remaining equations. Start by looking at solution 11. Complex solutions are possible so first look at the real and imaginary parts.

eqn11 = eqn2 /. sol2[[11]];
Plot3D[Re[eqn11], {c, 0, 0.1}, {ϕ, -100, 10},
PlotRange -> {All, All, {-0.1, 0.1}}]
Plot3D[Im[eqn11], {c, 0, 0.1}, {ϕ, -100, 10},
PlotRange -> {All, All, {-0.1, 0.1}}]


The real part is always greater than zero and the imaginary part always equals zero. So no solutions equalling zero here.

I suggest you go and look at the remaining 10 to 14 solutions and see if you get some possibilities there.

Hope that helps.

• Thank you so much! But can I ask what "First@" means in the first part of your code? – Qiaoke Zhang Jun 6 at 22:42
• Please look up in help. First gives the first element of a list. Solve gives a list of lists and in this case there is only one element which I want. – Hugh Jun 7 at 5:07