ListPlot3D of the results of FindRoot

Question

I have a system of 2 simultaneous equations where i want to calculate the values of Tc and Tc which ranging with the values of two variables which are wfc and V.

I tried to evaluate it by using FindRoot manage to solve the equation which gave me a table of values of two variables {{Tc-> values, Ta ->values}}. Now I would like to visualize the results by 3D plotting the values of Tc as a function of wfc and V as well as Ta as a function of wfc and V. There's 2 values in the solution of FindRoot and how do I separate them to make two different 3D plot with two different variables used to solve them. Here is my code for solving the values of Tc and Ta

qe = 1.6*10^-19;  (*electron charge*)
Ah = 0.0025;(*m^2*)
AL = 0.025; (*m^2*)
ALoss = 2.5*10^-4; (*m^2*) 
ATi = 2.5*10^-4; (*m^2*)
h = 6.626*10^-34; (*J.s*)
kB = 1.38*10^-23; (*J/K*)
Krd = 1.2*10^6; (*A/m^2.K^4*)
\[Epsilon]in = 0.18; 
\[Epsilon] = 0.9;
Th = 1500; (*K*)
Tl = 300; (*K*)
T0 = 300; (*K*)
Uh = 1.1022*10^4;  (*W/m.K^4*)
Ul = 2.2045*10^4; (*W/m.K^4*)
Uloss = 1.1022*10^4;  (*W/m.K^4*)
\[Sigma] = 5.67*10^-8; (*W/m.K^4*)
\[Epsilon]in = 0.18;

(*wfc=0.0; (*eV*)*)
(*V=0.0; (*V*)*)

k1 = (ATi*\[Sigma]*\[Epsilon]in)/(Uh*Ah) + (\[Epsilon]*\[Sigma])/Uh;
k2 = (2*kB*Krd*ATi)/(qe*Uh*Ah);
k3 = (Krd*ATi)/(qe*Uh*Ah);
k5 = (ATi*\[Sigma]*\[Epsilon]in)/(Uh*Ah);
k6 = Th + (\[Epsilon]*\[Sigma]*Th^4)/Uh;
k7 = (ATi*\[Sigma]*\[Epsilon]in)/(Ul*AL) + (\[Epsilon]*\[Sigma])/Ul;
k8 = (2*kB*Krd*ATi)/(qe*Ul*AL);
k9 = (Krd*ATi)/(qe*Ul*AL);
k11 = (ATi*\[Sigma]*\[Epsilon]in)/(Ul*AL);
k12 = Tl + (\[Epsilon]*\[Sigma]*Tl^4)/Ul;

eqn1[wfc_, V_] := 
  k1*Tc^4 + k2*Tc^3*Exp[-(wfc*qe + qe*V)/(kB*Tc)] + 
   k3*(wfc*qe + qe*V)*Tc^2*Exp[-(wfc*qe + qe*V)/(kB*Tc)] + 
   Tc - (k5*Ta^4) - k2*Ta^3*Exp[-(wfc*qe)/(kB*Ta)] - 
   k3*(wfc*qe + qe*V)*Ta^2*Exp[-(wfc*qe)/(kB*Ta)] - k6;
eqn2[wfc_, V_] := 
  k7*Ta^4 + k8*Ta^3*Exp[-(wfc*qe)/(kB*Ta)] + 
   k9*(wfc*qe)*Ta^2*Exp[-(wfc*qe)/(kB*Ta)] + Ta - k11*Tc^4 - 
   k8*Tc^3*Exp[-(wfc*qe + V*qe)/(kB*Tc)] - 
   k9*(wfc*qe)*Tc^2*Exp[-(wfc*qe + qe*V)/(kB*Tc)] - k12;

Table[FindRoot[{eqn1[wfc, V] == 0, 
   eqn2[wfc, V] == 0}, {{Tc, 1500}, {Ta, 300}}], {wfc, 0.1, 1.5, 
  0.05}, {V, 0.1, 1.5, 0.05}]

sol = Quiet[
   Flatten[Table[{wfc, V, Tc, Ta} /. 
      FindRoot[{eqn1[wfc, V] == 0, 
        eqn2[wfc, V] == 0}, {{Tc, 1500}, {Ta, 300}}], {wfc, 0.2, 1.5, 
      0.05}, {V, 0.2, 1.5, 0.05}], 1]];

ListPlot3D[sol[[All, {1, 2, 3}]], 
 AxesLabel -> {"\!\(\*SubscriptBox[\(\[Phi]\), \(C\)]\)(eV)", "V(V)", 
   "\!\(\*SubscriptBox[\(T\), \(emitter\)]\)(K) "}, 
 LabelStyle -> Bold, PlotRange -> All, ColorFunction -> "DarkRainbow"]

ListPlot3D[sol[[All, {1, 2, 4}]], 
 AxesLabel -> {"\!\(\*SubscriptBox[\(\[Phi]\), \(C\)]\)(eV)", "V(V)", 
   "\!\(\*SubscriptBox[\(T\), \(collector\)]\)(K)"}, 
 LabelStyle -> Bold, PlotRange -> All, ColorFunction -> "DarkRainbow"]

Many thanks!

Answer 1

One way to do this is to build a table of values rather than the solution rules:

(* build a table with each entry {wfc, V, Tc, Ta} *)
sol = Quiet[
   Flatten[Table[{wfc, V, Tc, Ta} /. 
      FindRoot[{eqn1[wfc, V] == 0, 
        eqn2[wfc, V] == 0}, {{Tc, 1500}, {Ta, 300}}], {wfc, 0.1, 1.5, 
      0.05}, {V, 0.1, 1.5, 0.05}], 1]];

ListPlot3D[sol[[All, {1, 2, 3}]], AxesLabel -> {"wfc", "V", "Tc "}, 
 LabelStyle -> Bold, PlotRange -> All]

ListPlot3D[sol[[All, {1, 2, 4}]], AxesLabel -> {"wfc", "V", "Ta "}, 
 LabelStyle -> Bold, PlotRange -> All]

Note that FindRoot is complaining about convergence. I suppressed the messages with Quiet. You might want to look into that.

Kind regards, David

EDIT:

Additional functional values can be added to the solution set by adding the expressions for them into the definition of the table being built. For example, this adds expressions for J1 and J2 from the comment below. The third expression for eff is difficult to make out from the text in the comment, but it can be done in the same manner. (If it is evaluated as an expression after j1 and j2 are defined then j1 and j2 will be replaced and the resulting form placed in the table definition.)

j1 = Krd Tc^2*Exp[-(wfc qe + qe V)/(kB Tc)];

j2 = Krd Ta^2*Exp[-(wfc qe)/(kB Ta)];

(* build a table with each entry {wfc, V, Tc, Ta, J1, J2} *)
sol2 = Quiet[
   Flatten[
    Table[{wfc, V, Tc, Ta, j1, j2} /.
      FindRoot[{eqn1[wfc, V] == 0, 
        eqn2[wfc, V] == 0}, {{Tc, 1500}, {Ta, 300}}],
     {wfc, 0.1, 1.5, 0.05}, {V, 0.1, 1.5, 0.05}],
    1]
   ];