CN107482673B - Economic dispatching method for multi-region fully-distributed active power distribution network - Google Patents

Economic dispatching method for multi-region fully-distributed active power distribution network Download PDF

Info

Publication number
CN107482673B
CN107482673B CN201710605843.8A CN201710605843A CN107482673B CN 107482673 B CN107482673 B CN 107482673B CN 201710605843 A CN201710605843 A CN 201710605843A CN 107482673 B CN107482673 B CN 107482673B
Authority
CN
China
Prior art keywords
distribution network
power distribution
region
active power
representing
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201710605843.8A
Other languages
Chinese (zh)
Other versions
CN107482673A (en
Inventor
吴文传
孙宏斌
许桐
郭庆来
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tsinghua Berkeley Shenzhen College Preparatory Office
Original Assignee
Tsinghua Berkeley Shenzhen College Preparatory Office
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Tsinghua Berkeley Shenzhen College Preparatory Office filed Critical Tsinghua Berkeley Shenzhen College Preparatory Office
Priority to CN201710605843.8A priority Critical patent/CN107482673B/en
Publication of CN107482673A publication Critical patent/CN107482673A/en
Application granted granted Critical
Publication of CN107482673B publication Critical patent/CN107482673B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • H02J3/383
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00Circuit arrangements for ac mains or ac distribution networks
    • H02J3/04Circuit arrangements for ac mains or ac distribution networks for connecting networks of the same frequency but supplied from different sources
    • H02J3/06Controlling transfer of power between connected networks; Controlling sharing of load between connected networks
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00Circuit arrangements for ac mains or ac distribution networks
    • H02J3/38Arrangements for parallely feeding a single network by two or more generators, converters or transformers
    • H02J3/46Controlling of the sharing of output between the generators, converters, or transformers
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2203/00Indexing scheme relating to details of circuit arrangements for AC mains or AC distribution networks
    • H02J2203/20Simulating, e g planning, reliability check, modelling or computer assisted design [CAD]
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02EREDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E10/00Energy generation through renewable energy sources
    • Y02E10/50Photovoltaic [PV] energy
    • Y02E10/56Power conversion systems, e.g. maximum power point trackers

Landscapes

  • Engineering & Computer Science (AREA)
  • Power Engineering (AREA)
  • Supply And Distribution Of Alternating Current (AREA)

Abstract

The invention provides an economic dispatching method for a multi-region fully-distributed active power distribution network, belonging to the technical field of operation and control of power systems, and comprising the following steps: 1) establishing an economic dispatching model of the active power distribution network considering the linearization network loss and the minimum light abandoning cost; 2) realizing each area through the communication of the boundary nodes among the neighbors: carrying out iteration initial value construction on a localized quasi-Newton matrix; performing localized quasi-Newton direction iterative computation; calculating a localization search step length; 3) after each area is converged, whether the whole situation is converged is judged through communication between neighbors. The method disclosed by the invention is aimed at the operation characteristics of the active power distribution network, fully utilizes the point-to-point communication technology, is based on a fully-distributed quasi-Newton method, is used for completing multi-region fully-distributed ultra-linear convergence calculation of economic dispatching of the active power distribution network, does not need a coordination center, improves the reliability, protects the privacy of each region, and realizes plug and play.

Description

