# Still having problem with solving high polynomials with exponential equations simultaneously

I need to find the value of Ta and Tc for two simultaneous equations to get a numerical answer using NSolve. I have tried using FindMinimum and NMinimize and still could not solve the problem. Here is my code:

qe = 1.6*10^-19;
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 = 1.0; (*eV*)
V = 1.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);
k4 = 1;
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);
k10 = 1;
k11 = (ATi*\[Sigma]*\[Epsilon]in)/(Ul*AL);
k12 = Tl + (\[Epsilon]*\[Sigma]*Tl^4)/Ul;

eqn1 = 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 -
k4*Ta^4 - k2*Ta^3*Exp[-(wfc*qe)/(kB*Ta)] -
k3*(wfc*qe + qe*V)*Ta^2*Exp[-(wfc*qe)/(kB*Ta)] - k5;
eqn2 = k6*Ta^4 + k7*Ta^3*Exp[-(wfc*qe)/(kB*Ta)] +
k8*(wfc*qe)*Ta^2*Exp[-(wfc*qe)/(kB*Ta)] + Ta - k9*Tc^4 -
k7*Tc^3*Exp[-(wfc*qe + V*qe)/(kB*Tc)] -
k8*(wfc*qe)*Tc^2*Exp[-(wfc*qe + qe*V)/(kB*Tc)] - k10;

sol = NSolve[{eqn1 == 0, eqn2 == 0}, {Ta, Tc}, Reals]


The value of Ta should be close to 300 and Tc is around 1500.

• The first line has a typo, should be qe = 1.6*10^-19;. – Roman Jun 11 at 21:22
• Edited. Thnks for noticing. The problem was still not resolved though – kamegheka Jun 11 at 21:28
• FindRoot is more similar to NSolve than FindMinimum. Have you tried it? – Chris K Jun 11 at 21:29
• Yes I have tried it but it says that "The function value is not a list of real numbers with dimension at (Ta,Tc)={1.,1.}" – kamegheka Jun 11 at 21:31
• Try Plot3D[eqn1, {Ta, 200, 400}, {Tc, 1000, 2000}] and Plot3D[eqn2, {Ta, 200, 400}, {Tc, 1000, 2000}] to see that neither eqn1 nor eqn2 are close to zero within this region. – Jack LaVigne Jun 12 at 1:56

If I try to rescale the problem substituting {Ta -> 1/\[CurlyEpsilon]a, Tc -> 1/\[CurlyEpsilon]c} I get

eqn={eqn1, eqn2} // Rationalize[#, 10^-20] & /. {Ta -> 1/\[CurlyEpsilon]a, Tc -> 1/\[CurlyEpsilon]c}


Assuming Ta>>0,Tc>>0 the solution range would be 0< \[CurlyEpsilon]a, \[CurlyEpsilon]< 1

NMinimize[{1, eqn\[CurlyEpsilon][] == 0,eqn\[CurlyEpsilon][] == 0, \[CurlyEpsilon]a >0, \[CurlyEpsilon]c > 0}, {\[CurlyEpsilon]a, \[CurlyEpsilon]c}]

{1/\[CurlyEpsilon]a, 1/\[CurlyEpsilon]c} /. %[] (* {Ta,Tc}*)
(* {343.695, 0.0501583} *)


Value Ta ~O whereas Tc<<1500