This is a spreadsheet giving a numerical integration calculated power output estimate for different steam engine designs.

## Numerical Calculations

The numerical integration performed by the spreadsheet is based off the idea of work diagram construction and the more modern approach of numerical integration. This section is meant to be a brief overview of the methodology employed in the calculations to attempt to approximate the performance of various engine designs. It should be noted that the calculations are meant to be of a comparative, rather than an absolute nature. That is, they should accurately reflect whether one bore size would be more efficient than some other bore size, but will not necessarily tell you the specific power output or steam consumption in absolute terms.

The calculation is modeled after the finite time difference method of numerical integration using the simple Euler's method to solve the state function. The spreadsheet is therefore organized like a computer program written in a traditional computer language like C or Java might use. Each subsection below will detail one portion of the calculation process, which can be thought of as analogous to a subroutine in the overall computer program. The reason for implementing the calculation in a spreadsheet rather than the analogous computer program was to make the calculations easily accessible to a much wider audience, at the cost of making the initial implementation a bit more time consuming.

### Parameters Sheet

The most important page from the users perspective is the Parameters sheet, the first sheet in the overall spreadsheet. This page allows the user to adjust the various engine parameters (stroke length, bore, valve timings, etc) and then see the results in the form of horsepower and steam consumption estimates, as well as graphs of the thermodynamic cycle undertaken by the engine. The parameters will not be documented in detail here as they are subject to change as the document evolves and are documented inside the spreadsheet itself. The cells with green background are input parameters, intended for users to change; the gray and white cells are labels or values calculated by the spreadsheet and are therefore not user editable. The spreadsheet is set to protected to prevent users from accidentally selecting or changing any cells other than the green ones, however there is no password required to turn off the protection (go to Tools->Protect Document and uncheck the 'Sheet' button) so anyone who really wants to can see/change how the calculation works.

The graphs at the bottom of the parameters page are the easiest way to see the performance of a given engine design. The first is the 'PV Diagram' for the cylinder pressure and volume; the computation of this curve is the primary purpose of the numerical integration. Since work done by the engine manifests itself as the enclosed area on this cycle, the PV diagram is the best diagnostic in determining the engine output. With the default numerical integration parameters, the diagram should show several full cycles of the engine with ~100 integration steps per cycle. This allows one to see how the engine will start up (from initially atmospheric pressure) to its steady state performance; allowing one to easily diagnose incomplete intake/exhaust cycles and their impact on overall performance. A properly functioning engine should show a single path around the PV cycle with all the simulated cycles following the same path. This means that the simulation has reached its final steady state behavior.