Economic dispatching method for multi-region fully-distributed active power distribution network
Technical Field
The invention relates to an economic dispatching method for a multi-region fully-distributed active power distribution network, and belongs to the technical field of operation and control of power systems.
Background
With the increasing addition of distributed power sources, controllable loads, electric vehicles and other distributed resources to the power distribution network, the traditional passive power distribution network is changed into a more active form, and each user in the power distribution network is not only a consumer but also a production consumer, namely an active power distribution network. The new power distribution network form requires that participants interact with the power grid effectively, distributed renewable energy sources can be consumed, and the plug-and-play function can be realized, so that the traditional centralized energy management mode is limited by the problems of increased computational burden and communication delay. Meanwhile, the power distribution side will present the situation of multi-benefit subjects (such as load aggregators) in the future, and the data of each subject will have valuable value, so the privacy protection problem is more prominent, and it is a very good solution idea to realize the collaborative optimization among multiple subsystems or areas through limited information communication among neighbors, but most of the mature algorithms nowadays only have linear or sub-linear convergence characteristics, the computing speed of a large system still needs to be improved, and it is very important to research the fully distributed optimization algorithm with the super-linear convergence property.
At present, a centralized or decomposition coordinated calculation method is mostly adopted for the power distribution network, because of the existence of a central controller, the central controller cannot meet the requirements of the future active power distribution network, and the traditional economic dispatching method can meet the bottlenecks of large calculation amount, high communication delay, poor information protection, low reliability and the like. Based on point-to-point communication technology, it is a trend and necessity to develop a fully distributed approach that does not rely on a central controller. Therefore, a new economic dispatching method of the active power distribution network with ultra-linear convergence and suitable for multiple areas needs to be researched.
Disclosure of Invention
The invention aims to provide a multi-region fully-distributed active power distribution network economic dispatching method, which overcomes the defects of the prior art, meets the aim of multi-region economic dispatching calculation of an active power distribution network by using a point-to-point communication technology and based on a fully-distributed quasi-Newton method, realizes that each region completes global optimization solution only through a small amount of information interaction between neighbors, can reduce the light abandon condition as far as possible, and maximizes the renewable energy utilization rate.
The invention provides an economic dispatching method for a multi-region fully-distributed active power distribution network, which comprises the following steps of:
(1) establishing an economic dispatching model of the active power distribution network, which comprises the minimum light abandoning cost:
the objective function F of the economic dispatching model of the active power distribution network is as follows:
wherein T represents the input to the active distribution networkThe total time of economic dispatch is passed, t represents any time in the total time of economic dispatch, i is the ith distributed conventional small generator in the power distribution network,representing the collection of all distributed conventional small generators in the distribution grid,representing a cost function of a conventional distributed small generator in a power distribution network, the cost function being calculated by the formula:wherein a isi、biAnd ciRespectively, the cost parameter, a, of the ith distributed conventional small generator in the power distribution networki、biAnd ciThe values of (A) are respectively 0.01-1, 20-20.5 and 600-1500,the planned active power of the ith distributed conventional small generator in the power distribution network at the moment t, and j is the jth distributed photovoltaic power generation equipment in the power distribution networkRepresenting a collection of distributed photovoltaic power plants in a power distribution network, omegajThe light abandonment penalty coefficient of the jth photovoltaic power generation equipment is 500-1000,representing the planned active power of the distributed photovoltaic power generation equipment in the power distribution network at time t,representing a predicted value of the distributed photovoltaic power generation equipment in the power distribution network at the time t;
the constraint conditions of the economic dispatching model of the active power distribution network comprise:
a) the load flow balance constraint condition considering the network loss of the active power distribution network is as follows:
wherein the content of the first and second substances,the active power of the line between the node m and the node h in the power distribution network at the moment t,representing the active power of the line between node n and node m in the distribution network at time t,representing the active power loss of the line between node n and node m in the distribution network at time t, andrespectively representing historical data points of the active power and the reactive power of a line between a node n and a node m in the power distribution network and the voltage of the node n,andr is obtained from historical data of a power distribution network scheduling systemmnRepresenting the line resistance between node n and node m in the distribution network,representing the active power output of node m in the distribution network at time t,
representing the active power of the generator at the node m in the distribution network at the moment t,is a variable to be solved for the method,representing the load active power demand of a node m in the power distribution network at the time t;
centralizing variables to be solved in the power flow balance constraint condition considering the network loss of the active power distribution network, and rewriting the variables to be solved into the following matrix form:
Bsx=b
wherein, BsRepresenting an extended adjacency matrix, Bs=[-Z B]Any element Z in the matrix ZopComprises the following steps:
arbitrary element B in matrix BnuComprises the following steps:
x represents the set of all variables to be solved, x is a vector,b represents the deviation of the active power of the nodes in the power distribution network, b is a vector,
b) the inequality constraint condition of the network where the active power distribution network is located is considered as follows:
the climbing constraint conditions of the conventional unit in the active power distribution network are as follows:
wherein, mudiAnd muiiRespectively representing the climbing coefficient and the descending coefficient of the ith conventional unit in the power distribution network, wherein delta t is the set economic dispatching time interval of the power distribution network;
the constraint conditions of the active power of the conventional unit in the active power distribution network are as follows:
wherein the content of the first and second substances,andrespectively representing the upper limit of active power and the lower limit of active power of the ith conventional unit in the power distribution network;
the constraint conditions of the line tide capacity in the active power distribution network are as follows:
wherein the content of the first and second substances,andrepresenting the upper limit/the lower limit of active power of a line power flow between a node n and a node m in the power distribution network;
the active power constraint conditions of the photovoltaic units in the active power distribution network are as follows:
writing inequality constraint conditions of the network where the active power distribution network is located into a matrix form as follows:
(2) and (3) converting the economic dispatching model of the active power distribution network into the following KKT condition equation by adopting a primal-dual interior point method:
wherein r ispriIs the original residual of the KKT conditional equation, rdualAnd rcentThe dual residual and the central residual of the KKT conditional equation, respectively, Dq represents the function gradient matrix,λ and ν are the dual multipliers of the above inequalities and equality constraints, respectively, γ is a calculation parameter,mu is a factor, the value of mu is 10, and the value of g is equal to twice of the total number Y of the inequality constraint conditions;
solving the KKT conditional equation by using a numerical iteration method, where y is (x, λ, v) and Δ y is (Δ x, Δ λ, Δ v), and then the iterative solution formula of the KKT conditional equation is:
yk+1=yk+dk△yk
where k is the number of iterations, ykRepresents the value of y at the k iteration, dkRepresents the value of step length at the kth iteration,. DELTA.ykRepresenting the search direction at the k-th iteration;
obtaining the search direction delta y according to a first-order Taylor expansion formula by using the following formulak
△yk=-(Drγ(yk))-1rγ(yk)
Wherein the content of the first and second substances, representing a quadratic gradient calculation, superscriptTRepresenting matrix transposition, diag representing an operation on a diagonal matrix, BsRepresents an extended adjacency matrix;
(3) solving y in the KKT conditional equation to realize economic dispatching of the multi-region fully-distributed active power distribution network, firstly defining the connection lines between the regions in the active power distribution network as follows: the method comprises the following steps of performing region allocation on the links between regions in the active power distribution network, wherein the links between any two regions belong to the region with the smaller number of nodes:
(3-1) setting the economic dispatching initial time t to be 0 during initialization;
(3-2) setting the initial iteration step number k of the active power distribution network multi-region fully-distributed economic scheduling calculation at the time t to be 0;
(3-3) solving the approximation matrix H, comprising the steps of:
(3-3-1) constructing an approximate matrix initial value H of a region N in an active power distribution networkN(0):
HN(0)=MN+LN
Wherein M isNIs a block diagonal matrix comprising an area N and other areas connected to the area N, each diagonal block in the block diagonal matrix being defined by the Dr of each area in the distribution networkγ(y0) Namely Drγ(y) initial value composition, LNIs a non-block diagonal symmetric matrix, the elements in the non-block diagonal symmetric matrix are formed by the initial value H of an approximate matrixN(0) The element composition in the B matrix of the corresponding position;
(3-3-2) dual residual errors and original residual errors r in KKT conditional equation of region N in active power distribution networkdualAnd rpriThe local dual residual r after correction is obtained by respectively correctingdualNAnd the original residual rpriN
Wherein, subscript N indicates that the element is an element in an area N in the active power distribution network, and vector MrdualNAny element in (1) is MrdualN(pN+n),
Region Q is connected to region N, pNIs the number of generators in the N region, vector MrpriNAny element in (1) is MrpriN(n),
pQIs the number of generators in the Q region, laNQIs the above BsThe values of the boundary elements at the junction of region N and region Q in the matrix,Nnand QmRespectively representing a boundary node N of an area N and a boundary node m of an area Q connected with the area N in the power distribution network;
(3-3-3) obtaining an area N approximate matrix H in the power distribution network of the (k + 1) th iteration by utilizing a quasi-Newton method and calculating according to the following formulaN(k+1):
Wherein the content of the first and second substances,tau is an adjustment parameter and takes the value of 10-3I is an identity matrix, subscriptnNRepresenting the extension variables of boundary nodes and corresponding lines which comprise an area N and an area connected with the area N in the power distribution network;
(3-4) calculating the extended search direction of the area N in the power distribution network by using the approximate matrix H obtained in the step (3-3) according to the following formula:
wherein, gamma is an adjusting parameter with the value of 10-4The extension variables representing the k-th iteration including the region N and the boundary nodes of the regions connected to the region N in the distribution network and the corresponding lines, the active power variables of the generators representing the region N and the boundary nodes of the regions connected with the region N in the power distribution network and the active power variable of the power flow of the corresponding lines,andrespectively representing the dual variables of an area N, an area boundary node connected with the area N and a corresponding line in the power distribution network;
(3-5) Using the extended search Direction obtained in the above step (3-4)The search direction for region N is calculated using the following equation:
wherein, wNRepresenting the total number of regions including region N and the regions contiguous to region N,the search direction of the variable in the area N corresponding to the expanded search direction of the distribution network area f in the k iteration is obtained;
(3-6) calculating the step length in the iterative solution formula of the KKT conditional equation of the region N in the power distribution network by using the following formulaThe method comprises the following steps:
(3-6-1) respectively calculating first middle step lengths of all areas in the power distribution network, wherein the first middle step length of the area N is recorded as
(3-6-2) exchanging the calculated first intermediate step length between all the areas and all the other connected areas;
(3-6-3) taking region N as an example, region N receives the first intermediate step size of all other connected regions, and calculates the second intermediate step size of region N by using the following formula
(3-6-4) according to the second intermediate step sizeRegion N calculates search step size for region N Using formulasFor search step lengthUpdating and utilizing according to the updated search step lengthUpdatingUsing updatedCalculating q (y)k+1) Up to q (y)k+1)≤0;
(3-6-5) according to the search step sizeUsing formulasFor search step lengthUpdating and utilizing according to the updated search step lengthUpdatingUsing updatedCalculate | | | rγ(yN k+1)||2Up toExchange with connected areaAnd updating the search step size Wherein, alpha and beta are backtracking parameters with the values of 0.01-0.1 and 0.3-0.8;
(3-6-6) search variable Deltay according to the region NN kAnd search step d of region NN kSolution of the above KKT conditional equation yN kUpdating:
yN k+1=yN k+dN k△yN k
(3-7) calculating residual total r of the region N in the power distribution networkfeasAnd proxy gap Setting a threshold value epsilon for convergence of iterative calculation of a regionfeasAccording to a threshold value epsilonfeasDetermining whether the region iteration calculation converges, wherefeasValue of 10-8~10-6
(3-8) based on the total residual rfeasAnd proxy gapSolution y to the above KKT conditional equationN kIs judged if r is the convergence state offeas≤εfeasAnd is andsolution y of the KKT conditional equation for region N in the distribution networkN kConvergence, setting region N convergence flag Conv in distribution networkNGo to step (3-9) if r is 1feasfeasAndone or both of them are satisfied, then the region N convergence flag Conv is assertedNGo to step (3-3-2) when k is 0 and k is k + 1;
(3-9) Conv according to the local convergence statusNAnd (3) performing global convergence state calculation:
the communication exchange local convergence state of the defined region N and the connected region is calculated as a global convergence mark:
(3-10) if the global convergence flag1, completing economic dispatching solution of the multi-region fully-distributed active power distribution network at the time t to obtain planned active power of each conventional generator in the power distribution network at the time tAnd planned active power of the photovoltaic power generation plantWhile setting t +1 to go to (3-11), ifMaking k equal to k +1, turning to the step (3-3-2), and continuing to perform the next iterative computation;
(3-11) if T is equal to T, finishing the calculation, and obtaining the planned active power of each conventional generator in the power distribution network at each time T in the T timeAnd planned active power of the photovoltaic power generation plantComplete moreEconomic dispatch of the regional fully distributed active distribution network if t<And T, turning to the step (3-2).
The economic dispatching method for the multi-region fully-distributed active power distribution network provided by the invention has the advantages that:
the economic dispatching method of the multi-region fully-distributed active power distribution network provided by the invention has the advantages that renewable energy sources can be utilized to the maximum extent in the multi-region economic dispatching of the active power distribution network, the light loss is greatly reduced, meanwhile, the economic dispatching method has the characteristics of high calculation efficiency and high convergence speed, and because a central controller does not exist, each region only exchanges limited boundary information by using a point-to-point communication technology, the communication traffic is small, the time delay is low, the privacy information is protected, the requirement of plug-and-play characteristics can be met, and in addition, the accurate dispatching can be realized under the condition of asynchronous communication or communication failure. In conclusion, the method and the device can play an important role in multi-area economic dispatching of the active power distribution network.
Drawings
Fig. 1 is a schematic diagram of an active power distribution network economic dispatching method according to the present invention, in which a region is divided into three regions.
Fig. 2 is a block diagram of an economic dispatching calculation process of a multi-region fully-distributed active power distribution network according to the present invention.
Fig. 3 is an exemplary schematic diagram of two regions for calculating an initial value of a KKT matrix according to the method of the present invention.
FIG. 4 is a schematic diagram illustrating two exemplary areas of search direction calculation information interaction according to the method of the present invention.
Detailed Description
The economic dispatching method for the multi-region fully-distributed active power distribution network, provided by the invention, considers an economic dispatching model of maximum renewable energy consumption and a fully-distributed multi-region optimization solving method. The fully-distributed multi-region optimization solving method completes the global optimization solving process through communication between neighbor regions, does not need a global central controller, is based on a quasi-Newton method with super-linear convergence, and can adapt to asynchronous communication conditions to a certain degree, and comprises the following steps:
(1) establishing an economic dispatching model of the active power distribution network, which comprises the minimum light abandoning cost:
the objective function F of the economic dispatching model of the active power distribution network is as follows:
wherein T represents the total time of economic dispatching on the active power distribution network, T represents any time in the total time of economic dispatching, i is the ith distributed conventional small generator in the power distribution network,representing the collection of all distributed conventional small generators in the distribution grid,representing a cost function of a conventional distributed small generator in a power distribution network, the cost function being calculated by the formula:wherein a isi、biAnd ciRespectively, the cost parameter, a, of the ith distributed conventional small generator in the power distribution networki、biAnd ciThe values of (a) are respectively 0.01-1, 20-20.5 and 600-1500, in one embodiment of the invention, ai、biAnd ciThe values of (A) are respectively 0.05, 20.3 and 800,j is the jth distributed photovoltaic power generation equipment (also can be distributed renewable energy generators in the power distribution network such as a small wind driven generator) in the power distribution network,representing a collection of distributed photovoltaic power plants in a power distribution network, omegajThe light abandonment penalty coefficient of the jth photovoltaic power generation equipment is 500-1000, and the photovoltaic power generation equipment is provided withIn one embodiment of the invention, ωjThe value of (a) is 600,representing the planned active power of the distributed photovoltaic power generation equipment in the power distribution network at the time t, is a variable to be solved by the method,representing a predicted value of the distributed photovoltaic power generation equipment in the power distribution network at the time t;
the constraint conditions of the economic dispatching model of the active power distribution network comprise:
a) the load flow balance constraint condition considering the network loss of the active power distribution network is as follows:
wherein the content of the first and second substances,the active power of the line between the node m and the node h in the power distribution network at the moment t,representing the active power of the line between node n and node m in the distribution network at time t,is a variable to be solved for the method,representing the active power loss of the line between node n and node m in the distribution network at time t,
andrespectively representing historical data points of the active power and the reactive power of a line between a node n and a node m in the power distribution network and the voltage of the node n,andr is obtained from historical data of a power distribution network scheduling systemmnRepresenting the line resistance between node n and node m in the distribution network,representing the active power output of node m in the distribution network at time t, representing the active power of the generator at the node m in the distribution network at the moment t,is a variable to be solved for the method,representing the load active power demand of a node m in the power distribution network at the time t;
centralizing variables to be solved in the power flow balance constraint condition considering the network loss of the active power distribution network, and rewriting the variables to be solved into the following matrix form:
Bsx=b
wherein, BsRepresenting an extended adjacency matrix, Bs=[-Z B]Any element Z in the matrix ZopComprises the following steps:
arbitrary element B in matrix BnuComprises the following steps:
x represents the set of all variables to be solved, x is a vector,b represents the deviation of the active power of the nodes in the power distribution network, b is a vector,
b) the inequality constraint condition of the network where the active power distribution network is located is considered as follows:
the climbing constraint conditions of the conventional unit in the active power distribution network are as follows:
wherein, mudiAnd muiiRespectively representing the climbing coefficient and the descending coefficient of the ith conventional unit in the power distribution network, wherein delta t is the set economic dispatching time interval of the power distribution network;
the constraint conditions of the active power of the conventional unit in the active power distribution network are as follows:
wherein the content of the first and second substances,andrespectively representing the upper limit of active power and the lower limit of active power of the ith conventional unit in the power distribution network;
the constraint conditions of the line tide capacity in the active power distribution network are as follows:
wherein the content of the first and second substances,andrepresenting the upper limit/the lower limit of active power of a line power flow between a node n and a node m in the power distribution network;
the active power constraint conditions of the photovoltaic units in the active power distribution network are as follows:
writing inequality constraint conditions of the network where the active power distribution network is located into a matrix form as follows:
(2) and (3) converting the economic dispatching model of the active power distribution network into a KKT (Karush-Kuhn-Tucker) condition equation as follows by adopting a primal-dual interior point method:
wherein r ispriIs the original residual of the KKT conditional equation, rdualAnd rcentThe dual residual and the central residual of the KKT conditional equation, respectively, Dq represents the function gradient matrix,λ and ν are the dual multipliers of the above inequalities and equality constraints, respectively, γ is a calculation parameter,mu is a factor, the value of mu is 10, and the value of g is equal to twice of the total number Y of the inequality constraint conditions;
solving the KKT conditional equation by using a numerical iteration method, where y is (x, λ, v) and Δ y is (Δ x, Δ λ, Δ v), and then the iterative solution formula of the KKT conditional equation is:
yk+1=yk+dk△yk
where k is the number of iterations, ykRepresents the value of y at the k iteration, dkRepresents the value of step length at the kth iteration,. DELTA.ykRepresenting the search direction at the k-th iteration;
obtaining the search direction delta y according to a first-order Taylor expansion formula by using the following formulak
△yk=-(Drγ(yk))-1rγ(yk)
Wherein the content of the first and second substances, representing a quadratic gradient calculation, superscriptTRepresenting matrix transposition, diag representing an operation on a diagonal matrix, BsRepresents an extended adjacency matrix;
(3) solving y in the KKT conditional equation to realize economic dispatching of the multi-region fully-distributed active power distribution network, firstly defining the connection lines between the regions in the active power distribution network as follows: area distribution is carried out on the links between areas in the active power distribution network, the links between any two areas belong to the areas with smaller number of nodes, as shown in fig. 1, the nodes n, m and h in the power distribution network respectively belong to the areas A, B, C and h<n<m, the region A passes through the connecting line lnm,lnhConnecting to the region B and the region C, the connecting line lnmBelonging to the area A, the tie line lnhBelonging to the region C. Comprising the following steps, as shown in figure 2:
(3-1) setting the economic dispatching initial time t to be 0 during initialization;
(3-2) setting the initial iteration step number k of the active power distribution network multi-region fully-distributed economic scheduling calculation at the time t to be 0;
(3-3) solving the approximation matrix H, comprising the steps of:
due to (Dr)γ(y))-1Direct full-distributed multi-area solution cannot be realized, and an approximate matrix H needs to be searched to replace DrγAnd (y), finally solving the KKT condition equation.
(3-3-1) constructing an approximate matrix initial value H of a region N in an active power distribution networkN(0):
The construction of an initial value matrix can be completed only through information interaction with a boundary node of a neighbor region, and the initial value H of an approximate matrix of a region NN(0) The device is composed of two parts:
HN(0)=MN+LN
wherein M isNIs a block diagonal matrix comprising an area N and other areas connected to the area N, each diagonal block in the block diagonal matrix being defined by the Dr of each area in the distribution networkγ(y0) Namely Drγ(y) initial value composition, LNIs a non-block diagonal symmetric matrix, the elements in the non-block diagonal symmetric matrix are formed by the initial value H of an approximate matrixN(0) The element composition in the B matrix of the corresponding position; as shown in fig. 3, for two regions HN(0) As an example of the initial value,Aintrespectively representing a boundary node n of an A area and an internal node of the A area in the power distribution network,representing a boundary node m in the distribution network,
(3-3-2) due to BsPresence of a network coupling element, rprirdualRealizing the residual localization calculation needs to realize certain correction with the information communication of the neighboring boundary region, rcentNo correction is required.
Dual residual and original residual r in KKT conditional equation of region N in active power distribution networkdualAnd rpriThe local dual residual r after correction is obtained by respectively correctingdualNAnd the original residual rpriN
Wherein, subscript N indicates that the element is an element in an area N in the active power distribution network, and vector MrdualNAny element in (1) is MrdualN(pN+n),
Region Q is connected to region N, pNIs the number of generators in the N region, vector MrpriNAny element in (1) is MrpriN(n),
pQIs the number of generators in the Q region, laNQIs the above BsThe values of the boundary elements at the junction of region N and region Q in the matrix,Nnand QmRespectively representing a boundary node N of an area N and a boundary node m of an area Q connected with the area N in the power distribution network;
(3-3-3) obtaining an area N approximate matrix H in the power distribution network of the (k + 1) th iteration by utilizing a quasi-Newton method and calculating according to the following formulaN(k+1):
Wherein the content of the first and second substances,tau is an adjustment parameter and takes the value of 10-3And I is an identity matrix, and compared with a traditional quasi-Newton method formula, the above formula can ensure the positive nature of an approximate matrix, thereby ensuring the convergence of the algorithm. SubscriptnNRepresenting the extension variables of boundary nodes and corresponding lines which comprise an area N and an area connected with the area N in the power distribution network;
(3-4) calculating the extended search direction of the area N in the power distribution network by using the approximate matrix H obtained in the step (3-3) according to the following formula:
wherein, gamma is an adjusting parameter with the value of 10-4The extension variables representing the k-th iteration including the region N and the boundary nodes of the regions connected to the region N in the distribution network and the corresponding lines, the active power variables of the generators representing the region N and the boundary nodes of the regions connected with the region N in the power distribution network and the active power variable of the power flow of the corresponding lines,andrespectively representing the dual variables of an area N, an area boundary node connected with the area N and a corresponding line in the power distribution network;
(3-5) Using the extended search Direction obtained in the above step (3-4)The search direction for region N is calculated using the following equation:
acquiring a corresponding local region partial value of the neighbor region extended search direction by utilizing communication, and simultaneously transmitting corresponding information to the neighbor region to realize the reassembly of the search direction, wherein finally the search direction of the N region can be obtained by the following formula:
wherein, wNRepresenting the total number of regions including region N and the regions contiguous to region N,the search direction of the variables in the area N corresponding to the extended search direction of the distribution network area f in the kth iteration is shown as an information interaction diagram of two areas AB in the distribution network in fig. 4, where the two areas AB only need to be exchangedAndrespectively representing a B area boundary variable part contained in the calculation of the A area expansion direction and an A area boundary variable part contained in the calculation of the B area expansion direction, wherein the rest variables are area internal interaction variables and do not need to communicate with the outside;
(3-6) calculating the step length in the iterative solution formula of the KKT conditional equation of the region N in the power distribution network by using the following formulaThe method comprises the following steps:
the step length is calculated by adopting the following fully distributed backtracking search mode:
(3-6-1) respectively calculating first middle step lengths of all areas in the power distribution network, wherein the first middle step length of the area N is recorded as
(3-6-2) exchanging the calculated first intermediate step length between all the areas and all the other connected areas;
(3-6-3) taking region N as an example, region N receives the first intermediate step size of all other connected regions, and calculates the second intermediate step size of region N by using the following formula
(3-6-4) according to the second intermediate step sizeRegion N calculates search step size for region N Using formulasFor search step lengthUpdating and utilizing according to the updated search step lengthUpdatingUsing updatedCalculating q (y)k+1) Up to q (y)k+1)≤0;
(3-6-5) according to the search step sizeUsing formulasFor search step lengthUpdating and utilizing according to the updated search step lengthUpdatingUsing updatedCalculate | | | rγ(yN k +1)||2Up toExchange with connected areaAnd updating the search step size Wherein, alpha and beta are backtracking parameters with the values of 0.01-0.1 and 0.3-0.8; in one embodiment of the invention, the values of α and β are 0.05 and 0.4, respectively;
(3-6-6) search variable Deltay according to the region NN kAnd search step d of region NN kSolution of the above KKT conditional equation yN kUpdating:
yN k+1=yN k+dN k△yN k
(3-7) calculating the residual sum of the region N in the power distribution networkQuantity rfeasAnd proxy gap Setting a threshold value epsilon for convergence of iterative calculation of a regionfeasAccording to a threshold value epsilonfeasDetermining whether the region iteration calculation converges, wherefeasValue of 10-8~10-6In one embodiment of the invention,. epsilon.feasIs taken as value of 10-6
(3-8) based on the total residual rfeasAnd proxy gapSolution y to the above KKT conditional equationN kIs judged if r is the convergence state offeas≤εfeasAnd is andsolution y of the KKT conditional equation for region N in the distribution networkN kConvergence, setting region N convergence flag Conv in distribution networkNGo to step (3-9) if r is 1feasfeasAndone or both of them are satisfied, then the region N convergence flag Conv is assertedNGo to step (3-3-2) when k is 0 and k is k + 1;
(3-9) Conv according to the local convergence statusNAnd (3) performing global convergence state calculation:
the communication exchange local convergence state of the defined region N and the connected region is calculated as a global convergence mark:
(3-10) if the global convergence flag1, completing economic dispatching solution of the multi-region fully-distributed active power distribution network at the time t to obtain planned active power of each conventional generator in the power distribution network at the time tAnd planned active power of the photovoltaic power generation plantWhile setting t +1 to go to (3-11), ifMaking k equal to k +1, turning to the step (3-3-2), and continuing to perform the next iterative computation;
(3-11) if T is equal to T, finishing the calculation, and obtaining the planned active power of each conventional generator in the power distribution network at each time T in the T timeAnd planned active power of the photovoltaic power generation plantCompleting economic dispatch of the multi-region fully-distributed active power distribution network if t<And T, turning to the step (3-2).

