Seat Vibration Controller Analysis

This case study analyzes the amount of vibration a passenger experiences for a vehicle traveling over a road disturbance (bump). We want to determine the amount of reduction in displacement and acceleration that can be achieved using two different controllers. We have test data across 12 vehicle speeds for three schemes:

Baseline Suspension System (uncontrolled)
Controller 1: Actively Controlled Seat + Suspension System
Controller 2: Low Cost Version of Controller 1

The primary difference between the two controller configurations are the implementation of the linear motor used to dampen oscillations in the seat. The dynamic range of the linear motors are:

Controller 1: +/- 100 N
Controller 2: +/- 10 N

Load data files
Filter data
Create visualization
Calculate Maximum Values
Visualize Maximum Values
Model Controller Acceleration
Load Survival Data and Visualize
ANOVA on Survival Data
Plot Survival Rate
Summary

Load data files

First, all of the Excel data files are imported using a for-loop. The data is stored in a dataset array (Stats Toolbox).

% Create a list of data files to import
filenames = dir('data\seatvib*.xls');
filenames = sort({filenames.name})';

nFiles = length(filenames);  % Number of data files
data   = cell(nFiles, 12); % Preallocate cell array for data storage

h = waitbar(0, 'Loading data files...'); % Progress indicator

% Loop through data files
for i = 1:nFiles

    name = filenames{i};

    % Update progress indicator
    waitbar(i/nFiles, h, ['Loading data files...' name]);

    % Read from Excel file, then arrange data into
    % a more convenient format
    x         = xlsread(name);
    data(i,:) = fcnStoreData(name, x);

end

% Convert cell array data to dataset object
ds = fcnConvert2Dataset(data);

close(h);  % Close waitbar window

Filter data

The data is a bit noisy, so we will filter out the high frequency component. We will use FDATOOL to design a digital filter, and the auto-generated M-file is used.

% Sampling interval
dt = diff(ds{1,'Time'});

% Sampling frequency
Fs = 1 / mean(dt);
% Cutoff frequency, [Hz]
Fc = 10;

% Create filter using auto-generated m-file.
Hd = myFilter(Fc, Fs);

h = waitbar(0, 'Applying filter...');  % Progress bar

% Loop through data
for i = 1:nFiles

    % Update progress indicator
    msg = ['Applying filter ' num2str(i) ' of ' num2str(nFiles)];
    waitbar(i/nFiles, h, msg);

    % Filter displacement data
    ds{i,'DispFilt'}  = filter(Hd, ds{i,'Displacement'});
    % Filter acceleration data
    ds{i,'AccelFilt'} = filter(Hd, ds{i,'Acceleration'});
end

close(h);

Create visualization

We will visualize the displacement and acceleration data (raw and filtered) using the auto-generated function. This function can be easily created using Plot Tool.

% Select a file to visualize.
obsname = 'seatvib-control2-speed040.xls';
fcnCreateFigure(ds{obsname, 'Time'        }, ...
                ds{obsname, 'Acceleration'}, ...
                ds{obsname, 'AccelFilt'   }, ...
                ds{obsname, 'Displacement'}, ...
                ds{obsname, 'DispFilt'    }, ...
                ds{obsname, 'CForce'      }, ...
                obsname);

Calculate Maximum Values

To evaluate the different configurations, we'll look at the maximum displacement, acceleration, and controller force (as % Full Range) observed at each speed.

h = waitbar(0, 'Calculating max values...');  % Progress indicator

% Loop through data
for i = 1:nFiles

    % Update progress indicator
    msg = ['Calculating max values ' num2str(i) ' of ' num2str(nFiles)];
    waitbar(i/nFiles, h, msg);

    % Determine max displacement
    ds{i, 'DispMax'}  = max(abs(ds{i, 'DispFilt'}));
    % Determine max acceleration
    ds{i, 'AccelMax'} = max(abs(ds{i, 'AccelFilt'}));
    % Determine max control force (percentage of max range)
    ds{i, 'ForceMax'} = max(abs(ds{i, 'CForce'}))/ds{i, 'MaxForceRange'}*100;
end

close(h);

Visualize Maximum Values

Plot the maximum values

dTypes     = {'DispMax'     , 'AccelMax'       , 'ForceMax'         };
labStrings = {'Max Disp (m)', 'Max Acc (m/s^2)', 'Max Force (% Max)'};
yScales    = {'log'         , 'log'            , 'linear'           };
colors     = {'red'         , 'green'          , 'black'            };
location   = {'northwest'   , 'center'         , 'southeast'        };

