Engineering 41

Lab 5: Characterization of a Gamma-Type Stirling Engine

Aron Dobos, Adem Kader, Brian Park, Brendan Bond

November 22, 2004

Abstract

In this laboratory, the performance of a small Gamma-type Stirling engine was characterized. Stirling engines rely on thermodynamic heat transfer in order to drive the engine. The engine was tested under load-less and numerous loaded conditions, at varying power inputs. The average no-load speed and torque under average power input (50W) were 0.0896m/s and 2.99W. For medium power (39.2W), the values were 0.0777m/s and 2.30W, and for low power (28.5W), the values were 0.0484m/s and 0.61W. The stirling engine was determined to be very inefficient due its design and specifications. Improvements and suggestions were developed to make the engine more efficient.

Introduction

The Stirling engine, originally known by its originator as the economizer, was invented in 1816 by Robert Stirling. It is an external combustion engine that operates by a sequence of pressure changes of a gas in a sealed enclosure that is converted into mechanical work. Because combustion is external, it is difficult to rapidly change the power output of a Stirling engine since it takes time for temperature changes to propagate throughout the engine. The result is that more of the Stirling engine’s energy is stored in the engine itself, meaning that Stirling engines tend to be larger for the same power output than internal combustion engines. While these facts make the Stirling engine a less attractive option for automobiles and other applications, the engine’s design results in a theoretical thermal efficiency equal to that of the Carnot cycle.

Each of the four stages of the ideal Stirling cycle is a fully reversible process, two of which involve heat transfer, one involving gas expansion, and the last the compression of the gas. More work is done in the expansion stage compared to the compression stage, and the difference is extracted from the system as work. The cycle operates between a maximum temperature at the hot end of the engine and the minimum temperature obtained at the cold end. The gas within is moved between the hot and cold regions by the displacer, resulting in a wave of pressure differences due to the temperature variation of the gas.

The first stage is the isothermal expansion stage (1 à 2). Heat from the external combustion is added to the system causing the working gas (air in this case) to expand and do work. The second stage is the isochoric displacement or cooling stage (2 à 3), in which the displacer piston shuttles the gas from the hot end to the cold end of the engine. Here, the working gas decreases in pressure and temperature. If a regenerator is implemented, the excess heat is absorbed by the regenerator to be used later in “pre-heating” the working gas in the final sequence. In the isothermal compression (3 à 4) stage, the cooled gas is compressed by the power piston in the compression region as heat is rejected to the cold region. During the last stage, known as isochoric displacement or heating (4 à 1), the gas is moved from the cold region to the hot region, thereby increasing the temperature and thus pressure of the working gas. A regenerative engine releases its stored heat from the second stage to help increase the temperature of the gas also. Thus, a fully regenerative Stirling engine can approach the efficiency of the Carnot cycle since all the heat-injection and heat-extraction sequences occur at isothermal places in the cycle. Heat loss stemming from imperfect regeneration is one reason Stirling engines cannot in practice achieve full Carnot efficiency. The P-v diagram for the ideal Stirling cycle is shown below.

Figure 1. Ideal Stirling Cycle P-v Diagram. [Haywood 3]

The ideal Stirling cycle is approximated in the design and mechanical operation of a simple gamma engine by a single continuous space with two cylinders. The power piston in the power cylinder converts the differential pressure waves into sinusoidal work output, while the displacer in the heat exchange unit shuttles the working gas between the hot and cold regions of the engine. The power piston is designed for full air-tight operation, while the displacer is laced with a 2mm radial clearance through which the working gas can moved around the engine. It is worth noting that ideally the working gas pressure at any given time is constant throughout the engine is constant. The displacer leads the power piston by 90 degrees phase angle, and due to the sinusoidal piston motion, fails to reach the “corners” of the ideal Stirling cycle in which stages 2 and 4 ideally occur at constant volume. The approximation is good enough, however, since the piston spends more than 50% of the cycle in the ideal isochoric positions of either bottom dead center (BDC) or top dead center (TDC). The machining parameters for the power piston area and stroke are 0.701 sq. in. and 1.5 in, respectively. The displacer also has a 1.5 in stroke, but a 1.16 in diameter.