Claims (1)

1. A multi-region full-distributed active power distribution network economic dispatching method is characterized by mainly comprising the following steps:
(1) establishing an economic dispatching model of the active power distribution network, which comprises the minimum light abandoning cost:
the objective function F of the economic dispatching model of the active power distribution network is as follows:
wherein T represents the total time of economic dispatching on the active power distribution network, T represents any time in the total time of economic dispatching, i is the ith distributed conventional small generator in the power distribution network, and GconRepresenting the collection of all distributed conventional small generators in the distribution grid,representing a cost function of a conventional distributed small generator in a power distribution network, the cost function being calculated by the formula:wherein a isi、biAnd ciRespectively, the cost parameter, a, of the ith distributed conventional small generator in the power distribution networki、biAnd ciThe values of (A) are respectively 0.01-1, 20-20.5 and 600-1500,the planned active power of the ith distributed conventional small generator in the power distribution network at the time t, j is the jth distributed photovoltaic power generation equipment in the power distribution network, GsolarRepresenting a collection of distributed photovoltaic power plants in a power distribution network, omegajThe light abandonment penalty coefficient of the jth photovoltaic power generation equipment is 500-1000,representing the planned active power of the distributed photovoltaic power generation equipment in the power distribution network at time t,representing a predicted value of the distributed photovoltaic power generation equipment in the power distribution network at the time t;
the constraint conditions of the economic dispatching model of the active power distribution network comprise:
a) the load flow balance constraint condition considering the network loss of the active power distribution network is as follows:
wherein the content of the first and second substances,the active power of the line between the node m and the node h in the power distribution network at the moment t,representing the active power of the line between node n and node m in the distribution network at time t,representing the active power loss of the line between node n and node m in the distribution network at time t, andrespectively representing historical data points of the active power and the reactive power of a line between a node n and a node m in the power distribution network and the voltage of the node n,andr is obtained from historical data of a power distribution network scheduling systemnmRepresenting the line resistance between node n and node m in the distribution network,representing the active power output of node m in the distribution network at time t, representing the active power of the generator at the node m in the distribution network at the moment t,is a variable to be solved for the method,representing the load active power demand of a node m in the power distribution network at the time t;
centralizing variables to be solved in the power flow balance constraint condition considering the network loss of the active power distribution network, and rewriting the variables to be solved into the following matrix form:
Bsx=b
wherein, BsRepresenting an extended adjacency matrix, Bs=[-Z B]Any element Z in the matrix ZopComprises the following steps:
arbitrary element B in matrix BnuComprises the following steps:x represents the set of all variables to be solved, x is a vector,b represents the deviation of the active power of the nodes in the power distribution network, b is a vector,
b) the inequality constraint condition of the network where the active power distribution network is located is considered as follows:
the climbing constraint conditions of the conventional unit in the active power distribution network are as follows:
wherein, mudiAnd muuiRespectively representing the climbing coefficient and the descending coefficient of the ith conventional unit in the power distribution network, wherein delta t is the set economic dispatching time interval of the power distribution network;
the constraint conditions of the active power of the conventional unit in the active power distribution network are as follows:
wherein the content of the first and second substances,andrespectively representing the upper limit of active power and the lower limit of active power of the ith conventional unit in the power distribution network;
the constraint conditions of the line tide capacity in the active power distribution network are as follows:
wherein the content of the first and second substances,andrepresenting the upper limit/the lower limit of active power of a line power flow between a node n and a node m in the power distribution network;
the active power constraint conditions of the photovoltaic units in the active power distribution network are as follows:
writing inequality constraint conditions of the network where the active power distribution network is located into a matrix form as follows:
(2) and (3) converting the economic dispatching model of the active power distribution network into the following KKT condition equation by adopting a primal-dual interior point method:
wherein r ispriIs the original residual of the KKT conditional equation, rdualAnd rcentThe dual residual and the central residual of the KKT conditional equation, respectively, Dq represents the function gradient matrix,λ and ν are the dual multipliers of the above inequalities and equality constraints, respectively, γ is a calculation parameter,mu is a factor, the value of mu is 10, and the value of g is equal to twice of the total number Y of the inequality constraint conditions;
using a numerical iteration method to solve the KKT conditional equation, where y is (x, λ, v), and Δ y is (Δ x, Δ λ, Δ v), then the iterative solution formula of the KKT conditional equation is:
yk+1=yk+dkΔyk
where k is the number of iterations, ykRepresents the value of y at the k iteration, dkRepresents the value of the step length at the kth iteration, Δ ykRepresenting the search direction at the k-th iteration;
obtaining the search direction delta y according to a first-order Taylor expansion formula by using the following formulak
Δyk=-(Drγ(yk))-1rγ(yk)
Wherein the content of the first and second substances, representing a quadratic gradient calculation, superscriptTRepresenting matrix transposition, diag representing an operation on a diagonal matrix, BsRepresents an extended adjacency matrix;
(3) solving y in the KKT conditional equation to realize economic dispatching of the multi-region fully-distributed active power distribution network, firstly defining the connection lines between the regions in the active power distribution network as follows: the method comprises the following steps of performing region allocation on the links between regions in the active power distribution network, wherein the links between any two regions belong to the region with the smaller number of nodes:
(3-1) setting the economic dispatching initial time t to be 0 during initialization;
(3-2) setting the initial iteration step number k of the active power distribution network multi-region fully-distributed economic scheduling calculation at the time t to be 0;
(3-3) solving the approximation matrix H, comprising the steps of:
(3-3-1) constructing an approximate matrix initial value H of a region N in an active power distribution networkN(0):
HN(0)=MN+LN
Wherein M isNIs a block diagonal matrix comprising an area N and other areas connected to the area N, each diagonal block in the block diagonal matrix being defined by the Dr of each area in the distribution networkγ(y0) Namely Drγ(y) initial value composition, LNIs a non-block diagonal symmetric matrix, the elements in the non-block diagonal symmetric matrix are formed by the initial value H of an approximate matrixN(0) The element composition in the B matrix of the corresponding position;
(3-3-2) KKT Condition for region N in active distribution networkDual residual and original residual r in the equationdualAnd rpriThe local dual residual r after correction is obtained by respectively correctingdualNAnd the original residual rpriN
Wherein, subscript N indicates that the element is an element in an area N in the active power distribution network, and vector MrdualNAny element in (1) is MrdualN(pN+n),Region Q is connected to region N, pNIs the number of generators in the N region, vector MrpriNAny element in (1) is MrpriN(n),pQIs the number of generators in the Q region, laNQIs the above BsThe values of the boundary elements at the junction of region N and region Q in the matrix,Nnand QmRespectively representing a boundary node N of an area N and a boundary node m of an area Q connected with the area N in the power distribution network;
(3-3-3) obtaining an area N approximate matrix H in the power distribution network of the (k + 1) th iteration by utilizing a quasi-Newton method and calculating according to the following formulaN(k+1):
Wherein the content of the first and second substances,tau is an adjustment parameter and takes the value of 10-3I is an identity matrix, subscriptnNRepresenting the extension variables of boundary nodes and corresponding lines which comprise an area N and an area connected with the area N in the power distribution network;
(3-4) calculating the extended search direction of the area N in the power distribution network by using the approximate matrix H obtained in the step (3-3) according to the following formula:
wherein, gamma is an adjusting parameter with the value of 10-4The extension variables representing the k-th iteration including the region N and the boundary nodes of the regions connected to the region N in the distribution network and the corresponding lines, the active power variables of the generators representing the region N and the boundary nodes of the regions connected with the region N in the power distribution network and the active power variable of the power flow of the corresponding lines,andrespectively representing the dual variables of an area N, an area boundary node connected with the area N and a corresponding line in the power distribution network;
(3-5) Using the extended search Direction obtained in the above step (3-4)The search direction for region N is calculated using the following equation:
wherein, wNRepresenting the total number of regions including region N and the regions contiguous to region N,the search direction of the variable in the area N corresponding to the expanded search direction of the distribution network area f in the k iteration is obtained;
(3-6) calculating the step length in the iterative solution formula of the KKT conditional equation of the region N in the power distribution network by using the following formulaThe method comprises the following steps:
(3-6-1) respectively calculating first middle step lengths of all areas in the power distribution network, wherein the first middle step length of the area N is recorded as
(3-6-2) exchanging the calculated first intermediate step length between all the areas and all the other connected areas;
(3-6-3) taking region N as an example, region N receives the first intermediate step size of all other connected regions, and calculates the second intermediate step size of region N by using the following formula
(3-6-4) according to the second intermediate step sizeRegion N calculates search step size for region NUsing formulasFor search step lengthUpdating and utilizing according to the updated search step lengthUpdatingUsing updatedCalculating q (y)k+1) Up to q (y)k+1)≤0;
(3-6-5) according to the search step sizeUsing formulasFor search step lengthUpdating and utilizing according to the updated search step lengthUpdatingUsing updatedCalculate | | | rγ(yN k+1)||2Up toExchange with connected areaAnd updating the search step size Wherein, alpha and beta are backtracking parameters with the values of 0.01-0.1 and 0.3-0.8;
(3-6-6) search variable Deltay according to the above region NN kAnd search step d of region NN kSolution of the above KKT conditional equation yN kUpdating:
yN k+1=yN k+dN kΔyN k
(3-7) calculating residual total r of the region N in the power distribution networkfeasAnd proxy gap Setting a threshold value epsilon for convergence of iterative calculation of a regionfeasAccording to a threshold value epsilonfeasDetermining whether the region iteration calculation converges, wherefeasValue of 10-8~10-6
(3-8) based on the total residual rfeasAnd proxy gapSolution y to the above KKT conditional equationN kIs judged if r is the convergence state offeas≤εfeasAnd is andarea in the distribution networkSolution y of KKT conditional equation of NN kConvergence, setting region N convergence flag Conv in distribution networkNGo to step (3-9) if r is 1feas>εfeasAndone or both of them are satisfied, then the region N convergence flag Conv is assertedNGo to step (3-3-2) when k is 0 and k is k + 1;
(3-9) Conv according to the local convergence statusNAnd (3) performing global convergence state calculation:
the communication exchange local convergence state of the defined region N and the connected region is calculated as a global convergence mark:
(3-10) if the global convergence flag1, completing economic dispatching solution of the multi-region fully-distributed active power distribution network at the time t to obtain planned active power of each conventional generator in the power distribution network at the time tAnd planned active power of the photovoltaic power generation plantWhile setting t +1 to go to (3-11), ifMaking k equal to k +1, turning to the step (3-3-2), and continuing to perform the next iterative computation;
(3-11) if T is equal to T, finishing the calculation, and obtaining the planned active power of each conventional generator in the power distribution network at each time T in the T timeAnd planned active power of the photovoltaic power generation plantAnd (4) completing economic dispatching of the multi-region fully-distributed active power distribution network, and if T is less than T, turning to the step (3-2).
CN201710605843.8A 2017-07-24 2017-07-24 Economic dispatching method for multi-region fully-distributed active power distribution network Expired - Fee Related CN107482673B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710605843.8A CN107482673B (en) 2017-07-24 2017-07-24 Economic dispatching method for multi-region fully-distributed active power distribution network

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710605843.8A CN107482673B (en) 2017-07-24 2017-07-24 Economic dispatching method for multi-region fully-distributed active power distribution network

