Solution 8.8.3.mtr3x The full MATLAB program used to nd the transfer function is: % %load the data % load mtr3step t=mtr3step(1:500,1);; y=mtr3step(1:500,2);; plot(t,y) print -deps sr883mtr3xa.eps y=y(51:500);; t=t(51:500);; y=y-y(1);; y(1) = 0;; t(1) = 0;; plot(t,y) print -deps sr883mtr3xb.eps y(180:450) = mean(y(180:450));; y(175:179) = mean(y(175:179));; y(170:174) = mean(y(170:174));; y(165:169) = mean(y(165:169));; y(160:164) = mean(y(160:164));; y(155:159) = mean(y(155:159));; y(150:154) = mean(y(150:154));; y(145:149) = mean(y(145:149));; y(140:144) = mean(y(140:144));; y(135:139) = mean(y(135:139));; y(130:134) = mean(y(130:134));; y(125:129) = mean(y(125:129));; y(120:124) = mean(y(120:124));; y(115:119) = mean(y(115:119));; y(110:114) = mean(y(110:114));; y(105:109) = mean(y(105:109));; y(105:109) = mean(y(105:109));; y(100:104) = mean(y(100:104));; y(95:99) = mean(y(95:99));; y(90:94) = mean(y(90:94));; y(85:89) = mean(y(85:89));; y(80:84) = mean(y(80:84));; y(75:79) = mean(y(75:79));; y(70:74) = mean(y(70:74));; y(65:69) = mean(y(65:69));; y(60:64) = mean(y(60:64));; 1 plot(t,y) print -deps sr883mtr3xc.eps % zero the summer % asum = 0.0;; % %Find the length of the array y(t) Note:this the step response. It has to be % pulled into Matlab before this program is run. % top = size(y);; top = top(1,1) % % Initialize counter. % k=1;; % % % Find maximum value of step response. % K0 = y(top) % % % Form function K0 - y(t). Note y(t) is the same as yu(t) in the write up. % ftemp = K0 - y;; % fsave1 = ftemp;; % % Reset counter % k=2;; % % Integrate K0 - y(t). and form the function y1(t). % while ( k <= top) l= k-1;; asum = asum + ( ( ( ftemp(k) + ftemp(l) ) /2 ) *( t(k) - t(l) ) );; t1(k) = t(k);; f1(k) = asum;; k=k+1;; end % % Save data for plotting. % 2 save time.dat t1 -ascii save f1.dat f1 -ascii save ftempa.dat ftemp -ascii % % Set K1 = to max(y1(t)), that is its final steady state value. % K1 = asum % % % form function: K1 - y1(t). % ftemp1 = K1 - f1;; % % Save function for later plotting % save ftempb.data ftemp1 -ascii % % reset summer % asum = 0.0;; % % Reset counter % k=2;; % % Integrate K1 - y1(t) to get y2(t) (called f2 in program) % while ( k <= top) l= k-1;; asum = asum + ( ( ( ftemp1(k) + ftemp1(l) ) /2 ) *( t(k) - t(l) ) );; f2(k) = asum;; k=k+1;; end % % Set K2 equal to maximum value of y2(t). % K2 = asum % % Save function y2(t) for later plotting. % save f2.data f2 -ascii % % Find coefficients and then poles of transfer function. % a1 = K1/K0 3 a2 = ( a1*K1 - K2 )/K0 cpol = [a2 a1 1] p =roots(cpol) % % Find gain of transfer function % Kplant = K0/a2 g=zpk([],[p(1) p(2)],Kplant) T=linspace(0,0.5,450);; u=ones(450,1);; [yhat,T] = lsim(g,u,T);; plot(t,y,'k-',T,yhat,'k--') print -deps sr883mtr3xd.eps load mtr3av topmav = size(mtr3av) topav = topm(1,1) tav = mtr3av(1:topav,1);; av = mtr3av(1:topav,2);; plot(tav,av) print -deps 883mtr3xav.eps p1 = 28;; p2 = 60;; yss = y(top);; KK = yss * p1*p2;; g=zpk([],[-p1 -p2],KK) [yhat1,T] = lsim(g,u,T);; plot(t,y,'k-',T,yhat1,'k--') The rst step is to get the data into MATLAB. That is done by the MAT- LAB statements: load mtr1step t=mtr3step(1:500,1);; y=mtr3step(1:500,2);; plot(t,y) print -deps sr883mtr3xa.eps The plot of the raw date is shown in Figure1. As can be seen, the plot does not start at t =0and the initial voltage is not zero. The MATLAB statements y=y(51:500);; t=t(51:500);; %y = y-y(1);; 4 -0.05 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0 0.5 1 1.5 2 2.5 3 3.5 4 Figure 1: Rawstep response data y(1) = 0;; t(1) = 0;; plot(t,y) print -deps sr883mtr3xb.eps produce the adjusted step respone shown in Figure 2 The rest of the MAT- LAB program is designed to do the successiveintegrations necessary to nd the twopoles. The MATLAB statements y(180:450) = mean(y(180:450));; y(175:179) = mean(y(175:179));; y(170:174) = mean(y(170:174));; y(165:169) = mean(y(165:169));; y(160:164) = mean(y(160:164));; y(155:159) = mean(y(155:159));; y(150:154) = mean(y(150:154));; y(145:149) = mean(y(145:149));; y(140:144) = mean(y(140:144));; y(135:139) = mean(y(135:139));; y(130:134) = mean(y(130:134));; 5 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0 0.5 1 1.5 2 2.5 3 3.5 4 Figure 2: Adjusted step response data y(125:129) = mean(y(125:129));; y(120:124) = mean(y(120:124));; y(115:119) = mean(y(115:119));; y(110:114) = mean(y(110:114));; y(105:109) = mean(y(105:109));; y(105:109) = mean(y(105:109));; y(100:104) = mean(y(100:104));; y(95:99) = mean(y(95:99));; y(90:94) = mean(y(90:94));; y(85:89) = mean(y(85:89));; y(80:84) = mean(y(80:84));; y(75:79) = mean(y(75:79));; y(70:74) = mean(y(70:74));; y(65:69) = mean(y(65:69));; y(60:64) = mean(y(60:64));; are a crude smoothing program designed to help the algorith converge. The smoothed response is shown in Figure 3 6 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0 0.5 1 1.5 2 2.5 3 3.5 4 Figure 3: Smoth step response data For this problem the alogrithm gives the answer G(s)= 4967:3 (s +36:576;7:8879)(s+36:576;7:8879) : This is unadjusted for the size of the input voltage, about 32 V, a shown inf Figure 4. So we nally arriveat G(s)= 152:23 (s +36:576;7:8879)(s+36:576;7:8879) : The step response of the model is compared to the actual recorded step response in Figure 5 The matchisprettygood. However, weknowfromthe physics of the motor that weshould havetwo real poles. The poles aboveareessentially critically damped. So what wehavefound is the approximate dominant time constantofthe motor. After a little searching, we arriveat G(s)= 186:27 (s +28)(s +60) ;; which produces the step response shown in Figure 6. 7 -0.05 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0 5 10 15 20 25 30 35 Figure 4: Armature voltage applied to motor 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 0.5 1 1.5 2 2.5 3 3.5 4 Figure 5: Comparison of model step response to actual step response 8 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 0.5 1 1.5 2 2.5 3 3.5 4 Figure 6: Comparison of model step response to actual step response 9