function [sys,x0,str,ts] = sfunmoment(t, x, u, flag,...
parameters),
% An S-function which computes the Mirror Moment
%
% The input is in the Form:
% u = [V_M(t), phi(t), V_L(t), V_R(t)]
%
% The parameter structure has the following entries:
% d0: gap width
% lm: mirror lenght
% bm: mirror width
%
% Heiko Wolfram 2005.
%
% $Id: sfunmoment.m,v 0.1 2005/08/12 12:57:04 hwol Exp $
% What is returned by SFUNC at a given point in time, T, depends on the
% value of the FLAG, the current state vector, X, and the current
% input vector, U.
%
% FLAG RESULT DESCRIPTION
% ----- ------ --------------------------------------------
% 0 [SIZES,X0,STR,TS] Initialization, return system sizes in SYS,
% initial state in X0, state ordering strings
% in STR, and sample times in TS.
% 1 DX Return continuous state derivatives in SYS.
% 2 DS Update discrete states SYS = X(n+1)
% 3 Y Return outputs in SYS.
% 4 TNEXT Return next time hit for variable step sample
% time in SYS.
% 5 Reserved for future (root finding).
% 9 [] Termination, perform any cleanup SYS=[].
switch flag,
%%%%%%%%%%%%%%%%%%
% Initialization %
%%%%%%%%%%%%%%%%%%
case 0,
[sys,x0,str,ts]=mdlInitializeSizes;
%%%%%%%%%%%
% Outputs %
%%%%%%%%%%%
case 3,
sys=mdlOutputs(t, x, u, parameters);
%%%%%%%%%%%%%%%%
% Unused flags %
%%%%%%%%%%%%%%%%
case { 1, 2, 4, 9 },
sys = [];
case 'moment',
sys = LocalMoment(parameters, u, x);
%%%%%%%%%%%%%%%%%%%%
% Unexpected flags %
%%%%%%%%%%%%%%%%%%%%
otherwise
error(['Unhandled flag = ',num2str(flag)]);
end % case flag
% end sfun_moment
%
%=============================================================================
% mdlInitializeSizes
% Return the sizes, initial conditions, and sample times for the S-function.
%=============================================================================
%
function [sys,x0,str,ts] = mdlInitializeSizes
% simsizes for a sizes structure
sizes = simsizes;
sizes.NumContStates = 0; % number of continuous states
sizes.NumDiscStates = 0; % number of discrete states
sizes.NumOutputs = 1; % number of system outputs (length of output y)
sizes.NumInputs = 4; % number of system inputs (length of input u)
sizes.DirFeedthrough = 1; % direct feedthrough flag; used to detect algebraic loops.
sizes.NumSampleTimes = 0; % at least one sample time is needed
sys = simsizes(sizes);
% initialize the initial conditions
x0 = [];
% str is always an empty matrix
str = [];
% initialize the array of sample times
ts = [];
% end mdlInitializeSizes
%
%=============================================================================
% mdlOutputs
% Return the block outputs.
%=============================================================================
%
function sys = mdlOutputs(t,x,u,params)
sys = LocalMoment(params, u(1) + u(4), u(2)) -...
LocalMoment(params, u(1) + u(3), -u(2));
% end mdlOutputs
%
%=============================================================================
% LocalMoment
% Calculate Moment for one Plate.
%=============================================================================
%
function yout = LocalMoment(params, u, phi)
epsilon0 = 8.8541978E-12;
if abs(phi) < sqrt(eps),
% Taylor Approximation
yout = (params.bm^2*params.lm*epsilon0*...
(-3*params.d0 + 2*params.bm*phi)*u^2)/(48*params.d0^3);
else;
yout = params.lm*epsilon0*u^2*cos(phi)/sin(phi)^2*...
( (params.bm*sin(phi))/(4*params.d0 + 2*params.bm*sin(phi)) - ...
atanh((params.bm*sin(phi))/(4*params.d0 + params.bm*sin(phi))) ...
);
end;
% end LocalMoment
%% End of File