function y=f(TB)

csf=0.013; hfg=2256e3; cp=4.215e3; Pr=1.76; mu=282.7e-6; sigma=58.917e-3;
g=9.807; rhol=958; rhog=0.6037; sigmarad=5.6697e-8; k=0.057; delta=5e-3;
TH=1273; Tsat=373; epsilon1=0.7; epsilon2=0.9; n=1.76; diam=0.2;

coeff=sqrt(g*(rhol-rhog)/sigma);
emissiv=1/(1/epsilon1+1/epsilon2-1);

y=-k*(TH-TB)/delta + sigmarad*(TH^4-TB^4)*emissiv - mu*hfg*coeff*((cp*(TB-Tsat))/(csf*hfg*Pr^n))^3;