Publications (2)

Publication Number Publication Date
CN107482673A CN107482673A (en) 2017-12-15
CN107482673B true CN107482673B (en) 2019-12-27

Family

ID=60595955

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710605843.8A Expired - Fee Related CN107482673B (en) 2017-07-24 2017-07-24 Economic dispatching method for multi-region fully-distributed active power distribution network

Country Status (1)

Country Link
CN (1) CN107482673B (en)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108711890A (en) * 2018-06-27 2018-10-26 广东电网有限责任公司 Ahead market goes out clearing method, system, device and computer readable storage medium
CN109120011B (en) * 2018-09-29 2019-12-13 清华大学 distributed power distribution network congestion scheduling method considering distributed power sources
CN110266038B (en) * 2019-05-28 2023-10-20 广东电网有限责任公司电力调度控制中心 Distributed coordination regulation and control method for multiple virtual power plants
CN110310014B (en) * 2019-06-05 2022-05-03 清华大学 Ubiquitous power Internet of things distributed economic dispatching method based on transition matrix
CN111525556B (en) * 2020-05-06 2023-03-10 华东交通大学 Multi-target optimal power flow calculation method considering wind power confidence risk
CN112926185B (en) * 2021-01-26 2022-11-22 吉林化工学院 Power economy emission scheduling method based on improved kernel search optimization
CN113472014B (en) * 2021-06-25 2023-09-26 国网山东省电力公司泗水县供电公司 Distribution network optimal scheduling method and system containing distributed power supply
CN115632404B (en) * 2022-12-22 2023-03-10 国网冀北电力有限公司 Distributed voltage control method and device based on digital twin power distribution system

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101699448A (en) * 2009-10-26 2010-04-28 清华大学 Transient stability distributed simulation method of electric power system
CN105119275A (en) * 2015-08-18 2015-12-02 河海大学 An algorithm for electric power system dynamic optimal power flows of a meter and a unified power flow controller

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101699448A (en) * 2009-10-26 2010-04-28 清华大学 Transient stability distributed simulation method of electric power system
CN105119275A (en) * 2015-08-18 2015-12-02 河海大学 An algorithm for electric power system dynamic optimal power flows of a meter and a unified power flow controller

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Optimal economic dispatch for power generation using artificial neural network;S. Panta et al.;《2007 International Power Engineering Conference (IPEC 2007)》;20080502;第1343-1348页 *
直角坐标下含零注入约束的电力系统状态估计修正牛顿法;郭烨等;《中国电机工程学报》;20120705;第32卷(第19期);第96-100、s14页 *

