标签：仿真 数值 k1 MATLAB delta Y1 omega 欧拉 out
1 发电机方程
有发电机的二阶微分方程如下：
其中δ为发电机转子角（rotor angle），ω为角速度（rotor angular speed），H为发电机惯性常数（inertia constant），D为发电机阻尼系数（damping constant）。
2 欧拉法
将微分方程离散化，在离散点处用差商代替导数。其中，根据不同的差商，欧拉法分为向前欧拉和向后欧拉。
2.1 向前欧拉法
2.2 向后欧拉法
向后欧拉是一种隐式格式，计算每一个迭代步时都需要求解一个非线性方程。
2.3 两点欧拉公式
两点欧拉公式需要用到两个初值，其中y0给定，y1一般通过向前欧拉算得。
3 欧拉预估-校正法
由于梯形公式是一种隐式格式，实际计算时不是很方便，因此将其修正成显示格式。
4 发电机方程的数值求解
向前欧拉法：
iter = 100;
h = 1/100;
x0 = [x0_af(1);x0_af(4)];
x(:,1) = x0;
omega = x(1,:);
delta = x(2,:);
Y1 = [ 1.0608 - 2.6305i 0.1082 + 0.6064i 0.1572 + 0.9459i];
K=[];
for i = 1:iter
k1 = f1(omega(i),delta(i),i,Y1);
k2 = f1(omega(i)+0.5*h*k1(1),delta(i)+0.5*h*k1(2),i,Y1);
k3 = f1(omega(i)+0.5*h*k2(1),delta(i)+0.5*h*k2(2),i,Y1);
k4 = f1(omega(i)+h*k3(1),delta(i)+h*k3(2),i,Y1);
k22 = f1(omega(i)+h*k1(1),delta(i)+h*k1(2),i,Y1);
%x(:,i+1) = x(:,i) + (h/6)*(k1+2*k2+2*k3+k4); %四阶龙格库塔方法
x(:,i+1) = x(:,i) + h*k1; %欧拉方法
% x(:,i+1) = x(:,i) + 0.5*h*(k1+k22); %改进的欧拉方法
omega = x(1,:);
delta = x(2,:);
K = [K,k1];
end
子函数
%% 4-order Runge-Kutta method 求k1 k2 k3 k4的函数
% 在power_grid ode45_RK4中被调用
%%
function k = f1(omega,delta,i,Y1)
k = zeros(2,1);
load Xaf_out
delta1 = [delta,Xaf_out(i+1,5),Xaf_out(i+1,6)];
E_pf =[1.0566,1.0502,1.0170];
clear sum
Pe1 = E_pf(1).*sum(E_pf.*(real(Y1).*cos(delta*ones(1,3)-delta1)+imag(Y1).*sin(delta*ones(1,3)-delta1)));
k(1) = (pi*60/23.6400)*(0.6611 - Pe1);
k(2) = omega - 2*pi*60;
end
5 ode45
ode45是MATLAB自带的解微分方程的函数，用的是四阶龙格库塔方法
[Taf_out , Xaf_out] = ode45(@(t,x)Gen_Fun(t,x,GEN,Yraf,E_pf),tspan,x0_af,options);
子函数：
function Xdf_out = Gen_Fun(t,x,GEN,Ynn,E_pf)
Xdf_out = zeros(6,1);
omega_r = 2*pi*60;
delta_ij_matrix = zeros(3);
for i = 1:3
for j = 1:3
delta_ij_matrix(i,j) = x(3+i) - x(3+j);
end
end
H = GEN(:,4);
Pm = GEN(:,5);
Pe = E_pf*E_pf'.*(real(Ynn).*cos(delta_ij_matrix) + imag(Ynn).*sin(delta_ij_matrix));
Pe = sum(Pe,2);
Xdf_out(1:3) = omega_r*(Pm - Pe)./(2*H);
Xdf_out(4:6) = x(1:3) - omega_r;
end
如上图是利用欧拉法求解，和ode45求解的对比。
标签：仿真,数值,k1,MATLAB,delta,Y1,omega,欧拉,out 来源： https://blog.csdn.net/weixin_42749944/article/details/115142299
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。