Changeset 905

Show
Ignore:
Timestamp:
04/14/10 11:24:04 (14 years ago)
Author:
vahalam
Message:
 
Files:
1 modified

Legend:

Unmodified
Added
Removed
  • applications/dual/IterativeLocal/pmsm_lqg.m

    r898 r905  
    3333e = DELTAt*ckp*cp*cp*cPSIpm/cJ; 
    3434 
    35 OMEGAt = 1.0015; 
     35OMEGAt = 1.15;%1.0015; 
    3636 
    3737%penalizace vstupu a rizeni 
     
    6868R = diag([0.0006, 0.0006]); 
    6969 
    70 x0 = [0 0 1-OMEGAt pi/2 OMEGAt]; 
     70x0 = [0 0 1.0-OMEGAt pi/2 OMEGAt]; 
    7171P = diag([0.01, 0.01, 0.01, 0.01, 0]); 
    7272 
     
    8181%zapinani a vypinani sumu, sumu v simulaci a generovani trajektorii s 
    8282%rozptylem 
    83 sum = 1;%0.01; 
     83% sum = 1;%0.01; 
    8484sumsim = 1;%0.01; 
    8585neznalost = 1; 
    8686 
    87  
     87errnans = 0; 
    8888 
    8989% vycisti kreslici okno 
     
    141141            xn(:, 1, n) = x00; 
    142142            for k = 1:Kt+K-1,                
    143                u(:, k) = L(:, :, k)*(xn(:, k, n));                
     143               u(:, k) = L(:, :, k)*(xn(:, k, n)); %tady se vyuzije(k) / nevyuzije(1) receding horizon               
    144144               xn(1, k+1, n) = a*xn(1, k, n) + b*(xn(3, k, n) + xn(5, k, n))*sin(xn(4, k, n)) + c*u(1, k) + sumsim*sqrt(Q(1, 1))*randn();   
    145145               xn(2, k+1, n) = a*xn(2, k, n) - b*(xn(3, k, n) + xn(5, k, n))*cos(xn(4, k, n)) + c*u(2, k) + sumsim*sqrt(Q(2, 2))*randn(); 
     
    148148               xn(5, k+1, n) = xn(5, k, n); 
    149149            end 
    150  
    151  
     150             
     151    %kontrola spatne L matice 
     152    lstore(:,:,n) = L(:,:,1); 
     153    if(isnan(sum(sum(L(:,:,1))))==1) 
     154        errnans = errnans + 1; 
     155    end 
     156     
    152157    %vykresleni 
    153158        subplot(2, 3, 1); 
     
    179184toc 
    180185 
    181 % %hlavni iteracni smycka 
    182 % for iterace = 1:It, 
    183 %     %generovani stavu 
    184 %     for n = 1:N, 
    185 %         xn(:, 1, n) = x0' + neznalost*sqrtm(P)*randn(4,1); 
    186 %         for k = 1:Kt-1, 
    187 %            tu =  L(:, :, k)*(xn(:, k, n) - [0 0 OMEGAt 0]'); 
    188 %            xn(1, k+1, n) = a*xn(1, k, n) + b*xn(3, k, n)*sin(xn(4, k, n)) + c*tu(1) + sumsim*sqrt(Q(1, 1))*randn();   
    189 %            xn(2, k+1, n) = a*xn(2, k, n) - b*xn(3, k, n)*cos(xn(4, k, n)) + c*tu(2) + sumsim*sqrt(Q(2, 2))*randn(); 
    190 %            xn(3, k+1, n) = d*xn(3, k, n) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sumsim*sqrt(Q(3, 3))*randn(); 
    191 %            xn(4, k+1, n) = xn(4, k, n) + xn(3, k, n)*DELTAt + sumsim*sqrt(Q(4, 4))*randn();            
    192 % %            xn(1, k+1, n) = a*xn(1, k, n) + b*(xn(3, k, n)+OMEGAt)*sin(xn(4, k, n)) + c*tu(1) + sum*sqrt(Q(1, 1))*randn();   
    193 % %            xn(2, k+1, n) = a*xn(2, k, n) - b*(xn(3, k, n)+OMEGAt)*cos(xn(4, k, n)) + c*tu(2) + sum*sqrt(Q(2, 2))*randn(); 
    194 % %            xn(3, k+1, n) = -OMEGAt + d*(xn(3, k, n)+OMEGAt) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sum*sqrt(Q(3, 3))*randn(); 
    195 % %            xn(4, k+1, n) = xn(4, k, n) + (xn(3, k, n)+OMEGAt)*DELTAt + sum*sqrt(Q(4, 4))*randn(); 
    196 %         end 
    197 %     end 
    198 %      
    199 %     %prumerny stav 
    200 %     xs = mean(xn, 3); 
    201 %      
    202 %     %napocteni S a L 
    203 %     for kt = 1:Kt-1, 
    204 % %     receding horizon 
    205 %      if((K-1+kt)<Kt) 
    206 %         S(:, :, K) = X; 
    207 %         for k = K-1+kt-1:-1:1+kt-1, 
    208 %             A(3, 1) = -e*sin(xs(4, k)); 
    209 %             A(3, 2) = e*cos(xs(4, k)); 
    210 %             A(1, 3) = b*sin(xs(4, k)); 
    211 %             A(2, 3) = -b*cos(xs(4, k)); 
    212 %             A(1, 4) = b*(xs(3, k))*cos(xs(4, k)); 
    213 %             A(2, 4) = b*(xs(3, k))*sin(xs(4, k)); 
    214 % %             A(1, 4) = b*(xs(3, k)+OMEGAt)*cos(xs(4, k)); 
    215 % %             A(2, 4) = b*(xs(3, k)+OMEGAt)*sin(xs(4, k)); 
    216 %             A(3, 4) = -e*(xs(2, k)*sin(xs(4, k) + xs(1,k)*cos(xs(4, k)))); 
    217 %             S(:, :, k-kt+1) = A'*(S(:, :, k-kt+2) - S(:, :, k-kt+2)*B*inv(B'*S(:, :, k-kt+2)*B + Y)*B'*S(:, :, k-kt+2))*A + X;         
    218 %         end   
    219 %         L(:, :, kt) = -inv(B'*S(:, :, 1)*B + Y)*B'*S(:, :, 1)*A; 
    220 %      else 
    221 %         L(:, :, kt) = L(:, :, kt-1); %kopiruje poslednich K kroku z Kt kde to nejde na K predpocitat 
    222 %      end 
    223 %          
    224 %     end 
    225 % end 
    226 %     toc 
    227 %  
    228 % %vysledky 
    229 %     clf 
    230 %     subplot(2, 3, 3); 
    231 %     hold all 
    232 %     plot(1:Kt, OMEGAt*ones(1,Kt)); 
    233 % % L(:,:,1) 
    234 %     for n = 1:N, 
    235 %         xn(:, 1, n) = x0' + neznalost*sqrtm(P)*randn(4,1); 
    236 %         for k = 1:Kt-1, 
    237 %            tu =  L(:, :, 1)*(xn(:, k, n) - [0 0 OMEGAt 0]'); 
    238 %            xn(1, k+1, n) = a*xn(1, k, n) + b*xn(3, k, n)*sin(xn(4, k, n)) + c*tu(1) + sumsim*sqrt(Q(1, 1))*randn();   
    239 %            xn(2, k+1, n) = a*xn(2, k, n) - b*xn(3, k, n)*cos(xn(4, k, n)) + c*tu(2) + sumsim*sqrt(Q(2, 2))*randn(); 
    240 %            xn(3, k+1, n) = d*xn(3, k, n) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sumsim*sqrt(Q(3, 3))*randn(); 
    241 %            xn(4, k+1, n) = xn(4, k, n) + xn(3, k, n)*DELTAt + sumsim*sqrt(Q(4, 4))*randn(); 
    242 % %            xn(1, k+1, n) = a*xn(1, k, n) + b*(xn(3, k, n)+OMEGAt)*sin(xn(4, k, n)) + c*tu(1) + sum*sqrt(Q(1, 1))*randn();   
    243 % %            xn(2, k+1, n) = a*xn(2, k, n) - b*(xn(3, k, n)+OMEGAt)*cos(xn(4, k, n)) + c*tu(2) + sum*sqrt(Q(2, 2))*randn(); 
    244 % %            xn(3, k+1, n) = -OMEGAt + d*(xn(3, k, n)+OMEGAt) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sum*sqrt(Q(3, 3))*randn(); 
    245 % %            xn(4, k+1, n) = xn(4, k, n) + (xn(3, k, n)+OMEGAt)*DELTAt + sum*sqrt(Q(4, 4))*randn(); 
    246 %             
    247 %            u(:, k) = tu;           
    248 %         end 
    249 %          
    250 % %         xn(3, :, n) = xn(3, :, n) + OMEGAt*ones(1, Kt); 
    251 %          
    252 %         subplot(2, 3, 1); 
    253 %         hold all 
    254 %         plot(1:Kt, xn(1, :, n)); 
    255 %         title('i_{\alpha}'); 
    256 %         subplot(2, 3, 2); 
    257 %         hold all 
    258 %         plot(1:Kt, xn(2, :, n)); 
    259 %         title('i_{\beta}'); 
    260 %         subplot(2, 3, 3); 
    261 %         hold all 
    262 %         plot(1:Kt, xn(3, :, n)); 
    263 %         title('\omega'); 
    264 %         subplot(2, 3, 4); 
    265 %         hold all 
    266 %         plot(1:Kt, xn(4, :, n)); 
    267 %         title('\theta'); 
    268 %         subplot(2, 3, 5); 
    269 %         hold all 
    270 %         plot(1:Kt, u(1, :)); 
    271 %         title('u_{\alpha}'); 
    272 %         subplot(2, 3, 6); 
    273 %         hold all 
    274 %         plot(1:Kt, u(2, :)); 
    275 %         title('u_{\beta}'); 
    276 %     end 
     186if(errnans > 0) 
     187        lstore 
     188        disp('L is NaN') 
     189        errnans 
     190end 
     191 
     192 
    277193 
    278194