Quiet@Remove["Global`*"];
$HistoryLength=2;
SetDirectory[NotebookDirectory[]];

kJ=1000;
EnergyOfOneBeer=(1850/2)kJ;
m=85;
g=9.81;
Solve[EnergyOfOneBeer==h*m*g]

Quiet@Remove["Global`*"];
cH2O=4186;
ΔT=1;
g=9.81;
m=1;
Solve[m*cH2O*ΔT==h*m*g]

Quiet@Remove["Global`*"];
sol=DSolve[{T''[x]==0,T[0]==T1,T[d]==T2},T[x],x][[1]];
temperature=T[x]/.sol;
Print["T[x]= ",temperature]
q=-λ*Grad[temperature,{x,y,z}]
qx=q[[1]]
Print["qx= ",qx]
P=S*qx;
Print["P= ",P]

(Sqrt[1000.])

Quiet@Remove["Global`*"];
(*convection, wire of 400kV overhead lines *)
Ts=70;
Ta=0;
cm=0.01;
d=5cm;
Dchar=d;
Tdef=((Ts+Ta)/2);
Print["Tdef= ",Tdef," °C"]
λ=29.99*(10^(-3));

(*https://www.engineersedge.com/heat_transfer/nusselt_number_13856.htm*)

v=20;
ν=16.45*(10^(-6));(*konematic viscosity*)
Rey=((v*Dchar)/ν);
Print["Re= ",Rey," -"];
Pr=0.7;
Nu=0.027*(Rey^0.805)*(Pr^(1/3));
Print["Nu= ",Nu," -"];
αForced=Nu*(λ/Dchar);
Print["αForced= ",αForced," W.(m^(-2)).(K^(-1))"]

NuNew[re_]:=0.3+((0.62*(re^0.5)*(Pr^(1/3)))/(((1+((((0.4/Pr)))^(2/3))))^(1/4)))*(((1+((((re/282000)))^(5/8))))^(4/5));
Print["NuNew= ",NuNew[Rey]," -"]

length=1;
S=Pi*d*length;
Pconv=S*αForced*(Ts-Ta);
Print["Power of convection from 1 meter of wire= ",Pconv," W"]

ϵ=0.8;
σ=5.67*(10^(-8));
Prad=S*ϵ*σ*((((Ts+273))^4)-(((Ta+273))^4));
Print["Power of radiation from 1 meter of wire= ",Prad" W"]

PtotalForced=Pconv+Prad;
Print["Total power from 1 meter of wire= ",PtotalForced" W"]

ρAl20=2.65*(10^(-8));
αR=0.00429;
ρAl=ρAl20*(1+αR*(Ts-20));
ScrossSection=((Pi*(d^2))/4);
R=((ρAl*length)/ScrossSection);
Print["R1m= ",(10^6)R," μΩ"]
electricalLoses=R*(current^2);
Solve[electricalLoses==PtotalForced]

(*start of free convection*)
Series[(1/(1+β*ΔT)),{ΔT,0,3}]

β=(1/(Tdef+273.));
g=9.81;
Gr=((β*(Ts-Ta)*g*(Dchar^3))/(ν^2));
Print["Gr= ",Gr," -"];
Print["Gr*Pr= ",Gr*Pr," -"]
NuFree=0.48*(((Gr*Pr))^0.25);
Print["NuFree= ",NuFree," -"]

αFree=NuFree*(λ/Dchar);
Print["αFree= ",αFree," W.(m^(-2)).(K^(-1))"]

PconvFree=S*αFree*(Ts-Ta);
Print["Power of free convection from 1 meter of wire= ",PconvFree," W"]

PtotalFree=PconvFree+Prad;
Print["Power of radiation from 1 meter of wire= ",Prad" W"]
Print["Total power from 1 meter of wire free= ",PtotalFree" W"]

Solve[electricalLoses==PtotalFree]

Quiet@Remove["Global`*"];
(*heater in a room*)
Ts=50;
Ta=20;
a=0.5;
b=1;
Dchar=a;
Tdef=((Ts+Ta)/2);
Print["Tdef= ",Tdef," °C"];
λ=26.99*(10^(-3));
ν=16.45*(10^(-6));(*konematic viscosity*)
β=(1/(Tdef+273.));
g=9.81;
Gr=((β*(Ts-Ta)*g*(Dchar^3))/(ν^2));
Pr=0.7;
Print["Gr= ",Gr," -"];
Print["Gr*Pr= ",Gr*Pr," -"]

Ra=Gr*Pr;
Nu=0.68+((0.67*(Ra^0.25))/(((1+((((0.492/Pr)))^(9/16))))^(4/9)));
α=Nu*(λ/Dchar);
NuMIT=(((0.825+((0.387*(Ra^(1/6)))/(((1+((((0.492/Pr)))^(9/16))))^(8/27)))))^2);
αMIT=NuMIT*(λ/Dchar);
{α,αMIT}

NuOld=0.135*(((Gr*Pr))^(1/3));
αOld=NuOld*(λ/Dchar);
{α,αMIT,αOld}

S=2*a*b;
Pconv=αMIT*S*(Ts-Ta);
Print["Power of convection= ",Pconv," W"]

ϵ=0.8;
σ=5.67*(10^(-8));
Prad=S*ϵ*σ*((((Ts+273))^4)-(((Ta+273))^4));
Print["Power of radiation from the heater= ",Prad" W"]

Print["Total power from the heater= ",Prad+Pconv," W"]

Quiet@Remove["Global`*"];
(*heat exchangers*)
Plot[(1-((2*(ξ-1))/((ξ+1)*Log[ξ])))*100,{ξ,0.5,2},AxesLabel->{"ξ[-]","δ[%]"},GridLines->Full,PlotStyle->{Red,Thickness[0.006]}];
Limit[((Δ1-Δ2)/(Log[((Δ1)/(Δ2))])),Δ2->Δ1]

(*parallel or counter flow?*)
volume=2.5*5*11;
Print["Volume= ",volume," (m^3)"];
Tinside=20;
Toutside=-12;
cAir=1000;
density=(349./(Tinside+273));
Print["ρAir= ",density," (kg/(m^3))"]
timeOfAirExchange=3*3600;
massOfAir=volume*density;
Print["mAir= ",massOfAir," kg"]
energy=massOfAir*cAir*(Tinside-Toutside);
Print["energy for one exchange= ",(energy/(10^6))," MJ"]
power=(energy/timeOfAirExchange);
Print["power for ventilation= ",power," W"]
mDot=(massOfAir/timeOfAirExchange);
Print["mDot= ",mDot," kg/s"]

(*counter flow*)
T1c=Tinside;
T3c=Toutside;
powerC[kS_]:=Module[{eqns,Δ1,Δ2,P,T4,T2},
Δ1=T1c-T4;
Δ2=T2-T3c;
eqns={P==kS*((Δ2+Δ1)/2),P==mDot*cAir*(T1c-T2),
mDot*cAir*T1c+mDot*cAir*T3c==mDot*cAir*T2+mDot*cAir*T4
};
P/.Solve[eqns][[1]]
];
(*parallel flow*)
powerP[kS_]:=Module[{eqnsLinear,Δ1,Δ2,P,T4,T2,solLin,forFindRoot,eqns},
Δ1=T1c-T3c;
Δ2=T2-T4;
eqnsLinear={P==kS*((Δ2+Δ1)/2),P==mDot*cAir*(T1c-T2),
mDot*cAir*T1c+mDot*cAir*T3c==mDot*cAir*T2+mDot*cAir*T4
};
solLin=Solve[eqnsLinear][[1]];
forFindRoot=solLin/.Rule->List;
eqns={P==kS*((Δ1-Δ2)/(Log[((Δ1)/(Δ2))])),P==mDot*cAir*(T1c-T2),
mDot*cAir*T1c+mDot*cAir*T3c==mDot*cAir*T2+mDot*cAir*T4
};
Quiet[Chop[P/.FindRoot[eqns,forFindRoot]]]
];

dataC=Table[{kS,powerC[kS]},{kS,0.1,100,1}];
dataP=Table[{kS,powerP[kS]},{kS,0.1,100,1}];
plC=ListPlot[dataC,AxesLabel->{"kS [W/K]","P [W]"},Joined->True,PlotStyle->{Red,Thickness[0.006]},GridLines->Full];
plP=ListPlot[dataP,AxesLabel->{"kS [W/K]","P [W]"},Joined->True,PlotStyle->{Blue,Thickness[0.006],Dashed},GridLines->Full,PlotRange->All];
Show[plC,plP,PlotRange->{Automatic,{0,power}}]

Quiet@Remove["Global`*"];
T1=80;
T2=40;
c=4186.;
T3=10;
Δ1=T1-T4;
Δ2=T2-T3;
m1=0.1;
m2=0.08;
eqns={P==kS*((Δ1-Δ2)/(Log[((Δ1)/(Δ2))])),P==m1*c*(T1-T2),m1*c*T1+m2*c*T3==m1*c*T2+m2*c*T4};
Quiet@Solve[eqns][[1]]

Quiet@Remove["Global`*"];
(*training for exam*)
cm=0.01;
d=3cm;
v=0.5;
Ts=0;
Ta=4;
Dchar=d;
Tdef=((Ts+Ta)/2);
Print["Tdef= ",Tdef," °C"]

dataλ={{0.01,555.75},{10,578.64}};
λ=(10^(-3))Fit[dataλ,{1,T},T]/.T->Tdef;
Print["λ= ",λ," W.(m^(-1)).(K^(-1))"]

ν=0.0000016761;
Rey=((v*Dchar)/ν);
Print["Re= ",Rey," -"]

dataPr={{0,13.6},{5,11.2}};
Pr=Fit[dataPr,{1,T},T]/.T->Tdef;
Print["Pr= ",Pr," -"]

Nu=0.683*(Rey^0.466)*(Pr^(1/3));
Print["Nu= ",Nu," -"]

α=Nu*(λ/Dchar);
Print["α= ",α," W.(m^(-2)).(K^(-1))"]

S1m=Pi*d*1;
P1m=α*S1m*(Ta-Ts);
Print["maximal power from 1m of tube is ",P1m," W"]

(*wire in oil, free convection*)
Quiet@Remove["Global`*"];
Ts=90;
Ta=60;
mm=0.001;
d=3mm;
Dchar=d;
Tdef=((Ts+Ta)/2);
Print["Tdef= ",Tdef," °C"]

ν=30*(10^(-6));
λ=0.31;
Pr=795;
β=0.0007;
(*for air easily β=(1/(Tdef+273))*)
g=9.81;
Gr=((β*(Ts-Ta)*g*(Dchar^3))/(ν^2));
Print["Gr= ",Gr," -"]

Ra=Gr*Pr;
Nu=0.68+((0.67*(Ra^0.25))/(((1+((((0.492/Pr)))^(9/16))))^(4/9)));
α=Nu*(λ/Dchar);
NuMIT=(((0.825+((0.387*(Ra^(1/6)))/(((1+((((0.492/Pr)))^(9/16))))^(8/27)))))^2);
αMIT=NuMIT*(λ/Dchar);
{α,αMIT}

(*https://ahtt.mit.edu/*)

length=1;
S1m=Pi*d*length;
power=αMIT*S1m*(Ts-Ta);
elconductivity=50*(10^6);
resistance=(length/(elconductivity*((Pi*(d^2))/4)));
Print["R= ",resistance," Ω"];
Solve[resistance*(i^2)==power]

Quiet@Remove["Global`*"];
$HistoryLength=2;
SetDirectory[NotebookDirectory[]];Off[Solve::"ifun",General::"spell1"]

p1[t_]:=Exp[((1500.3+23.5*t)/(234+t))];
t=100;
{p1[t],p2[t]}

Plot[p1[T],{T,20,100},AxesLabel->{"T[°C]","p[Pa]"},GridLines->Full]

(*heating the air*)
Quiet@Remove["Global`*"];
TA=0;
fiA=0.8;
xA=3.0592;(*grams H2O per 1 kg of dry air*)
hA=7.6508;(* enhtalpy, kJ/kg*)
xB=xA;
fiB=0.174;
hB=31.01;
TB=23.;
neededHeat=hB-hA;
Print["needed heat fo heating= ",neededHeat," kJ/kgda"];
V=6*3*11;
density=(349./(TB+273));
Print["Volume= ",V," (m^3),  density= ",density," kg/(m^3)"]
mass=V*density;
Print["there is ",mass," kg of dry air in the room"]
timeOfAirExchange=3*3600.;
energy=mass*neededHeat;
Print["needed energy for one air exchange= ",(energy/1000)," MJ"]

power=(energy/timeOfAirExchange);
Print["power for heating= ",power," kW"]

Quiet@Remove["Global`*"];
(*drying the air*)
fiA=0.8;
TA=25;
xA=16.1866;
hA=66.4757;
fiC=0.5;
TC=25;
xC=10.0189;
hC=50.7671;
xB=xC;
fiB=1;
TB=13.865;
hB=39.3155;
heatFromAir=hA-hB;
Print["energy of cooling= ",heatFromAir," kJ/kgda"];
energyIntoAir=hC-hB;
Print["energy of heating= ",energyIntoAir," kJ/kgda"];

Quiet@Remove["Global`*"];
(*increasing relative moisture*)
TA=15;
fiA=0.4;
xA=4.2737;
hA=25.9561;
TB=23;
fiB=0.7;
xB=12.4871;
hB=54.9874;
hC=hB;
xC=xA;
TC=43.5218;
energy=hC-hA;
Print["energy of heating= ",energy," kJ/kgda"];
water=xB-xA;
Print["we need ",water," g of H2O for each kg of dry air"]

Quiet@Remove["Global`*"];
(*mixing and heating*)
TA=25;
fiA=0.8;
xA=16.1866;
hA=66.4757;
TB=10;
fiB=0.7;
xB=5.3946;
hB=23.6906;
TC=20;
fiC=0.6;
xC=8.8559;
hC=42.6736;
xD=xC;
solmA=Solve[xA*mA+xB*(1-mA)==xD][[1]]
hD=mA*hA+(1-mA)*hB/.solmA;
Print["enthalpy after mixing= hD= ",hD," kJ/kgda"];
neededHeat=hC-hD;
Print["energy of heating= ",neededHeat," kJ/kgda"];

Quiet@Remove["Global`*"];
(*drying the air*)
TA=25;
fiA=80;
xA=16.1866;
hA=66.4757;
TB=20;
fiB=60;
xB=8.8559;
hB=42.6736;
TC=12;
hC=34.4717;
heatFromAir=hA-hC;
Print["energy of cooling= ",heatFromAir," kJ/kgda"];
energyIntoAir=hB-hC;(*kJ/kgda*)
Print["energy of heating= ",energyIntoAir," kJ/kgda"];
Δx=xA-xB;
Print["from each kg of dry air ",Δx," grams of water condensates"]
volume=10000;
density=(349./(TB+273));
mass=volume*density;
totalEnergyForCooling=mass*heatFromAir;
Print["energy for cooling= ",(10^(-6))*1000*totalEnergyForCooling," MJ"]
time=1*3600;
power=((totalEnergyForCooling*1000)/time);
Print["cooling power= ",power," W"]