不可能啊 ，
你重新安装一下matlab吧

这是已知矩阵，
 
A =
 
[ 418.1625,  12.4299,     0,       0,         0,         0   ]
[ 12.4299,   35.4397,     0,       0,         0,         0   ]
[    0,         0,        0,       0,         0,         0   ]
[    0,         0,        0,   12.5800,       0,         0   ]
[    0,         0,        0,       0,    3.4000*C_55,    0   ]
[    0,         0,        0,       0,         0,      18.3600]


 
D =
 
[ 402.8299,  11.9742,  0,    0,     0,         0    ]
[ 11.9742,   34.1403,  0,    0,     0,         0    ]
[   0,         0,      0,    0,     0,         0    ]
[   0,         0,      0,  12.1187, 0,         0    ]
[   0,         0,      0,    0,  3.2753*C_55,  0    ]
[   0,         0,      0,    0,     0,       17.6868]

下面是程序，


 i=sqrt(-1);
 rho=0.0000015; 
 h=3.4;  
 f = (1:22);                % Frequencies for which experimental results of cg are available
%---------------- -------------------------------------------------------
    ka1=sqrt(pi^2/12);             % shear correction factor 
%--------------------------------------------------------------------------
    I0= rho*2*h/2;                 % intergal of I0 with (-z,z) interval.
    I2= rho*2*(h/2)^3/3;           % intergal of I2 with (-z,z) interval.
%--------------------------------------------------------------------------
% Symbolic expresion used to perform Algebraic calculation
   
    
    phi = 0;                       % Propagation wave angle
    
   for phi_i=phi:phi+1
 

    syms C_55;
    syms k ;
    syms omega;
    
         
    c=cos(phi_i/360*2*pi );
    s=sin(phi_i/360*2*pi );
    
    
     if phi_i==0
   
    kx=k*round(c);                 % wave number in x axis
    ky=k*round(s);                 % wave number in y axis    
     end
   
   
  
    kx=k*(c);                      % wave number in x axis
    ky=k*(s);                      % wave number in y axis
     

 
%--------------------------------------------------------------------------- 
    L33 = -ka1^2*A(5,5)*kx^2-2*ka1^2*A(4,5)*kx*ky-ka1^2*A(4,4)*ky^2+omega.^2.*I0;
    L34 = i*ka1*(ka1*A(5,5)*kx+ka1*A(4,5)*ky);
    L35 = i*ka1*(ka1*A(4,5)*kx+ka1*A(4,4)*ky);
    L44 = D(1,1)*kx^2+2*D(1,6)*kx*ky+D(6,6)*ky^2+ka1^2*A(5,5)-omega.^2.*I2;
    L45 = D(1,6)*kx^2+(D(1,2)+D(6,6))*kx*ky+D(2,6)*ky^2+i*ka1^2*A(4,5);
    L55 = D(6,6)*kx^2+2*D(2,6)*kx*ky+D(2,2)*ky^2+ka1^2*A(5,5)-omega.^2.*I2;
    

%--------------------------------------------------------------------------        
%Antisymmetric waves matrix:
  
   M   =  [L33  L34  L35;
           L34  L44  L45;
           L35  L45  L55];
   
   
   M=vpa(M);

   

%-------------------------------------------------------------------------
%Solve the determenent of the matrix and derive the roots

    det_M=vpa(det(M));
    
    Sol= solve(det_M,k);

    k_A0=sym(Sol(1,1));

    omega=2*pi*f;

    digits(5); 
    

    k_A0=sym((subs(k_A0,omega)));
    
    
     for p=1:length(f)
        
        Cp_A0((phi_i+1),p)=omega(p)/real(k_A0(p));
     end
%--------------------------------------------------------------------------
% Calcualtion of Group velocity values using and analytical differention of
% omega with respect to wavelength
    
    syms omega;
    
    diff_w=diff(Sol(1,1),omega);
        
    omega=2*pi*f;

    digits(5); 

    Cg_A0=subs(diff_w,omega);
    

   end
% --------------------------------------------------------------------------     
    for m=1:length(f) 
        
        Cg_A0_1(:,m)=1/Cg_A0(m);
        
    end
    
        Cg_A0=Cg_A0_1;
    
    end

%Group velocity calculation of Anti-symmetric waves.
   
      Cg2_A0=((Cp_A0(2,:)-Cp_A0(1,:))/(1/360*2*pi));
      Cgx_A0=(Cg_A0*cos(phi))-(Cg2_A0.*sin(phi));
      Cgy_A0=(Cg_A0*sin(phi))+(Cg2_A0.*cos(phi)); 
      C_55=3.0;
      Cgx_A0=subs(Cgx_A0,C_55);
      Cgy_A0=subs(Cgy_A0,C_55);
      Cg_A0_0=sqrt(Cgx_A0.^2+Cgy_A0.^2);
      save Cg_A0_0.mat Cg_A0_0;
     我觉得程序没什么问题呀，但是每次运行都超出内存范围，请高手指点一下，万分感谢！