Also Published As

Publication number Publication date
CN107482673A (en) 2017-12-15

Similar Documents

Publication Publication Date Title
CN107482673B (en) Economic dispatching method for multi-region fully-distributed active power distribution network
CN110266038B (en) Distributed coordination regulation and control method for multiple virtual power plants
CN110460036B (en) Distributed optimization method for alternating current-direct current power distribution network considering wind power uncertainty
CN111416356B (en) Transmission and distribution network linkage optimization method based on alternate direction multiplier method and optimal power flow
CN110690702B (en) Active power distribution network optimal scheduling and operation method considering comprehensive bearing capacity
CN106712076A (en) Power transmission system optimization method on offshore wind farm cluster scale
CN107069814A (en) The Fuzzy Chance Constrained Programming method and system that distribution distributed power source capacity is layouted
CN108075487A (en) The hierarchical control method for the isolated island micro-capacitance sensor that adaptive sagging and uniformity is combined
CN104967126A (en) Interbasin hydropower station group multiple power grid combination adjusting peak method facing regional power grid
CN103514374B (en) The discrimination method of infeasible transmission profile constraints in the online rolling scheduling of power system
Hu et al. Hierarchical distributed scheme for demand estimation and power reallocation in a future power grid
CN111553544B (en) Industrial park distributed comprehensive demand response method based on consistency algorithm
CN115809282A (en) Transformer substation carbon emission monitoring method and system
CN116345564A (en) Multi-time-scale distributed collaborative optimization scheduling method and system for comprehensive energy system
CN104659812A (en) Multi-microgrid coordination control method based on predictive control
CN106611966B (en) Multi-inverter type exchanges micro-capacitance sensor distribution economy Automatic Generation Control algorithm
CN113344283A (en) Energy internet new energy consumption capacity assessment method based on edge intelligence
CN110380407B (en) Power distribution network operation optimization method considering agricultural electric irrigation and drainage loads
CN110535123B (en) Rapid analytic distributed multi-target multi-microgrid economic scheduling optimization method
CN112671047A (en) Active power distribution network reconstruction and reactive power joint robust optimization method considering limit scene
CN109494750B (en) Hierarchical distributed voltage optimization control method for high and medium voltage distribution network
CN113673912B (en) Distribution-gas network distributed collaborative planning method and system considering influence of power transmission network
CN116707023A (en) Active power distribution network layering and partitioning comprehensive optimization method based on source-load correlation clustering
CN107370182B (en) Distributed power supply access planning method for active power distribution network ternary planning system
CN115833105A (en) Power distribution network planning method based on cluster division

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20191227

Termination date: 20210724

CF01 Termination of patent right due to non-payment of annual fee