The next graph is the Temperature Pressure diagram. Like the PV diagram it shows another view of the thermodynamic cycle undertaken by the engine, however less effort goes into determining the values on this calculation. The reason for this is that the only thing the diagram is really used for is as a check on the final engine design that the steam will remain energetic enough to not condense in the cylinder. To calculate this the temperature is assumed to fall from its starting temperature (the superheated temperature parameter) and undergo cooling only due to adiabatic expansion. The key feature to observe on the graph is that the decompression line does not cross over to the left side of the non-moving water vapor pressure line (i.e. the temperature at which steam will re-condense back into water. The temperature calculation is done post-hoc to the actual work calculations used to arrive at the PV diagram since the cylinder temperature does not affect the work produced (unless it crosses the critical temperature to condense and we are building the engine to avoid that). After a given engine design has been decided as a final candidate, the user can increase/decrease the superheating temperature to see how hot the steam needs to be for that engine to avoid condensation. Note that since the adiabatic cooling starts from the user supplied temperature, the actual steam used to feed the engine would have to be slightly hotter as some mixing with the leftover exhaust gases in the cylinder will lower the starting temperature. For now the exact details of this mixing induced cooling is beyond the scope of the calculation for this spreadsheet.

The final graph shows the time derivatives of the work output and the steam intake, normalized to the engine cycle phase (which goes from 0.0 to 1.0 on the x-axis). The blue graph (the derivative of the work function) is shown so that you can see how much work is produced during the decompression stroke (the positive spike from 0.5 to 1.0), verses how much is consumed during the compression stroke (the negative spike from 0.0 to 0.5). Also shown is the derivative of steam consumption (i.e. the flow rate of steam into the cylinder). For now, only the steam intake is graphed (since the engine should fully exhaust all of this at steady state) but it is shown as the derivative of the flow through the valve so it can be negative for some of its operation indicating the cylinder is pushing high pressure steam backwards into the steam chest. (Note that the two curves are drawn on different axis so the 0 points are different).

### Steam Tables Sheet

This is a sheet containing tabular data as well as analytical approximations of the thermodynamic properties of water/steam. It is used by the other pages to work out the boiler pressure, as well as the density of the vapor for steam usage calculations and flow rate calculations. Since the data is the physical behavior of water, the user should never have to interact with this sheet as its contents are fixed by the laws of physics.

### Integration Parameters Sheet

This is a small page of parameters that are specific to the numerical integration, rather than the engine design. The most important thing on this sheet is the number of simulation steps per cycle. Since the spreadsheet uses a fixed number of rows in the integration, more steps per cycle will mean that you end up simulating fewer cycles. The other option you can set is the adiabatic flow rate constant: at 1.0 the valves work normally, at 2.0 they flow twice as fast, etc. Setting the adiabatic constant term to 0.0 will disable the flow through the valves so you can examine the engine behavior during compression without letting in any steam.

### Integration Sheet

The last sheet is the numerical integration calculation itself. The whole calculation is far to complicated to do in a single step so it is broken down into multiple sub calculations to determine various easy to compute things, like the cylinder position, or whether the valve should be open or not given the current control scheme. Once these baseline things are established in the first few columns, they are used as inputs to other, more sophisticated calculations like pressures and flow rates in the subsequent columns. The list below contains the "pseudo code" for the algorithm. Each bulleted item corresponds to a calculation performed in one (or a few) columns, with the major computation steps being shown as groups of columns separated in between by empty columns. As a general rule I try to use only the columns in one small cluster, or the clusters to the left of a given column as inputs for the calculation in that column. This means that the leftmost columns are the easiest things to compute and are assumed not to depend on anything to the right in this sheet. Once these are established as absolute, they are used as inputs in the columns to the right... and so on. The graphs, and many of the horsepower and steam consumption figures on the Parameters sheet are mostly just graphs of various columns in this large numerical calculation, or are individual figures found by doing row lookups (such as how much steam has been consumed in total at the end of the 1st cycle, 2nd cycle...).

• Start with a column stating the angular position ($\theta$) of the crankshaft/flywheel and let it increase from 0 to $2 \pi n$ where $n$ is the number of engine cycles to be simulated. Since the flywheel position increases steadily, the effect of feedback torque on the engine speed is ignored (i.e. we are simulating the steady state behavior of a machine with a sufficiently heavy flywheel).
• From the angular position of the crankshaft calculate the position and velocity of the piston in the cylinder starting at BDC and moving up from there. Currently the position is just given as $x(t) = S*\sin(\theta)$, however future versions will implement the true motion of a piston driven by a crasnkshaft.
• Now that piston position is established, on a bump valve design the valve states are also fully determined so their values are computed and stored, both in terms of the linear distance the valve is "open" and the cross sectional area of the valve available for flow at the given time step.
• Now using the volume of the cylinder and bootstrapping the pressure with $P(t=0) = Exhaust press.$ to set the initial conditions. Once this is done the adiabatic constant defined as $C = PV^\gamma$ with $\gamma \approx 1.31$ for steam is computed. From then on the, if the valves are closed the adiabatic constant remains unchanged from the previous step and the new cylinder volume is used to compute the pressure using $P(t) = \frac{C(t-dt)}{V^\gamma}$ to work out the new pressure.
• If the valves are open, an additional calculation takes place which adds or subtracts from the current adiabatic constant, based on how much gas flows into or out of the cylinder. To compute $\Delta C$ the flow rate is calculated (see next bullet point for details) and the pressure and volume of the resulting flow during the next time step are used to compute $C$ as if the flowing material were a closed system separate from the engine simulation using the same adiabatic constant equation used to calculate the adiabatic constant in the first time step. This new adiabatic constant for the flow is then added to (for intake flows) or subtracted from (for exhaust flows) the current adiabatic constant for the cylinder to produce the constant for the next time step. That constant is then used to calculate the next cylinder pressure and so on...
• To compute the flow volume through the intake and exhaust valves they are modeled as an opening in a thin sheet of a given cross sectional area and a pressure difference across the opening using orifice plate theory. The pressure difference if found by taking the difference of the current cylinder pressure (as above) and the intake or exhaust pressure taken from the parameters sheet. The valve area is taken from the leftmost columns as it is one of the first things computed. Steam density is assumed to be constant during the cycle and its value is taken to be the density of steam at the operating pressure and temperature of the boiler. This is roughly correct for the steam flow through the intake valves, but due to the expansion stroke this value is too high for the one in the exhaust valve flow (meaning the exhaust valves can be smaller than the size estimated here). The reason for not using a better steam density model is simply the author's lack of an adequate analytical approximation for $\rho_{steam}(P, T)$. A better model for this would remove some of the uncertainty in the model and help to establish absolute quantitative numbers for power output and steam consumption, rather than the qualitative results produced by the present model.
• Finally, now that the pressure of the cylinder has been established by adiabatic compression and corrected for fluid flow, the working pressure and change in cylinder position during this time step are used to calculate $dw$ the "infinitesimal" amount of work done during this time step. This work function can then be plotted as on the parameters page, and also summed to find the total work done during each of the complete operating cycles. Knowing the work per cycle and multiplying by the engine speed in RPM gives us the HP rating of the engine.