S函数运行报错

偶编写了个S函数，用于求解偏微分方程，运行时出现报错：
Output returned by S-function 'cyjzibian' in 'untitled/S-Function' during flag=3 call must be a real vector of length 1.
不知道错误处在哪里，亟盼高手指教，谢谢哈

你应该把你的代码贴出来
要不然叫别人怎么帮你。。

你的初始化有没有问题．

我就是找不出问题在哪里，以前好像出过一次这样的情况，后来不记得怎么弄好了
我把源代码粘过来，亟盼各位大侠指点
function [sys,x0,str,ts]=cyjzibian(t,x,u,flag,n,NUM)
switch flag
  case 0         
    [sys,x0,str,ts] = mdlInitializeSizes(n,NUM);
  case {1,2,4,9}
    sys = [];
  case 3
    sys = mdlOutputs(t,x,u,n,NUM);
    otherwise
    error(['unhandled flag = ',num2str(flag)]);
end
function [sys,x0,str,ts] = mdlInitializeSizes(n,NUM)
sizes = simsizes;
sizes.NumContStates  = 0;
sizes.NumDiscStates  =0;
sizes.NumOutputs     = 1;
sizes.NumInputs      = 1;
sizes.DirFeedthrough = 0;
sizes.NumSampleTimes = 1;
sys = simsizes(sizes);
str = [];
x0  = [];
ts  = [-1 0];  
function sys = mdlOutputs(t,x,u,n,NUM)
nH=720;f=50;p=4;n0=60*f/p;n=6;iMB=n0/n;lmdk=2.5;PH=45;omgg0=2*pi*n0/60;SH=(n0-nH)/n0;Sm=SH*(lmdk+sqrt(lmdk^2-1)); MH=9550*PH/nH;
Bw=0;l=1.025;m1=1440+1350*NUM;j1=1/3*m1*l*l;j2=2.12;Mc=l*m1*9.81;ycl=0.8702;tal=12/180*pi;H=6.1;G=2.6;I=3.2;K=sqrt((H-G)^2+I^2);
R=1.1;C=2.5;A1=4.4;P=3.6;theta0=0/180*pi;af=asin(I/K);
faimax=acos((C^2+K^2-(R+P)^2)/(2*C*K));faimin=acos((C^2+K^2-(P-R)^2)/(2*C*K));S=(faimax-faimin)*A1;
Er=2.12e+11;Dr=0.0256;Ar=pi*(Dr/2)^2;Lr=900;pr=7850;wr=pr*Ar*Lr*9.8;c=sqrt(Er/pr);
w=100e-3;Dc=0.05;Dti=0.076;v=(12*pi*w/(pr*Ar))*(Dr/(Dti-Dr))*((0.20+0.39*Dr/Dti)+2.1970*10e4/25*(Dc/Dti-0.3810)^2.57*((Dc^2-Dr^2)/(Lr*Dr)));
deltx=Lr/30;q=2/3;omg0=2*pi*nH/(60*iMB);bj=70;dym=700;ym=950;bjm=bj*bj*pi/4/1000000;ygjm=Dti*Dti*pi/4;
lmdt=bjm*dym*ym*9.8*Lr/Er/Ar;lmdg=bjm*dym*ym*9.8*Lr/Er/(ygjm);lmdtg=lmdt+lmdg;
x0=0.2;sp=S-lmdtg+x0+lmdg;W0=ym*9.8*(dym*ygjm-Ar*Lr);
load -ascii t1.txt;
load -ascii thetad.txt;
load -ascii thetd1.txt;
load -ascii thetd2.txt;
for i=1:3
    theta(i)=thetad(i);
    thet1(i)=thetd1(i);
    thet2(i)=thetd2(i);
end
deltt=t-t1;r=deltt/deltx;
for j=1:1:29           
    Ek=q*c^2*r^2;Ck=1+2*c^2*r^2*q+deltt*v*0.5;     
    A(j+1,j)=-Ek;A(j+1,j+2)=-Ek;A(j+1,j+1)=Ck;
end
    A(1,1)=1; A(31,29)=1/2*(Er*Ar);A(31,30)=-2*Er*Ar;A(31,31)=3/2*(Er*Ar);
  u=zeros(31,1000);up=zeros(1,1000);up1=zeros(1,1000);ppj=zeros(1,6000);
  u(:,1)=0;u(:,2)=0;