Figure 2. Gamma-Type Stirling Engine

(http://www.ent.ohiou.edu/~urieli/stirling/engines/gamma.html)

As seen in the above schematic, the gamma engine imposes a full separation between the heating/cooling units and the power unit. This implementation results in a potentially large dead volume that may introduce significant losses due to gas flow around sharp corners and heat loss through the additional surface areas. The specifications give 4.7 cm^3 for the dead volume, and 1/8^th in. for the internal dead volume piping diameter. The Stirling engine characterized in the laboratory proved entirely bereft of a regenerative mechanism, thereby immediately reducing its maximum theoretical thermal efficiency. Heating was provided by an electric heating coil, and cooling simply by a properly thermally greased heat sink to the surroundings. The rated hot and cold temperatures are 400 K and 310 K respectively.

Laboratory Procedures

Heat input to the gamma-type Stirling engine was accomplished by an electric heating coil whose power was readily calculated knowing the input voltage and current. A temperature probe was also inserted into the hot end to record the actual operating temperature for each power input level. The cooling side involved nothing more than the standard heatsink at the ambient temperature. The internal pressure inside the engine was measured using a basic pressure transducer connected to a data acquisition system. The engine shaft speed was determined by an optical tachometer that recorded a voltage spike when the power piston reached top dead center.

To determine the hoisting power of the engine, various weights were loaded in different combinations resulting in weight differentials of 46, 71, 89, 121, 143, and 347 g. The pvrealtime.m MATLAB script recorded real-time pressure and tachometer readings for later analysis. Each weight differential was hoisted at three different input power levels, after allowing the engine to achieve equilibrium steady state operation.

Selected P-V diagrams for zero load conditions at the three power levels are included below.

Power Piston and Displacer Phase Considerations

The various power piston and displacer position relationships are shown below at several 45 degree intervals around one cycle.

Power Piston	Displacer	System V	Power Piston	Fluid
-45 à+ 45	45 à 135	no net change	first doing then absorbing work	cooled down
+45à +135	135 à -135	decrease	absorbing work
+135à-135	-135 à -45	no net change	First absorbing then doing work	heated up
- 135à -45	-45 à 45	increase	doing work

Looking at the P-V graph above and Table 1, the line from 1 to 2 represents work done by the piston and expansion of the fluid; this corresponds to the section where the power piston goes from -135 to -45. In this section the piston does work and the fluid expands. (Isothermal expansion)

The line from 2 to 3 in the above graph corresponds to the section where the power piston goes from -45 to 45. Here, the power piston first does a little bit of work and then absorbs some work. In theory, the main thing happening in this section is the cooling down of the fluid but in practice, this section is very short because actually, the parts from 1 to 2 and 2 to 3 are combined in that the fluid starts cooling down while the work is being done. This is why in the measured P-V graph, we observe that the upper part of the graph is concave up and that we can hardly observe the 2 to 3 section. (Cooling stage)

The line from 3 to 4 corresponds to the section where the power piston goes from 45 to 135. In this section the system volume decreases and the piston absorbs work. (Isothermal compression)

In the last part (3 to 4), the fluid in the system is heated up and this corresponds to the section where the power piston goes from 135 to -135 in Table 1. In theory there is no change in volume, however in practice the volume increases a little bit and then decreases as illustrated in the figure below. Again, in the real application of the system, the last two parts are combined; the fluid starts heating up while the volume is decreasing. (Heating

stage)

Varying the phase lead of the Displacer relative to the Power Piston the system is a possible way to increase the efficiency of the system; however it is not very smart to do it.

As we mentioned before, in practice the isothermal expansion and the cooling stages overlap and this is also the case with the isothermal compression and the heating stages. We can take advantage of this overlapping and heat up the fluid to a higher temperature while in the compression stage which would then push the piston with more force, resulting in a greater torque. The down-side of this change is the early deceleration of the power piston when entering the cooling stage. Namely, this change would make the flywheel turn faster at some stages of the rotation and slower in other stages and this is why we wouldn’t want to make such a change in real applications. Imagine a car wheel having a non-uniform angular velocity over a period of rotation; this would make the passengers very dizzy. Even in our experimental settings, this change would result in a very inefficient mechanism because the loads would be pulled with varying speeds and acceleration. It could also stall the engine easier because when the flywheel jerks the load up, then the load would jump a little bit (we have just wasted some of the energy by not keeping the load at the max height it jumped to), and when it goes back down, it would actually pull the flywheel down. Even as it is, the velocity of the flywheel is not constant throughout the cycle and making the velocity vary even more would be disastrous.

Results

The pressure and tachometer data recorded using the matlab script pvrealtime.m was analyzed using the matlab script analyze.m (see appendix B). Analyze.m calculated the integral of the p-V diagram and the number of revolutions per minute for a few selected “clean” sterling cycles. Knowing the revolutions per minute, RPM, the radius of the shaft, r, and the load, L, the torque, τ, and hoisting power, P, can be calculated using the following equations:

(Torque Eqn)

where g is the gravitational constant

(Hoisting Power Eqn)

The equations used to calculate the thermal and mechanical efficiencies of the system are as follows:

(Thermal Eff.)

Where, Q_in is unknown, but can be approximated by the W_in to the heating elements

And W_out is equal to the integral of the p-V diagram multiplied by the cycles per second of the system.

(Mechanical Eff.)

The results of the previously mentioned calculations are shown in tables R1 through R3 for the varying watt inputs and varying loads. Figures R1 through R3 contains Power-Torque versus shaft speed and Hoisting Power- Mechanical Efficiency versus shaft speed plots for 28.5, 39.2 and 50.2 watt inputs respectively.

Table R1. Results from Varying Loads with 28.5 Watt Input

Load (g)	Shaft Speed (m/s)	p-V Integral (W)	Thermal Efficiency of Sterling Engine, ηth	Torque (N-m)	Hoisting Power (W)	Mechanical Efficiency, ηm
0	0.0484	0.61	1.22%	0.0000	0.0000	0.00%
46	0.0431	0.52	1.03%	0.0011	0.0194	3.75%
71	0.0419	0.50	0.99%	0.0017	0.0291	5.88%
89	0.0394	0.47	0.94%	0.0021	0.0344	7.30%

Table R2. Results from Varying Loads with 39.2 Watt Input

Load (g)	Shaft Speed (m/s)	p-V Integral (W)	Thermal Efficiency of Sterling Engine, ηth	Torque (N-m)	Hoisting Power (W)	Mechanical Efficiency, ηm
0	0.0777	2.30	5.9%	0.0000	0	0.0%
71	0.0683	1.21	3.1%	0.0017	0.0476	4.4%
89	0.0662	1.08	2.8%	0.0021	0.0578	5.6%
121	0.064	1.01	2.6%	0.0028	0.0759	7.9%
143	0.0616	0.81	2.1%	0.0033	0.0864	10.9%
374	0.044	0.40	1.0%	0.0088	0.1617	41.2%

Table R3. Results from Varying Loads with 50.2 Watt Input

Load (g)	Shaft Speed (m/s)	p-V Integral (W)	Thermal Efficiency of Sterling Engine, ηth	Torque (N-m)	Hoisting Power (W)	Mechanical Efficiency, ηm
0	0.0896	2.99	5.95%	0.0000	0.0000	0.00%
46	0.0708	1.51	3.01%	0.0011	0.0319	2.11%
71	0.0854	2.71	5.39%	0.0017	0.0594	2.19%
89	0.0855	2.82	5.61%	0.0021	0.0746	2.65%
121	0.0818	2.56	5.11%	0.0028	0.0970	3.78%
143	0.0751	1.64	3.27%	0.0033	0.1052	6.41%
374	0.0593	0.62	1.24%	0.0088	0.2175	34.97%

Figure R1. Power-Torque and Power-Efficiency versus Shaft Speed for 28.5 Watt Input

Figure R2. Power-Torque and Power-Efficiency versus Shaft Speed for 39.2 Watt Input

Figure R3. Power-Torque and Power-Efficiency versus Shaft Speed for 50.2 Watt Input

The ideal thermodynamic efficiency of the stirling engine can be calculated using the following equation:

(Ideal Thermal Eff.)

The resulting estimations of the ideal thermal efficiency of the stirling engine are displayed in table R4, below.

Table R4. Theoretical Thermodynamic Efficiency at Each Input

Power Input (W)	TC (Celsius)	TH (Celsius)	Ideal Thermal Efficiency
50.2	22	254	44.0%
39.2	22	223	40.5%
28.5	22	177.5	34.5%

Higher power inputs resulted in higher shaft speeds at each load. In general, the shaft speed decreased as the load applied to the shaft was increased. Higher loads produce a greater torque force that acts against the rotation of the engine and thereby decreases the engine rpm. An inverse relationship existed between the output power (pV-integral) and the hoisting power. The mechanical efficiency, which is a measurement of this relationship, increased with higher loads. Since higher loads decrease rpm, less work is needed to maintain the engine at a lower rpm and more work is needed to move the load.

Torque, power, and efficiency all followed the same relationship to the shaft speed, as can be seen from the figures above. Ideally the torque at the max power would be determined through curve fits of the power and torque. However, the stirling engine stalls at low rpm with high loads and no definite conclusions about the performance of the engine.

The stirling engine is far from ideal thermal efficiency. The difference can be a result of a variety of factors with the design and layout of the model engine. The measured temperatures at each power input were not accurate assessments of the actual temperatures within the engine. A thermocouple measured the outside glass cylinder temperature as T_h rather than the internal high temperature surrounding the displacer. The actual values of T_h are lower than the measured temperatures. In addition, T_c was assumed to be constant room temperature. The actual values of T_c are higher than the room temperature assumption. As a result, the actual ideal thermal efficiencies are lower than the ideal thermal efficiencies in Table R4.

Thermal losses due to the engine design and model contribute to the overall low performance efficiencies. Ideally, the displacer cylinder is perfectly insulated but in our model, the displacer cylinder was poorly insulated, allowing heat to leave the system, as it was encased in glass tubing for pedagogic viewing purposes. Additionally, the heat diffuser was not effectively removing heat from the engine, due to poor fin design. Therefore, the temperature differences across the displacer piston were lowered and a lower amount of work was done by the power piston. The engine also has large dead volumes which introduce significant pressure losses due to gas flow around sharp corners and heat loss through the additional surface areas.

Error during the recording of data also had significant effects on the analysis of the sterling engine’s performance. The data that was recorded by the matlab script was not limited to the period during which the sterling engine was carrying the load. As a result data points on either end of the period of load bearing may have been included in the analysis. Error also occurred within the period of loading as the p-V plot estimations of the work out of the sterling engine had variances as follows in table D1 (Note the high variances for the 39.2 watt input with loads of 0 and 71 grams, and the high variances for the 50.2 watt input with loads of 0, 46, 71, and 121 grams). Error in the loading of the sterling engine occurred as the weights used to load the engine often swayed back and forth and on occasion collided.

Table D1. Variance Calculations

Load (g)	Variance of p-V integral (28.5W)	Variance of p-V integral (39.2W)	Variance of p-V integral (50.2W)
0	0.00074	0.20391	0.08036
46	0.00086	No Data	0.2927
71	0.00086	0.19623	0.1610
89	0.00024	0.07609	0.0420
121	No Data	0.07988	0.1342
143	No Data	0.02654	0.4592
374	No Data	0.00217	0.0012

Potential Engine Improvements

It would definitely be helpful to drill the extra holes in the crank wheel to facilitate observing engine efficiency with displacer lead angles other than 90 degrees. As it is, it is somewhat difficult to assess immediately in a conceptual way how changing the phase shift might affect engine performance.

As for a regenerator, it might be worthwhile to make the displacer piston replaceable with other ones of varying shapes and materials. This may or may not be feasible, depending on how easily the end of the heating unit can be removed from the displacer cylinder. The end of the displacer rod could be threaded to easily accept a variety of displacer piston heads, some with regenerators and some without. A simple metal mesh could work as a regenerator, although frictional losses might limit any noticeable performance gains.

The ability of the flywheel to maintain relatively constant rpm considering the single-acting power cylinder could be improved by shifting more of the flywheel’s mass to the outside. The larger moment of inertia would increase the time constant of the mechanical system’s response to the ‘pulsing’ power stroke, and would thus help smooth out the shaft speed variation. It is probably not necessary to replace the flywheel with a larger one, as the increased total mass would probably load down the engine too much. Removing part of the center on a metal lathe might be sufficient to increase the moment of inertia. It seems plausible to include as part of the laboratory a few different flywheels that could be switched onto the engine, although this would be moving away somewhat from the thermodynamic focus of the lab. Also, better instrumentation for measure the engine rpm would be required, involving possibly a disk with numerous evenly spaced slits through which a light could be shined onto a phototransistor (or photoresistor) to record the spike data.

Another potential improvement could involve drilling a small hole into the cooling end’s heatsink as close to the cylinder as possible to accept a temperature probe. In this way the actual T_Low could be measured to assist with efficiency calculations; that is, how much of the temperature difference is actually converted to work. Also, it might be instructive to assess engine performance when a small fan is directed over the heatsink to assist the heat rejection. A 12V computer processor fan could used to this end. An extreme approach would be to apply ice or cold water to the heatsink for serious cooling, and to observe the effects on performance. Since ice melts at a consistent temperature, the approximate temperature of the heatsink would be known with better accuracy, would not require additional instrumentation, and would result in a greater temperature differential.

An additional way to measure work output could be to attach a small dynamo to the shaft and measure the generated DC current. This could be compared to the work computed from the P-v diagrams, and an approximate measure of the losses involved could be achieved, noting of course that the addition of the generator might introduce significant losses itself. This sort of small electricity generating Stirling engine model could be helpful in gaining insight to the efficiencies, capabilities, and utility of a Stirling engine in remote parts of the world where various fuel sources are available.

References

Haywood, David. An Introduction to Stirling Cycle Machines. http://www.mech.canterbury.ac.nz/documents/sc_intro.pdf. Linked Nov 19 2004.

Stirling Cycle Machine Analysis. http://www.ent.ohiou.edu/~urieli/stirling/me422.html. Linked Nov 19 2004

Stirling Engine History. http://web.mit.edu/2.670/www/Stirl/stirl.html. Linked Nov 19 2004.

How Stirling Engines Work. http://travel.howstuffworks.com/stirling-engine.htm. Linked Nov 19 2004.

Stirling Engine. http://www.infoplease.com/ce6/sci/A0909803.html. Linked Nov 19 2004.

Cengel, Yunus, and Michael Boles. Thermodynamics: An Engineering Approach. 4th ed. New York: McGraw Hill, 2002.

Ilak, Milos, and Hartigan, Jesse. Engineering 90 Report. 2003.

Appendix A1. p-V Plots for 28.5 Watt Input

Appendix A2. p-V Plots for 39.2 Watt Input

Appendix A3. p-V Plots for 50.2 Watt Input

Appendix B1. Matlab Scripts

analyze.m

function [ ] = analyze(file_name, description, flabel, power_in)

load(file_name);

%% analyze

%% calculate P-V work

%% all pressures in kPa

%% all volumes in m^3

pressure_volt = 1000 * ourdata(:,1); % put pressure transducer data in mV

pressures = 188.29*pressure_volt-332.22; % result is in pascals (N/m2)

pressures = pressures/1000; % put in kPa

tach_volt = ourdata(:,2);

V_dead = 3.7925e-6; % all in m^3

V_E = 2.1284e-5;

V_C = 1.7240e-5;

% tachometer data: find spike locations by looking for those points

% which dip below 90% of the maximum voltage (black tape so CdS

% photocell voltage dips when tape mark is seen)

minspikes = min(tach_volt);

maxspikes = max(tach_volt);

trigspikes = (maxspikes - minspikes)*0.8 + minspikes; % 80% above min is trigger level

spikes = find(tach_volt<trigspikes); % find spike locations

%using 9VDC supply changed output to 0.76V bright 0.91V dark

%spikes = find(tach_volt<2.3);

if isempty(spikes) == 1

spikes = 0;

else

cleandata = zeros(size(time));

cleandata(spikes) = 1; % creates unit step at spike locations

spikestart = 1; % cleans spikestart

if spikes(1) ~= 1 % if/else statement takes care of the situation where the graph begins with only part of a spike

spikestart(1) = spikes(1);

i = 2;

else

i = 1;

end

for ctr = 2:size(spikes)

if spikes(ctr - 1) ~= spikes(ctr) - 1 % find locations where spikes begin

spikestart(i) = spikes(ctr); % spikestart stores spike locations

if i >1

period = time(spikestart(i)) - time(spikestart(i-1));

pressure_array = pressures(spikestart(i-1):spikestart(i));

time_array = [0:1:(length(pressure_array)-1)]*0.001;

rpm(i) = 60/period;

phi = 2*pi/period*time_array;

volume = 1/2*V_C*(1+cos(phi-pi)) + V_dead + V_E;

%*********** now plot a PV diagram ******

%the measurement is 9.2 kPa above the real pressure due

%to suction from vacuum pump; so we subtract 9.2 kPa to get

%real values

plotdata = [volume', pressure_array];

plotdata = sortrows(plotdata, 1);

% integrate and plot the stuff

integral(i) = do_pv_integral(plotdata(:,1), plotdata(:,2), description, flabel, power_in, rpm(i), i);

figure(i)

end

i = i + 1;

end

end

end

Appendix B2. Matlab Scripts

do_pv_integral.m

function [integral] = do_pv_integral(volumes, pressures, graph_label, file_label, power_in, rpm, index)

data_points = length(volumes);

high_pressures = zeros(data_points)-1;

low_pressures = zeros(data_points)-1;

high_data_X = [ ];

high_data_Y = [ ];

low_data_X = [ ];

low_data_Y = [ ];

low_index = 1;

high_index = 1;

for M = 1 : data_points

average = 0;

count = 0;

for K = M : M+5

if K <= data_points

count = count + 1;

average = average + pressures(K);

end

end

if count > 0

average = average / count;

end

%average

if pressures(M) > average

high_data_Y(high_index) = pressures(M);

high_data_X(high_index) = volumes(M);

high_index = high_index+1;

else

low_data_X(low_index) = volumes(M);

low_data_Y(low_index) = pressures(M);

low_index = low_index+1;

end

end

area_high = trapz(high_data_X, high_data_Y);

area_low = trapz(low_data_X, low_data_Y);

integral = area_high - area_low;

plot(low_data_X, low_data_Y, high_data_X, high_data_Y, volumes, pressures);

title(sprintf('Stirling Engine P-V diagram (%s, %.1fW input, %.1f rpm) Power: %.2f W', graph_label, power_in, rpm, 1000*rpm*integral/60));

xlabel('Volume (m^3)');

ylabel('Pressure (kPa)');

filename = sprintf('graph_%s_%.1fW_%d.bmp', file_label, power_in, index);

if rpm > 30 & rpm < 400

filename

print('-dbitmap', filename);

end

Appendix B3. Matlab Scripts

pvrealtime.m

% program pvrealtime.m

% based on a program by Milos Ilak and Jesse Hartigan for their 2004 E90 project.

% This program takes input voltages from a pressure transducer (Ch.1) and tachometer

% (Ch.2) and plots real-time pressure-vs-volume diagrams. It also calculates area

% inside the p-v loop, which is the total work done by the system on its surrounding.

clc

clear all

close all

ai = analoginput('nidaq',1); % get parameters for NI-DAQ board

ai.SampleRate=1000; % sample each channel 1 kHz

ai.SamplesPerTrigger=2000; % two seconds of data;

% NOTE: If rpm is lower then 150 increase ai.SamplesPerTrigger.

% Program will not function correctly if only one spike is seen.

ai.InputType = 'SingleEnded'; % default is differential inputs

ai.TriggerType='Immediate'; % do it on start(ai)

chan = addchannel(ai,[1 2]); % note that channel names are chan(1) and chan(2)

total_time = 10; % this is how long the program will run before printing summary info

secs_per_interval = ai.SamplesPerTrigger/ai.SampleRate;

num_intervals = total_time/secs_per_interval;

%input main filename

filename = input('enter main file: ','s');

for kk = 1:num_intervals

start(ai)

[ourdata,time]=getdata(ai); % grab one data chunk from memory

%save ourdata and time to file

fileNameTemp = [filename,num2str(kk-1)];

pressure_volt = 1000 * ourdata(:,1); % put pressure transducer data in mV

pressures = 188.29*pressure_volt-332.22; % result is in pascals (N/m2)

pressures = pressures/1000; % put in kPa

tach_volt = ourdata(:,2);

ourdata(:,3) = pressures;

save(fileNameTemp,'ourdata','time');

%Define the following engine volumes, as determined from running 'stirling_power_apr10.m'

%Dead volume V_dead, expansion space volume V_E (swept volume of displacer),

%and compression space volume V_C (swept volume of power piston)

V_dead = 3.7925e-6; % all in m^3

V_E = 2.1284e-5;

V_C = 1.7240e-5;

% tachometer data: find spike locations by looking for those points

% which dip below 90% of the maximum voltage (black tape so CdS

% photocell voltage dips when tape mark is seen)

minspikes = min(tach_volt);

maxspikes = max(tach_volt);

trigspikes = (maxspikes - minspikes)*0.8 + minspikes; % 80% above min is trigger level

spikes = find(tach_volt<trigspikes); % find spike locations

%using 9VDC supply changed output to 0.76V bright 0.91V dark

%spikes = find(tach_volt<2.3);

if isempty(spikes) == 1

spikes = 0;

else

cleandata = zeros(size(time));

cleandata(spikes) = 1; % creates unit step at spike locations

spikestart = 1; % cleans spikestart

if spikes(1) ~= 1 % if/else statement takes care of the situation where the graph begins with only part of a spike

spikestart(1) = spikes(1);

i = 2;

else

i = 1;

end

for ctr = 2:size(spikes)

if spikes(ctr - 1) ~= spikes(ctr) - 1 % find locations where spikes begin

spikestart(i) = spikes(ctr); % spikestart stores spike locations

if i >1

period = time(spikestart(i)) - time(spikestart(i-1));

pressure_array = pressures(spikestart(i-1):spikestart(i));

time_array = [0:1:(length(pressure_array)-1)]*0.001;

rpm(i) = 60/period;

phi = 2*pi/period*time_array;

volume = 1/2*V_C*(1+cos(phi-pi)) + V_dead + V_E;

%*********** now plot a PV diagram ******

%the measurement is 9.2 kPa above the real pressure due

%to suction from vacuum pump; so we subtract 9.2 kPa to get

%real values

plot(volume,pressure_array)

xlabel('volume (m^3)')

ylabel('pressure (kPa)')

title('P-V Diagram')

pause(0.1)

end

i = i + 1;

end

end

current_mean_rpm = mean(rpm(2:length(rpm)))

end

end

delete(ai)