campell图求解

admin

[mw_shl_code=applescript,true]function Crisp = campell(K1,K2,numOflines,min_freq,criticalSpd_noLimit)

fid = fopen('F:\123\2.txt','wt');
str1 = ['K1=',num2str(K1),'\n'];
str2 = ['K2=',num2str(K2),'\n'];
fprintf(fid,str1);
fprintf(fid,str2);
fclose(fid);
fid = fopen('F:\123\4.txt','wt');
fid1 = fopen('F:\123\1.txt','rt');
data1 = fread(fid1);
fwrite(fid,data1);
fid2 = fopen('F:\123\2.txt','rt');
data2 = fread(fid2);
fwrite(fid,data2);
fid3 = fopen('F:\123\3.txt','rt');
data3 = fread(fid3);
fwrite(fid,data3);
fclose(fid);
fclose(fid1);
fclose(fid2);
fclose(fid3);

dos('"C:\Program Files\ANSYS Inc\v130\ansys\bin\winx64\ANSYS130.exe" -b -i "F:\123\4.txt" -o "G:\123\123.txt"');
a = textread('123.txt');
a = a';

N = numOflines;

line_width = 1.5;



w = 0:14;
w = 5000*w*2*pi/60;
freq_line = w/2/pi;


max_freq = max(N*freq_line);
[tempi,tempj] = find(a(:,1)<=max_freq);
a = a(unique(tempi),;
clear tempi tempj
[tempi,tempj] = find(a(:,1)>=min_freq);
a = a(unique(tempi),;
clear tempi tempj
[m,n] = size(a);
tempi = [];
for ii = 1:m
    if abs(a(ii,1)-a(ii,n))<=5
        tempi = [tempi ii];
    end
end
a(tempi, = [];
clear tempi tempj

slp = (a(:,end) - a(:,1))/(w(end) - w(1));
[tempi,tempj] = find(slp>0);
Crisp = a(:,1)./(1/2/pi-slp)/2/pi;
Crisp = Crisp(tempi);
clear tempi tempj

[nmod,solnum] = size(a);

figure
hold on
for ii = 1:nmod
    slope(ii) = a(ii,end) - a(ii,1);
    if slope(ii)>=0
        col = 'r-';
        h1 = plot(w,a(ii,,col,'LineWidth',line_width);
    else
        col = 'k-';
        h2 = plot(w,a(ii,,col,'LineWidth',line_width);
    end
   
end
slope;
for jj = 1:N
    h3 = plot(w,jj*freq_line,'--g','LineWidth',line_width);
end

legend([h1 h2 h3],'正进动','反进动','等转速线',2)
fig_title = ['K1=',num2str(round(K1)),'Campell图'];
title(fig_title)
axis tight
xlabel('转速 rad/s')
ylabel('频率 Hz')
grid on
hold off
Crisp = Crisp(1:criticalSpd_noLimit);
return
end
[/mw_shl_code]