Toko. 2020-10-21 15:16 采纳率: 14.3%
浏览 158

MATLAB索引超出数组范围

以下是程序:

clear;
clc;
E=20
L=3.3e-3
C=1e-3
R=19;
T=4e-4;
a = 1/(2*R*C);
w = sqrt(1/L/C-1/4/R^2/C^2);
A2 = [0,-1/L;1/C,-1/R*C];
A1 = A2;
A3=[0,0;0,-1/C*R];
B1=[1/L;0];
B2=[0;0];
B3=B2;
I=[1,0;0,1];
for m=1:5001
    Iref=0.1+0.01*(m-1);
    i(1)=0;
    u(1)=0;

    for n=1:5000

       X = [ i( n ); u( n ) ];   
        p1=(exp(-a*T)/w*((a*sin(w*T))+(w*cos(w*T))));
        q1=exp(-a*T)/w*sin(w*T);
        Ib1=(Iref/p1)+(q1*u(n)/L/p1)-(((1-p1)/R)+(q1/L))*E/q1;
        if i(n)<=Ib1
            X(:,n+1)=exp(A1*T)*X+(exp(A1*T)-I)*A1^-1*B1*E;
        else
         syms t1;
         p3=exp(A1*t1)*X+(exp(A1*t1)-I)*A1^-1*B1*E;
         q3=p3(1,:);
            t1=solve(Iref==q3,t1);
           t1=double(t1);  
           p2=(exp(-a*t1)/w*(a*sin(w*t1)+w*cos(w*t1)));
           q2=exp(-a*t1)/w*sin(w*t1);
            Ib2 = Iref/p2+q2*u(n)/L/p2-((1-p2)/R+q2/L)*E/p2;
            if i(n) < Ib2
                t2=T-t1;
                X(:,n+1)=exp(A2*t2)*exp(A1*t1)*X+exp(A2*t2).*(exp(A1*t1)-I)*(A1^-1)*B1*E;


            else
                syms t2 ut1;
                S=[Iref;ut1];
               p4=exp(A2*t2)*S;
               q4=p4(1,:);
               t2=solve(0==q4,t2);
               t2=double(t2);
                t3=T-t1-t2;
                X(:,n+1)=exp(A3*t3)*exp(A2*t2)*exp(A1*t1)*X+exp(A3*t3)*exp(A2*t2)*(exp(A1*t1)-I)*A1^-1*B1*E;
            end
        end
    end
    i1(m,:) = i(end-150:end);

%   uo(m,:) = u(end-80:end);
    H= plot(Iref*ones(1,length(i1(m,:))),i1(m,:));
%  set(H,'linestyle','none','marker','.','markersize',3);
    hold on
    drawnow
end

程序是从X(:,n+1)=exp(A3*t3)*exp(A2*t2)*exp(A1*t1)*X+exp(A3*t3)*exp(A2*t2)*(exp(A1*t1)-I)*A1^-1*B1*E;循环一次后开始第二次循环时在 X = [ i( n ); u( n ) ]; 发生错误,我认为是初始化不能放在n循环里,但是X又跟n循环有关,又不能直接放在外面,恳请各位大佬们教教新手,如何把初始化那句 X = [ i( n ); u( n ) ]; 放在n循环外面?或者让这句话只在n循环中执行一次?谢谢各位!!!!!!

  • 写回答

1条回答 默认 最新

  • 你知我知皆知 2024-07-25 01:38
    关注

    以下回答参考 皆我百晓生券券喵儿 等免费微信小程序相关内容作答,并由本人整理回复。

    你的问题主要在于 X 的初始化。在 MATLAB 中,当你使用 X = [ i( n ); u( n ) ]; 这样的语法时,你实际上是在定义了一个新的矩阵 X,而不是对当前的 iu 列进行修改。

    如果你想要将这个初始化操作移到 n 循环之外,你可以这样写:

    clear;
    clc;
    E=20
    L=3.3e-3
    C=1e-3
    R=19;
    T=4e-4;
    a = 1/(2*R*C);
    w = sqrt(1/L/C-1/4/R^2/C^2);
    A2 = [0,-1/L;1/C,-1/R*C];
    A1 = A2;
    A3=[0,0;0,-1/C*R];
    B1=[1/L;0];
    B2=[0;0];
    B3=B2;
    I=[1,0;0,1];
    
    % 初始化步骤
    i = zeros(5001,1);
    u = zeros(5001,1);
    
    for m=1:5001
        i(1)=0;
        u(1)=0;
    
        for n=1:5000
    
           X = [ i( n ); u( n ) ];   
            p1=(exp(-a*T)/w*((a*sin(w*T))+(w*cos(w*T))));
            q1=exp(-a*T)/w*sin(w*T);
            Ib1=(Iref/p1)+(q1*u(n)/L/p1)-(((1-p1)/R)+(q1/L))*E/q1;
            if i(n)<=Ib1
                X(:,n+1)=exp(A1*T)*X+(exp(A1*T)-I)*A1^-1*B1*E;
            else
             syms t1;
             p3=exp(A1*t1)*X+(exp(A1*t1)-I)*A1^-1*B1*E;
             q3=p3(1,:);
                t1=solve(Iref==q3,t1);
               t1=double(t1);  
               p2=(exp(-a*t1)/w*(a*sin(w*t1)+w*cos(w*t1)));
               q2=exp(-a*t1)/w*sin(w*t1);
                Ib2 = Iref/p2+q2*u(n)/L/p2-((1-p2)/R+q2/L)*E/p2;
                if i(n) < Ib2
                    t2=T-t1;
                    X(:,n+1)=exp(A2*t2)*exp(A1*t1)*X+exp(A2*t2).*(exp(A1*t1)-I)*(A1^-1)*B1*E;
    
    
                else
                    syms t2 ut1;
                    S=[Iref;ut1];
                   p4=exp(A2*t2)*S;
                   q4=p4(1,:);
                   t2=solve(0==q4,t2);
                   t2=double(t2);
                    t3=T-t1-t2;
                    X(:,n+1)=exp(A3*t3)*exp(A2*t2)*exp(A1*t1)*X+exp(A3*t3)*exp(A2*t2)*(exp(A1*t1)-I)*A1^-1*B1*E;
                end
            end
        end
        i1(m,:) = i(end-150:end);
    
    %   uo(m,:) = u(end-80:end);
        H= plot(Iref*ones(1,length(i1(m,:))),i1(m,:));
    %  set(H,'linestyle','none','marker','.','markersize',3);
        hold on
        drawnow
    end
    

    这段代码首先定义了两个向量 iu,然后在每次循环前都初始化它们。这样,即使在计算过程中改变了 iu,也不会影响到 X 的初始值。

    评论

报告相同问题?

悬赏问题

  • ¥60 如何批量获取json的url
  • ¥15 对法兰连接元件所承受的表面载荷等效转化为法兰开孔接触面上的等效表面载荷?
  • ¥15 comsol仿真压阻传感器
  • ¥15 Python线性规划函数optimize.linprog求解为整数
  • ¥15 llama3中文版微调
  • ¥15 在win系统Django虚拟环境下载mysqlclient报错
  • ¥15 pg数据库导入数据序列重复
  • ¥15 三分类机器学习模型可视化分析
  • ¥15 本地测试网站127.0.0.1 已拒绝连接,如何解决?(标签-ubuntu)
  • ¥50 Qt在release捕获异常并跟踪堆栈(有Demo,跑一下环境再回答)