Assumptions:
Under these assumptions we can estimate the temperature of the board by simply equating heat flows. The heat coming in per unit time is from Joule heating from the current running through the copper and is given by
$$
q_{in}=I^2R.
$$
The heat flowing out has two escape mechanisms; radiative heat transfer to the surroundings which is given by
$$
q_{rad}=\epsilon\sigma A(T^4-T_s^4)
$$
and convective heat transfer to the air which is given by
$$
q_{conv}=hA(T-T_s).
$$
Now we just equate the heat flows
$$
\begin{align}
q_{in}&=q_{out}\\
I^2R&=A\left[\sigma(T^4-T_s^4)(\epsilon_{cu}+\epsilon_{pcb})+(T-T_s)(h_{up}+h_{down})\right]
\end{align}
$$
we can rearrange this to look like a quartic equation
$$
\sigma(\epsilon_{cu}+\epsilon_{pcb})T^4+(h_u+h_d)T -
\left[\frac{I^2R}{A}+\sigma(\epsilon_{cu}+\epsilon_{pcb})T_s^4+(h_u+h_d)T_s\right]=0.
$$
This isn't easy to solve analytically, but Mathematica would have no problem. I put it into Python and numerically found the minimum of the absolute value. I assumed a total resistance of $R=1\Omega$ so your actual results may vary. The results are shown below.
Here is the Python code also:
import scipy.optimize as opt
import numpy as np
import matplotlib.pyplot as plt
# Define the function with some extra variables
def tempOpt( t, ts, i):
sigma = 5.67e-8
ec = 0.78
ep = 0.50
hu = 7.25
hd = 3.63
r = 1
a = 0.0103
e = ec + ep
h = hu + hd
out = sigma*e*t**4 + h*t - i**2*r/a - sigma*e*ts**4 - h*ts
return out
# Decide the ranges for the current and temperature
currents = np.linspace( 0.1, 20, 50)
temps = [0, 20, 40, 60, 80, 100]
# Calculate the value
resDict = dict()
for temp in temps:
resVec = np.zeros( np.shape( currents))
cnt = 0
for current in currents:
# Define a new function of one variable and a minimum at zero
# Don't forget to convert between Kelvin and Celcius
def tempNow( t):
return abs( tempOpt( t+273.15, temp+273.15, current))
# Find the minimum with a reasonable guess
guess = temp + 30
resVec[cnt] = opt.fmin( tempNow, 30)
cnt += 1
# Store results in dictionary
resDict[temp] = resVec
# Plot
plt.figure(1)
plt.clf()
for temp in temps:
plt.plot( currents, resDict[temp], lw=2)
plt.xlabel('Current (A)')
plt.ylabel('Board Temp ($^\circ$C)')
leg = plt.legend( temps, loc=2)
leg.set_title( '$T_s$ ($^\circ$C)')
plt.title('Board Temperature (R=1$\Omega$)')
