|
楼主 |
发表于 2010-6-10 09:07
|
显示全部楼层
这里问题主要体现在 上,两种不同的算法得出的结果完全不一样:
1. 用quadl计算
建立一个“k12.m”文件
%求解电磁刚度系数K12
function y1=k12(alpha)
Rg=6.24; %发电机定子内圆半径m
L1=2.1; %转子有效长度m
k_mu=1.102; %饱和度
mu_0=4*pi*10^(-7); %空气导磁系数
delta_0=18e-003; %均匀气隙大小m
f=50; %电网电流频率Hz
t=0; %计算时间,假设取个定值,也可以是其他值,如果是初始状态的话,那么时间t应该为0
theta=0.17*pi; %内功率角,用pi来表示
phi=pi/6.22; %功率因数角
p=40; %磁极对数
Fsm=19210; %发电机定子绕组三相基波磁势幅值At
Fjm=24214; %发电机转子绕组基波磁势幅值At
Lambda_0=mu_0./(k_mu*delta_0); %发电机均匀气隙磁导
sigma=k_mu*delta_0; %饱和度与均匀气隙大小的乘积
Wf=2*pi*f./p; %发电机同步转速
y1=(Rg.*L1.*Lambda_0/(2.*sigma^2))*sin(2.*alpha).*((Fsm.*cos(Wf.*t-p.*alpha)+Fjm*cos(Wf.*t-p.*alpha+theta+phi+pi./2)).^2);
%积分里面的东西涉及到乘法或除法运算时,“*”或“/”最好换成“.*”和“./”,否则有时计算会出现问题!
然后在主窗口输入K12=quadl(@k12,0,2*pi),得到
>> K12=quadl(@k12,0,2*pi)
Warning: Maximum function count exceeded; singularity likely.
> In quadl at 104
K12 =
-1.9148e+007
2用字符型积分int
首先定义一个参数文件“coefficient_elec_stiffness.m”,里面包括计算电磁刚度K12的参数
%%% 计算电磁刚度会用到的参数 %%%
Rg=6.24; %发电机定子内圆半径m
L1=2.1; %转子有效长度m
k_u=1.102; %饱和度
mu_0=4*pi*10^(-7); %空气导磁系数
delta_0=18e-003; %均匀气隙大小m
f=50; %电网电流频率Hz
t=0; %计算时间,假设取个定值,也可以是其他值,如果是初始状态的话,那么时间t应该为0
% t=2;
theta=30.64/180*pi; %内功率角,用pi来表示
phi=acos(0.875); %功率因数角
p=40; %磁极对数
Fsm=19210; %发电机定子绕组三相基波磁势幅值At
Fjm=24214; %发电机转子绕组基波磁势幅值At
Lambda_0=mu_0./(k_u*delta_0); %发电机均匀气隙磁导
sigma=k_u*delta_0; %饱和度与均匀气隙大小的乘积
omega_f=2*pi*f./p; %发电机同步转速
然后在主窗口中输入如下内容,或者定义为另一个“solution_elec_stiffness.m”文件运行,
clc
clear
syms alpha %气隙宽度等于delta的周向位置与x轴之间的夹角,积分项
coefficient_elec_stiffness %%% 求解电磁刚度时采用的参数 %%%%
digits(8) %设置计算结果精度,保留8位有效数字
%%%%% 小写字母代表被积分内容,大写字母代表积分结果 %%%%%%
integral_same=(Fsm.*cos(omega_f.*t-p.*alpha)+Fjm*cos(omega_f.*t-p.*alpha+theta+phi+pi./2)).^2; %%%积分相同项
k12=vpa((Rg.*L1.*Lambda_0/(2.*sigma^2))*(sin(2*alpha)).*integral_same)
K12=vpa(int(k12,alpha,0,2*pi))
得到结果是
k12 =
1.0549247*sin(2.0*alpha)*(19210.0*cos(40.0*alpha) + 24214.0*cos(2.6109257 - 40.0*alpha))^2
K12=0
因为这两个积分的表达式是完全相同的,只不过前面的都用数值代替了,后面的用符号算,但是结果相差巨大 |
|