up(1)=0.0005;up(2)=0.001;up(3)=0.002,up1(1)=0;up1(2)=0;up1(3)=1.1;ppj(1)=0;ppj(2)=50;
      for i=1:3
        theta2(i)=2*pi-theta(i)+af;L(i)=sqrt(R^2+K^2-2*R*K*cos(theta2(i)));
        belta(i)=asin(R*sin(theta2(i))/L(i));theta3(i)=acos((P^2+L(i)^2-C^2)/(2*P*L(i)))-belta(i);
        theta4(i)=acos((P^2-L(i)^2-C^2)/(2*C*L(i)))-belta(i);x(i)=acos((C^2+L(i)^2-P^2)/(2*C*L(i)));fai(i)=x(i)+belta(i);
        PR(i)=(faimax-fai(i))*A1;u1(i)=(fai(i)-faimin)*A1;
        TF(i)=-R*(A1*sin(theta3(i)-theta2(i)))/(C*sin(theta3(i)-theta4(i)));
        theta21(i)=-thet1(i); theta31(i)=R*theta21(i)/P*sin(theta4(i)-theta2(i))/sin(theta3(i)-theta4(i));
        theta41(i)=R*theta21(i)/C*sin(theta3(i)-theta2(i))/sin(theta3(i)-theta4(i));
        vA(i)=theta41(i)*A1; theta22(i)=-thet2(i);
        theta42(i)=theta41(i)*(theta22(i)/theta21(i)-(theta31(i)-theta41(i))*cot(theta3(i)-theta4(i))+(theta21(i)-theta31(i))*cot(theta2(i)-theta3(i)));
        aA(i)=theta42(i)*A1;  
      end
    ii=3;
     D(1,1)=0;        
  for j=2:1:30
      D(j,1)=(2-2*c^2*r^2*(1-2*q))*u(j,ii-1)+(1-2*q)*c^2*r^2*(u(j+1,ii-1)+u(j-1,ii-1))+c^2*r^2*q*(u(j+1,ii-2)+u(j-1,ii-2))-(1+2*c^2*r^2*q-deltt*v*0.5)*u(j,ii-2)-((u1(ii)-2*u1(ii-1)+u1(ii-2))/(deltt^2)+v*(u1(ii)-u1(ii-1))/(deltt))*(deltt^2);
end
              
           
if up1(ii)>=0&up(ii)<=lmdt
    ppj(ii)=W0*((up(ii))/lmdt);
end
if up1(ii)>0&up(ii)>lmdt
    ppj(ii)=W0;
end
if up1(ii)<0&up(ii)>(sp-lmdt)
   ppj(ii)=W0*(up(ii)-(sp-lmdt))/lmdt;
end
if up1(ii)<0&up(ii)<=(sp-lmdt)
   ppj(ii)=0;
end
D(31,1)=deltx*ppj(ii);  
u(:,ii)=A\D;            
up(ii+1)=sp-u(31,ii)-u1(ii);up1(ii+1)=-(u(31,ii)-u(31,ii-1))/deltt-(u1(ii)-u1(ii-1))/deltt;
PRL(ii)=wr+Er*Ar/deltx*(-1.5*u(1,ii)+2*u(2,ii)-0.5*u(3,ii));
Mm(ii)=(2*lmdk*MH*Sm*omgg0*(omgg0-iMB*thet1(ii)))/(Sm^2*omgg0^2+(omgg0-iMB*thet1(ii))^2);
       if Mm(ii)>0
           Med(ii)=iMB./u(1)*(2*lmdk*MH*Sm*omgg0*(omgg0-iMB*thet1(ii)))/(Sm^2*omgg0^2+(omgg0-iMB*thet1(ii))^2);
       else
           Med(ii)=iMB*u(1)*(2*lmdk*MH*Sm*omgg0*(omgg0-iMB*thet1(ii)))/(Sm^2*omgg0^2+(omgg0-iMB*thet1(ii))^2);
       end
       P2=Mm(ii)*thet1(ii)*iMB;
             if vA(ii)>0
               Mef(ii)=TF(ii)*(PRL(ii)-Bw)*ycl^(-1)-Mc*sin(theta0+theta(ii)-tal);
            else
                Mef(ii)=TF(ii)*(PRL(ii)-Bw)*ycl^(1)-Mc*sin(theta0+theta(ii)-tal);
             end
                              m2=PRL(ii)/9.8;
                              Je(ii-1)=j1+j2*iMB*iMB+TF(ii-1)^2*m2;
                              Je(ii)=j1+j2*iMB*iMB+TF(ii)^2*m2;
                              dJe=(Je(ii)-Je(ii-1))/(theta(ii)-theta(ii-1));         
           theta(ii+1)=theta(ii)+thet1(ii)*deltt;
           thet1(ii+1)=thet1(ii)+thet2(ii)*deltt;
           thet2(ii+1)=1/Je(ii)*((Med(ii)-Mef(ii))-0.5*thet1(ii)*thet1(ii)*(dJe));
           if theta(ii+1)>=2*pi
              theta(ii+1)=theta(ii+1)-2*pi;
           end
     for i=1:3
         thetad(i)=theta(i+1);
         thetd1(i)=thet1(i+1);
         thetd2(i)=thet2(i+1);
     end
     save t1.txt t -ascii;
     save thetad.txt thetad -ascii;
     save thetd1.txt thetd1 -ascii;
     save thetd2.txt thetd2 -ascii;
     sys=P2;

:L 上述源代码中的两个笑脸处原为P，粘过来居然成了这样子

解决了么？ 可以联系我  QQ  ： 754962190

没解决呢啊，帮帮忙吧:)

出现这个报错有很多可能性：
１。 输出P２ 按照你的设置sizes.NumOutputs= 1 应为标量，而运算结果为向量
２。 输出应为实数，结果为复数，或者为空量 ［ ］。

检查方法，在P２ 那一行P2=Mm(ii)*thet1(ii)*iMB 的下一行
             if vA(ii)>0 加一个 breakpoint
运行模型，然后查看该P２的值是多少。 如果为空量，则你前面读取 输入量u() 可能有问题.

这个在用sfunction做仿真时遇到这种问题很常见，这个肯定是你的输出出现了NAN或者inf的数。

蝶起梦观涛 发表于 2009-7-13 21:08
我就是找不出问题在哪里，以前好像出过一次这样的情况，后来不记得怎么弄好了
我把源代码粘过来，亟盼各位 ...

学习一下啊