for i = 1:length(dTypes)

      figure('WindowStyle','docked');

      tt = unique(ds.TestType);
      for j = 1:length(tt)
            id = ds.TestType == tt(j);
            x  = ds.Speed(id);
            y  = ds.(dTypes{i})(id);
            hLine(j) = line(x, y, 'Marker', 'o',       ...
                                  'Color',  colors{j}); %#ok
      end

      % Format the axes and add markup
      box on;
      set(gca,'YScale', yScales{i});
      xlabel('Speed (mph)');
      ylabel(labStrings{i});
      legend(hLine, char(tt),                      ...
                    'Location',    'SouthOutside', ...
                    'Orientation', 'Horizontal');
end

Looking at displacement and acceleration, we see that Controller 1 and 2 are different primarily at low speeds and are similar at high speeds. We need to determine if this difference is significant and not a result of measurement error or noise.

Model Controller Acceleration

We want to determine if the difference in the two controllers are significant in terms of the acceleration observed. The acceleration curves for both controllers appear to be similar and are significantly different at the low speed range. Let's determine if a single model captures both trends with only a change in one parameter. The model we will fit is:

$y = ax_1^{b+cx_2}+d$

where x1 is speed and x2 is a binary variable representing the different configurations (0 for Controller 1, 1 for Controller 2).

% Functional representation of above equation
f = @(a,x) a(1)*x(:,1).^(a(2)+a(3)*x(:,2))+a(4);

fcnModelController(ds,f);

    'Coef'    'Estimate'    'Lower'        'Upper'   
    'a'       [  2.8174]    [   2.1349]    [     3.5]
    'b'       [ -3.5529]    [  -3.9698]    [  -3.136]
    'c'       [  1.5559]    [   1.2894]    [  1.8223]
    'd'       [0.017304]    [0.0091098]    [0.025497]

Note that the confidence interval for 'c' does not include 0, so the difference between Controllers is significant and not a result of noise or measurement error. At this point, Controller 1 has the best overall performance. However, if we are only concerned with the performance at moderate to high speed (>10 m/s) either configuration will be acceptable. Controller 2 has considerable cost advantage and we don't want to rule it out at this point.

Review of the controller force shows that Controller 2 is saturated (max force) over the low speeds. This is a concern from a reliability point of view. With a saturated controller, we may run the risk of early failure and need to quantify if this will have an impact on our selection of controllers. For our intended application, we need to have >90% survive at least 1,000 hours.

Load Survival Data and Visualize

To gain insight into the reliability of our two controller configurations, we tested 20 of each until failure. The test consisted of the entire seat and controller assembly and cycled under worst case conditions. The data was saved in an Excel Spreadsheet as the number of hours the assembly was tested before failure.

lifeData = dataset('XLSFile', 'data\lifeData.xls');

%avgLife = grpstats(lifeData)     %#ok
avgLife = mean(double(lifeData))

figure;
plot(lifeData,'o-')
xlabel('Observation');
ylabel('Survival Time (hours)');
legend('Controller 1', 'Controller 2');

avgLife =
       1951.7       1716.7

The average values of the two configurations are similar in magnitude. In order to determine if the difference is real, and not a result of our sample size, we'll perform an Anslysis of Variance test to determine if the average values are statistically different.

ANOVA on Survival Data

ANOVA is an analysis technique that allows us to test for differences in means across groups of data. We'll use ANOVA to test if the difference in means is actually different given our small sample size.

%p = kruskalwallis(double(lifeData));
p = anova1(double(lifeData));
xlabel('Controller');
ylabel('Survival Time (hours)');

Plot Survival Rate

Survival data commonly follows a Weibull distribution. We'll check that our data follows a Weibull distribution and use the fitted distribution to estimate survival rates for our two configurations. The Statistics Toolbox has a graphical tool (dfittool) that provides a central access point for fitting distribution to data and performing related calculations from the fitted distributions. We'll use this tool to perform our Weibull Analysis on Controller 1, generate M-files to capture our analysis, and apply the analysis to Controller 2.

% Survivor Function for Controller 1
figure;
fcnPlotSurvivalRate(lifeData.control1, gca);
title('Controller 1');
movegui('west');

% Survivor Function for Controller 2
figure;
fcnPlotSurvivalRate(lifeData.control2, gca);
title('Controller 2')
movegui('east');

The Survivor Function Plots show the fraction of the population that are likely to fail over time. Controller 1 has >95% survival rate at 1,000 hours. This meets our goal. Controller 2, on the other hand, has a survival rate of around 80% at 1,000 hours. This does not meet our goal. Therefore, Controller 1 is the selected configuration for our application.

Summary

We have seen that the two controller configurations testted significantly reduce the amount of displacement and acceleration a passenger is exposed to. Controller 1 meets our reliability criteria as well as had the best overall performance. Controller 2, while cheaper, does not meet our reliability requirement of >90% for >1,000 hours. Thus, Controller 1 is the configuration selected for further development.

close('all');