Error using daeic12. This DAE appears to be of index greater than 1.

9 次查看（过去 30 天）
Lazar
Lazar 2024-5-5
评论： Lazar 2024-5-6
I'm trying to simulate the response of a power system to a three-phase fault at one of the power lines. The system in question can be described by a set of differential-algebraic equations (44 differential and 34 algebraic equations). The response consists of three intervals:
  • normal operation (before the fault),
  • faulted operation (during the three-phase fault at one of the power lines),
  • post-fault operation (after the disconnection of the faulted power line).
The fault and line disconnection are simulated by modifying the bus-admittance matrix which is used in evaluating the system algebraic equations. The differential equations remain unchanged during all of the three intervals.
The script for solving the DAEs over the three intervals is as follows:
%% Initialization
% Prepare the environment
clear variables
close all
clc
% Load the test system
load('Two-area system.mat','system') ;
% Perform the power flow to initialize the system and extract the initial
% conditions for state and algebraic variables
system = system.setInitialConditions() ;
z = [system.x0
system.y0] ;
u = system.u0 ;
% Determine the number of state and algebraic variables
numberOfStateVars = length(system.x0) ;
numberOfAlgVars = length(system.y0);
numberOfStateAndAlgVars = numberOfStateVars + numberOfAlgVars ;
% Creating a mass matrix which weights the derivatives of the state and
% algebraic variables
M = eye(numberOfStateAndAlgVars) ;
M(numberOfStateVars+1:end,numberOfStateVars+1:end) = 0 ;
% Define the options for ODE solver
options = odeset('Mass',M,'RelTol',1e-6,'AbsTol',1e-8) ;
%% Disturbance and interval definition
% Add the disturbances to the system
threePhaseFault = Fault(1,1.1,"3-phase fault",3) ;
lineOutage = Outage(1.1,inf,"Branch outage",5) ;
system = system.addDisturbance(threePhaseFault) ;
system = system.addDisturbance(lineOutage) ;
% Total simulation time
tFinal = 10 ;
% Start by adding the initial time
timePoints = [0
[system.disturbances.tStart]'
[system.disturbances.tEnd]'
tFinal] ;
% Remove duplicates and sort the time points
timePoints = unique(timePoints) ;
timePoints(isinf(timePoints)) = [] ;
%% Main loop
% Loop through each interval and solve the DAEs
tOutput = [] ;
zOutput = [] ;
for i = 1 : length(timePoints) - 1
% Define the current time span
tSpan = [timePoints(i) timePoints(i+1)-1e-6] ;
% Call the ODE solver
[tSol,zSol] = ode15s(@(t,z)system.evaluateDAEs(t,z,u),tSpan,z,options) ;
% Append the current solution
tOutput = [tOutput
tSol] ;
zOutput = [zOutput
zSol] ;
% Extract the values of state and algebraic variables at the end of the
% current interval
x = zOutput(end,1:numberOfStateVars) ;
y = zOutput(end,numberOfStateVars+1:end) ;
% While the state variables can't change immediately, the algebraic
% variables can. Determine the initial values for the algebraic
% variables after the fault is applied. The initial values of the state
% variables for the next interval remain unchanged.
y = fsolve(@(y) system.evaluateAlgEqns(x,y,u,timePoints(i+1)),y) ;
% Create a new initial vector
z = [x y]' ;
end
The code properly solves the DAEs during the normal operation. The initial conditions for the faulted interval are also properly determined (the algebraic variable corresponding to the voltage of the faulted bus is 0, which is expected). However, when the integration for the faulted interval starts, I receive an error:
Error using daeic12
This DAE appears to be of index greater than 1.
Error in ode15s (line 298)
[y,yp,f0,dfdy,nFE,nPD,Jfac] = daeic12(ode,odeArgs,t,ICtype,Mt,y,yp0,f0,...)
Can someone provide any insight in how to debug this issue?
  6 个评论
显示 4更早的评论
Lazar
Lazar 2024-5-6
Sorry, I got confused. dz values are finite and nonzero for some of the differential equations and zero for algebraic equations at the beginning of the faulted interval.
The differential equations for one machine are (4 machines x 11 differential equations = 44 equations):
The complex algebraic equation for one machine, later decoupled to two real equations, is (4 machines x 2 equations = 8 algebraic equations):
The complex algebraic equations, later decoupled to two real equations, for each bus in the system, are (13 buses x 2 equations = 26 algebraic equations):
At the start of the faulted interval, the only thing that changes is Y_ii value for the faulted bus i. Therefore, if the fault is at or near bus 3, Y_33 will be placed at a high value (I used 1e6 as recommended in one of the books) which effectively reduces the voltage V_3 to zero (this is consistent with the results of fsolve). Do you think that removing the two algebraic equations for the faulted bus during the faulted interval would solve the issue?
Lazar
Lazar 2024-5-6
@Torsten, I have managed to solve this error by evaluating dz using the newly obtained values of algebraic variables and providing dz as 'InitialSlope' to the ode15s solver at the beginning of the new interval. Now the solver successfully solves the DAEs for the faulted system.
When solving the DAEs for the post-fault system, solver failed several times and displayed an error:
Failure at t=1.167596e+00. Unable to meet integration tolerances without reducing the step size below thesmallest value allowed (3.552714e-15) at time t. > In ode15s (line 690)
I solved this error by slightly reducing the relative tolerance to 1e-3 and absolute tolerance to 1e-5.
In addition, I have found that the DAE system is extremely sensitive to the initial values of the algebraic variables y0 determined by fsolve. Since the algebraic equations are nonlinear, I need to provide a relatively accurate initial guess for y0 so that the solver works properly.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

产品


版本

R2023b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by