# Analyze the Stirling Cycle with Real-World Materials

Examine the Stirling cycle for 10 moles of air.

 In[1]:= XmAir = 10 ChemicalData["Air", "MolarMass"]
 Out[1]=

Calculate the isothermal expansion and compression curves at 200 °C and 100 °C, respectively.

 In[2]:= XisothermalExpansion = Transpose[{ ThermodynamicData["Air", "SpecificVolume", {"Temperature" -> Quantity[200, "DegreesCelcius"], "Pressure" -> Quantity[Range[100, 200] 1000, "Pascals"]}] mAir, Quantity[Range[100, 200] 1000, "Pascals"]}];
 In[3]:= XisothermalCompression = Transpose[{ ThermodynamicData["Air", "SpecificVolume", {"Temperature" -> Quantity[120, "DegreesCelcius"], "Pressure" -> Quantity[Range[100, 200] 1000, "Pascals"]}] mAir, Quantity[Range[100, 200] 1000, "Pascals"]}];

Create the diagram.

 In[4]:= XListLinePlot[{isothermalExpansion, isothermalCompression}, Axes -> False, PlotRange -> All, FrameLabel -> Automatic, Frame -> True, GridLines -> {{0.2, 0.3}, {}}]
 Out[4]=

Calculate mechanical work through .

 In[5]:= XW = Quantity[ NIntegrate[ Evaluate[ Interpolation[QuantityMagnitude[isothermalExpansion]][V]], {V, 0.2, 0.3}] - NIntegrate[ Evaluate[ Interpolation[QuantityMagnitude[isothermalCompression]][V]], {V, 0.2, 0.3}], "Joules"]
 Out[5]=

Find maximal and minimal pressures for expansion and compression.

 In[6]:= X{pMax1, pMin1} = Quantity[Interpolation[QuantityMagnitude[isothermalExpansion]][{0.2, 0.3}], "Pascals"]
 Out[6]=
 In[7]:= X{pMax2, pMin2} = Quantity[ Interpolation[QuantityMagnitude[isothermalCompression]][{0.2, 0.3}], "Pascals"]
 Out[7]=

Calculate the change in internal energy.

 In[8]:= XU[T_, p_] := ThermodynamicData["Air", "InternalEnergy", {"Temperature" -> T, "Pressure" -> p}]*mAir; \[CapitalDelta]U = U[ Quantity[200, "DegreesCelcius"], pMax1 ] - U[ Quantity[120, "DegreesCelcius"] , pMin2]
 Out[8]=
 Out[8]=

Using work and the change in internal energy, calculate efficiency.

 In[9]:= XW/\[CapitalDelta]U
 Out[9]=

Compare with maximal possible Carnot efficiency.

 In[10]:= X1. - UnitConvert[Quantity[120, "DegreesCelcius"], "Kelvins"]/ UnitConvert[Quantity[200, "DegreesCelcius"], "Kelvins"]
 Out[10]=

## Mathematica

Questions? Comments? Contact a Wolfram expert »