EP1852820A1

EP1852820A1 - Method of automatically optimising a natural gas transport network

Info

Publication number: EP1852820A1
Application number: EP07107550A
Authority: EP
Inventors: Eglantine Peureux; Benoît Casoetto; Tony Pillay; Marie Benoit
Original assignee: Gaz de France SA
Current assignee: Engie SA
Priority date: 2006-05-05
Filing date: 2007-05-04
Publication date: 2007-11-07
Also published as: RU2007116343A; FR2900753B1; US20070260333A1; US7561928B2; CA2587070A1; FR2900753A1

Abstract

The method involves determining values for continuous variables e.g. pressure of natural gas, and discrete variables e.g. startup of compressors. An initial state of optimization of value intervals and value assemblies for the continuous and discrete variables, respectively, is selected. Combinations of values are explored for the variables by constructing a tree with branches that is connected to interconnection points (1.1-1.13) e.g. nodes. The values of required sizes considered as optimal are evaluated, when preset constraints are not transgressed and a preset criterion is minimized.

Description

La présente invention a pour objet un procédé d'optimisation automatique d'un réseau de transport de gaz naturel en régime permanent, le réseau de transport de gaz naturel comprenant à la fois un ensemble d'ouvrages passifs tels que des canalisations ou des résistances, et un ensemble d'ouvrages actifs comprenant des vannes de régulation, des vannes d'isolement, des stations de compression avec chacune au moins un compresseur, des dispositifs de stockage ou d'alimentation, des dispositifs de consommation, des éléments de dérivation des stations de compression et des éléments de dérivation des vannes de régulation, les ouvrages passifs et les ouvrages actifs étant reliés entre eux par des jonctions, le procédé d'optimisation comprenant la détermination de valeurs pour des variables continues telles que la pression et le débit du gaz naturel en tout point du réseau de transport, et la détermination de valeurs pour des variables discrètes telles que l'état de démarrage des compresseurs, l'état d'ouverture des stations de compression, l'état d'ouverture des vannes de régulation, l'état des éléments de dérivation des stations de compression, l'état des éléments de dérivation des vannes de régulation, l'orientation des stations de compression et l'orientation des vannes de régulation.The present invention relates to a method of automatic optimization of a natural gas transmission network in steady state, the natural gas transmission network comprising both a set of passive structures such as pipes or resistors, and a set of active structures comprising control valves, isolation valves, compression stations with each at least one compressor, storage or supply devices, consumption devices, station diversion elements compressors and diverting elements of the regulating valves, the passive structures and the active structures being interconnected by junctions, the optimization method comprising the determination of values for continuous variables such as the pressure and the flow rate of the gas in all points of the transport network, and the determination of values for discrete variables such as compressor run-up, the opening state of the compressor stations, the opening state of the control valves, the state of the bypass elements of the compressor stations, the state of the bypass elements of the control valves, the orientation of the compressor stations and the orientation of the control valves.

La présente invention vise à permettre de déterminer notamment les valeurs optimales de pression et de débit en tout point d'un réseau de transport de gaz naturel en régime permanent. L'invention vise également à permettre de déterminer de façon optimale et automatique non seulement des variables continues, telles que le débit, qui peuvent prendre toutes les valeurs comprises dans un intervalle, mais également des variables discrètes qui ne peuvent prendre qu'un nombre fini de valeurs.The present invention aims to determine in particular the optimum values of pressure and flow at any point of a natural gas transmission network in steady state. The invention also aims to allow optimal and automatic determination of not only continuous variables, such as the flow rate, which can take all the values included in an interval, but also discrete variables which can take only a finite number of values.

A titre d'exemple, l'ouverture d'une vanne est une variable discrète, dès lors que cette vanne ne peut être qu'ouverte (ce qui peut être représenté par exemple par un 1) ou fermée (ce qui peut alors être représenté par un 0).By way of example, the opening of a valve is a discrete variable, since this valve can only be open (which can be represented for example by a 1) or closed (which can then be represented by a 0).

Le procédé selon l'invention vise ainsi à permettre de déterminer de façon automatique et de manière optimale notamment des facteurs tels que l'ouverture des vannes, le démarrage des compresseurs, l'orientation des ouvrages actifs (station de compression et vannes de régulation), l'état des éléments de dérivation (by-pass) de ces ouvrages actifs, voire l'adaptation série ou parallèle de certains compresseurs.The method according to the invention thus aims to allow to determine automatically and optimally including factors such as the opening of the valves, the start of the compressors, the orientation of the active structures (compressor station and control valves) , the state of the bypass elements of these active structures, or the serial or parallel adaptation of some compressors.

Pour déterminer par le calcul les caractéristiques d'un réseau de transport de gaz, quelle que soit la modélisation physique retenue, on prend en compte traditionnellement la loi des noeuds et la loi des mailles (encore dénommées lois de Kirchhoff du fait de leur emprunt à la théorie des circuits électriques).To determine by calculation the characteristics of a gas transport network, whatever the physical model chosen, we take into account traditionally the law of the nodes and the law of the meshes (still called laws of Kirchhoff because of their borrowing at the theory of electrical circuits).

Un réseau de transport de gaz peut être représenté sous la forme d'un graphe composé de noeuds (sommets) et d'arcs qui établissent une relation orientée entre deux noeuds. Les arcs possèdent un attribut « ETAT » qui indique si l'arc est activé ou désactivé.A gas transport network can be represented as a graph consisting of nodes (vertices) and arcs that establish a directed relationship between two nodes. Arcs have a "STATUS" attribute that indicates whether the arc is enabled or disabled.

Selon la loi des noeuds, il y a pour tous les noeuds du réseau égalité entre la quantité de gaz entrant dans un noeud et la quantité de gaz sortant de ce noeud et globalement tout ce qui entre dans le réseau doit en sortir.According to the law of the nodes, there is for all the nodes of the network equality between the quantity of gas entering a node and the quantity of gas leaving this node and globally everything which enters the network must leave it.

De manière synthétique, selon la loi des noeuds, on obtient le système d'équations linéaires suivant : $B . E_{arcs} = E_{consommations} + E_{ressources} + C_{stations}$

avec

B : matrice d'incidence du réseau traduisant la correspondance entre les arcs et les noeuds du réseau,
E_arcs : vecteur des quantités s'écoulant en chaque arc,
E_{consommations} : vecteur des quantités livrées aux consommations,
E_ressources : vecteur des quantités émises ou injectées aux ressources (dispositifs de stockage ou d'alimentation),
C_station : vecteur des quantités de gaz carburant consommées par les stations de compression.

In a synthetic way, according to the law of the nodes, one obtains the system of linear equations following:

B . E_{bows} = E_{consumption} + E_{resources} + {VS}_{stations}

with

B: matrix of incidence of the network translating the correspondence between the arcs and the nodes of the network,
E _arcs : vector of quantities flowing in each arc,
E _consumption : vector of quantities delivered to consumptions,
E _resources : vector of quantities emitted or injected into resources (storage or feeding devices),
C _station : vector of the quantities of fuel gas consumed by the compressor stations.

La loi des noeuds permet de définir ainsi un système d'équations linéaires.The law of the nodes makes it possible to define thus a system of linear equations.

Toutes les quantités entrant ou sortant sont algébriques et leur signe est défini par le choix d'une convention. On peut considérer que ce qui entre dans un noeud est positif tandis que ce qui en sort est négatif.All quantities entering or leaving are algebraic and their sign is defined by the choice of a convention. We can consider that what enters a node is positive while what comes out is negative.

Selon la loi des mailles, la somme algébrique, le long d'une maille, des différences de pression du gaz entre deux noeuds consécutifs est nulle. La loi des mailles permet ainsi de définir un système d'équations : $\sum_{maille} ΔP = 0$

avec ΔP : différence des pressions entre deux noeuds consécutifs d'une maille.According to the law of the meshes, the algebraic sum, along a mesh, of the differences of pressure of the gas between two consecutive nodes is null. The law of meshes thus makes it possible to define a system of equations:

\underset{stitch}{Σ} .DELTA.P = 0

with ΔP: difference in pressures between two consecutive knots of a mesh.

Parce que l'on connaît les formules de perte de charge dans les canalisations selon la forme suivante : $P_{1}^{2} - P_{2}^{2} = α \times Q \times |Q|,$

on peut exprimer également de façon équivalente la loi des mailles à l'aide de différences de pression au carré :

\sum_{maille} Δ P^{2} = 0

avec ΔP² : différence des pressions quadratiques entre deux noeuds consécutifs d'une maille.Because we know the pressure loss formulas in the pipes according to the following form:

P_{1}^{2} - P_{2}^{2} = α \times Q \times |Q|,

the mesh law can also be expressed equally by squared pressure differences:

\underset{stitch}{Σ} Δ P^{2} = 0

with ΔP ² : difference of the quadratic pressures between two consecutive nodes of a mesh.

La loi des mailles permet ainsi de définir un système d'équations non linéaires.The law of meshes thus makes it possible to define a system of nonlinear equations.

On connaît déjà des procédés de calcul de réseau qui abordent le problème en supposant celui-ci parfaitement déterminé, c'est-à-dire en supposant que le nombre d'inconnues est égal au nombre d'équations.Network computation methods are already known which approach the problem by assuming it perfectly determined, that is to say assuming that the number of unknowns is equal to the number of equations.

Si l'on considère un réseau de N noeuds et M mailles, on en déduit le nombre d'arcs égal à N + M - 1, auxquels correspondent autant d'équations indépendantes, à savoir N - 1 équations selon la loi des noeuds et M équations selon la loi des mailles.If we consider a network of N nodes and M meshes, we deduce the number of arcs equal to N + M - 1, to which correspond as many independent equations, namely N - 1 equations according to the law of the nodes and M equations according to the law of meshes.

Les deux lois de Kirchhoff permettent de déterminer des débits (dans la mesure où la loi des mailles remplace les différences de pressions quadratiques par leur expression équivalente fonction du débit, en général de la forme ΔP² = α × Q² où α est considéré constant.The two laws of Kirchhoff make it possible to determine flow rates (insofar as the law of meshes replaces the differences of quadratic pressures by their equivalent expression function of the flow, in general of the form ΔP ² = α × Q ² where α is considered constant .

Lorsque le système d'équations de ces deux lois est résolu, les débits sont partout connus et l'imposition d'une pression particulière en un noeud quelconque du réseau permet de connaître les pressions en tous les noeuds.When the system of equations of these two laws is solved, the flows are everywhere known and the imposition of a particular pressure in any node of the network makes it possible to know the pressures in all the nodes.

De manière traditionnelle, les procédés de simulation visant à déterminer les variables continues en tout point d'un réseau comprennent une première phase de résolution des deux lois de Kirchhoff et d'obtention des débits partout et une deuxième phase d'imposition d'une pression en un noeud particulier et d'obtention des pressions partout.In a traditional way, the simulation processes for determining continuous variables at any point in a network include a first phase of solving both Kirchhoff laws and obtaining flow rates everywhere and a second phase of imposing a pressure. in a particular knot and getting pressures everywhere.

Généralement, le processus itère plusieurs fois entre la phase n° 1 et la phase n°2 car les coefficients α intervenant dans les relations de la loi des mailles ne sont pas parfaitement constants et dépendent très légèrement des pressions et des débits.Generally, the process iterates several times between the phase n ° 1 and the phase n ° 2 because the coefficients α intervening in the relations of the law of the meshes are not perfectly constant and depend very slightly on the pressures and flow rates.

Cette démarche impose deux fortes restrictions. La première restriction est qu'elle ne s'applique que sur des réseaux ne comportant que des canalisations ou, de façon plus générale, des ouvrages passifs. En effet, les ouvrages passifs présentent une relation entre la différence des pressions en amont et en aval de l'ouvrage et son débit. Cette relation est l'équation de perte de charge proprement dite. Disposant de cette relation, il est toujours possible de remplacer les différences de pressions par leur expression fonction du débit. En revanche, un ouvrage actif, tel qu'une vanne de régulation ou une station de compression, ne présente pas nécessairement une telle relation ou du moins, si cette équation existe, elle contient au moins une inconnue supplémentaire.This approach imposes two strong restrictions. The first restriction is that it only applies to networks with only pipelines or, more generally, passive works. In fact, passive structures have a relationship between the difference in pressures upstream and downstream of the structure and its flow. This relationship is the actual pressure loss equation. Having this relationship, it is always possible to replace pressure differences by their flow-dependent expression. On the other hand, an active structure, such as a control valve or a compressor station, does not necessarily have such a relationship or at least, if this equation exists, it contains at least one additional unknown.

Les ouvrages actifs constituent des organes de commande du réseau en introduisant des inconnues supplémentaires telles que, par exemple, le degré d'ouverture d'une vanne de régulation. Connaissant le degré d'ouverture et considérant un certain nombre de coefficients caractéristiques fournis par le constructeur, on relie les pressions en amont, en aval et le débit à ce pourcentage d'ouverture.The active structures constitute control members of the network by introducing additional unknowns such as, for example, the degree of opening of a control valve. Knowing the degree of opening and considering a number of characteristic coefficients provided by the manufacturer, the pressures upstream, downstream and the flow rate are related to this percentage of opening.

Pour les stations de compression, l'inconnue introduite est la puissance motrice de compression (puissance dépensée pour la compression) puisque cette dernière est reliée au débit et au taux de compression (rapport entre sa pression en aval et sa pression en amont).For compression stations, the unknown introduced is the driving power of compression (power spent for compression) since the latter is related to the flow rate and the compression ratio (ratio between its downstream pressure and its upstream pressure).

Généralement, les procédés de calcul de réseau permettant la simulation d'ouvrages actifs imposent à l'utilisateur de fixer lui-même la valeur de ces inconnues. Implicitement, les ouvrages actifs ne le sont alors plus car ils présentent une véritable équation de perte de charge (ou de gain s'il s'agit de compression). Typiquement, le contournement proposé par ces procédés consiste à demander à l'utilisateur d'imposer soit la puissance de compression s'il s'agit d'une station de compression, soit le degré d'ouverture de la vanne s'il s'agit d'une détente, etc. L'imposition de ces grandeurs établit un lien entre le débit de l'ouvrage et ses pressions amont et aval. Disposant alors d'une telle relation, il est donc possible de résoudre la seconde loi de Kirchhoff.Generally, network calculation methods for simulating active structures require the user to set the value of these unknowns himself. Implicitly, active works are no longer so because they present a real loss of load equation (or gain if it is compression). Typically, the workaround proposed by these methods consists in asking the user to impose either the compression power if it is a compression station, or the degree of opening of the valve if it is acts of a relaxation, etc. The imposition of these quantities establishes a link between the flow of the structure and its upstream and downstream pressures. Having then such a relation, it is thus possible to solve the second law of Kirchhoff.

Toute la difficulté consiste à déterminer quelle puissance des stations de compression ou quel degré d'ouverture des vannes de régulation imposer. Il n'est pas toujours possible, du moins en un temps raisonnable, de trouver manuellement selon une approche essais/erreurs un jeu de valeurs convenables notamment pour un réseau complexe où les mailles sont interconnectées les unes avec les autres.The difficulty lies in determining the power of the compressor stations or the degree of opening of the control valves. It is not always possible, at least in a reasonable time, to find manually, according to a trial / error approach, a set of suitable values, in particular for a complex network where the meshes are interconnected with each other.

La seconde restriction est la nécessité d'imposer une pression en un noeud particulier de la phase n° 2. Du fait de la première restriction, le réseau est supposé composé uniquement d'ouvrages passifs. En imposant cette pression particulière et après avoir résolu les deux lois de Kirchhoff, les pressions peuvent être connues partout.The second restriction is the need to impose pressure on a particular node of phase 2. Due to the first restriction, the network is assumed to consist solely of passive works. By imposing this particular pressure and after having solved Kirchhoff's two laws, the pressures can be known everywhere.

Si le réseau ne comporte qu'une unique source, il paraît naturel d'imposer la pression au noeud particulier qu'est le noeud de cette source. En général, on impose la pression la plus haute possible en ce point et l'ensemble des pressions en tous les noeuds constitue alors le régime de pression maximal. Une autre approche est de choisir au noeud source une pression la plus basse possible dans la limite où les pressions en tous les noeuds ne soient pas inférieures à un seuil fixé. L'ensemble des pressions en tous les noeuds constitue alors le régime de pression minimal.If the network has only one source, it seems natural to impose pressure on the particular node that is the node of this source. In general, it imposes the highest pressure possible at this point and the set of pressures in all nodes is then the maximum pressure regime. Another approach is to choose at the source node a pressure as low as possible in the limit where the pressures in all the nodes are not lower than a fixed threshold. The set of pressures at all nodes then constitutes the minimum pressure regime.

Si le régime de pression minimal présente des pressions supérieures au régime de pression maximal, cela signifie que l'on ne peut trouver une pression au noeud source qui, à la fois :

soit inférieure à la pression maximale de ce noeud,
soit supérieure à une valeur limite qui permet de satisfaire tous les seuils de pression minimale en tous les noeuds.

If the minimum pressure regime has pressures greater than the maximum pressure regime, it means that a pressure can not be found at the source node which

is less than the maximum pressure of this node,
is greater than a limit value which makes it possible to satisfy all the thresholds of minimum pressure in all the nodes.

Le réseau est dit saturé.The network is said to be saturated.

Dans le cas d'un réseau ne comportant qu'une seule source, le débit de celle-ci injecté sur le réseau est parfaitement déterminé par la loi des noeuds. Ce n'est plus le cas si une seconde source est présente sur le réseau car le nombre de noeuds n'est pas modifié (donc pas d'équation supplémentaire) et une inconnue additionnelle correspondant au débit de cette seconde source est introduite dans le problème.In the case of a network having only one source, the flow rate thereof injected into the network is perfectly determined by the law of the nodes. This is no longer the case if a second source is present on the network because the number of nodes is not modified (thus no additional equation) and an additional unknown corresponding to the flow of this second source is introduced in the problem. .

Les procédés traditionnels de calcul de réseau contournent ce cas en créant une maille fictive entre les deux sources. Cette maille est qualifiée de fictive car l'on suppose que les deux sources sont reliées par une canalisation de longueur nulle et de diamètre très grand. L'introduction de cette nouvelle maille pose, dans le système d'équations, l'équation supplémentaire manquante. L'équilibre entre le nombre d'inconnues du problème et le nombre d'équations est alors rétabli. En général, la canalisation fictive est construite de telle façon qu'elle présente une loi de perte de charge constante (ΔP² = C^ste). Avec cette canalisation fictive, l'imposition d'une pression en une seule des deux sources du réseau permet de connaître les pressions en tous les noeuds du réseau.Traditional network computing methods bypass this case by creating a fictitious mesh between the two sources. This mesh is called fictional because it is assumed that the two sources are connected by a pipeline of zero length and very large diameter. The introduction of this new mesh poses, in the system of equations, the missing additional equation. The balance between the number of unknowns problem and the number of equations is then restored. In general, the dummy duct is constructed in such a way that it has a constant law of pressure loss (ΔP ² = C ^ste ). With this fictitious pipeline, the imposition of a pressure in only one of the two sources of the network makes it possible to know the pressures in all the nodes of the network.

Cette méthode présente toutefois l'inconvénient que si la constante de la loi de perte de charge de la canalisation fictive est trop importante, alors la résolution de la seconde loi de Kirchhoff conduit à trouver un débit qui sort du réseau pour la seconde source, ce qui peut ne pas être souhaitable puisqu'il s'agit d'une source, autrement dit d'une entrée de gaz sur le réseau.However, this method has the disadvantage that if the constant of the law of loss of load of the fictitious pipeline is too important, then the resolution of the second law of Kirchhoff leads to find a flow which leaves the network for the second source, this which may not be desirable since it is a source, in other words a gas inlet on the network.

L'approche du calcul de réseau dans toute sa généralité par simulation n'est donc pas satisfaisante puisque la recherche des valeurs optimales des puissances et des pressions et des débits à imposer doit être prise en charge manuellement.The approach of network calculation in all its generality by simulation is therefore not satisfactory since the search for the optimal values of the powers and the pressures and the flows to be imposed must be taken care of manually.

Pour remédier à ces inconvénients, on a déjà proposé de disposer d'un nombre d'inconnues supérieur au nombre d'équations, de sorte qu'il existe plusieurs solutions au problème posé et que le choix d'une solution particulière va s'opérer selon un critère donné, qui détermine une optimisation.To overcome these drawbacks, it has already been proposed to have a number of unknowns greater than the number of equations, so that there are several solutions to the problem and that the choice of a particular solution will take place. according to a given criterion, which determines an optimization.

Certains procédés connus sont toutefois conçus pour le calcul de réseaux en régime dynamique plutôt qu'en régime permanent.Some known methods, however, are designed for calculating networks in dynamic rather than steady state.

D'autres procédés d'optimisation pour le calcul de réseaux en régime permanent ou dynamique imposent des conditions et contraintes particulières qui rendent ces procédés incomplets ou peu souples.Other optimization methods for calculating steady-state or dynamic networks impose particular conditions and constraints that make these processes incomplete or not very flexible.

La présente invention vise à remédier aux inconvénients précités et à permettre de déterminer automatiquement de façon optimale tous les degrés de liberté d'un réseau de transport de gaz en régime permanent, avec minimisation d'un critère économique et non transgression des contraintes, ou transgression des contraintes a minima.The present invention aims to remedy the above-mentioned drawbacks and to make it possible to automatically determine optimally all the degrees of freedom of a gas transport network in steady state, with minimization of an economic criterion and non-transgression of the constraints, or transgression minimal constraints.

L'invention vise plus particulièrement à effectuer une hybridation d'une méthode d'optimisation combinatoire et continue afin de déterminer les valeurs de l'ensemble des variables discrètes et continues, de façon entièrement automatique.The object of the invention is more particularly to perform a hybridization of a combinatorial and continuous optimization method in order to determine the values of the set of discrete and continuous variables, in a fully automatic manner.

Ces buts sont atteints, conformément à l'invention, grâce à un procédé d'optimisation automatique d'un réseau de transport de gaz naturel en régime permanent, le réseau de transport de gaz naturel comprenant à la fois un ensemble d'ouvrages passifs tels que des canalisations ou des résistances, et un ensemble d'ouvrages actifs comprenant des vannes de régulation, des vannes d'isolement, des stations de compression avec chacune au moins un compresseur, des dispositifs de stockage ou d'alimentation, des dispositifs de consommation, des éléments de dérivation des stations de compression et des éléments de dérivation des vannes de régulation, les ouvrages passifs et les ouvrages actifs étant reliés entre eux par des jonctions, le procédé d'optimisation comprenant la détermination de valeurs pour des variables continues telles que la pression et le débit du gaz naturel en tout point du réseau de transport, et la détermination de valeurs pour des variables discrètes telles que l'état de démarrage des compresseurs, l'état d'ouverture des stations de compression, l'état d'ouverture des vannes de régulation, l'état des éléments de dérivation des stations de compression, l'état des éléments de dérivation des vannes de régulation, l'orientation des stations de compression et l'orientation des vannes de régulation,
caractérisé en ce que l'on choisit comme état initial de l'optimisation des intervalles de valeurs pour les variables continues et des ensembles de valeurs pour les variables discrètes, en ce que l'on explore les possibilités de valeurs pour les variables en construisant au fur et à mesure un arbre avec des branches reliées à des noeuds décrivant les combinaisons de valeurs envisagées en utilisant une technique de séparation des variables, c'est-à-dire de découpage conduisant à la génération de nouveaux noeuds dans l'arbre, et d'évaluation, c'est-à-dire de détermination avec une probabilité forte des branches de l'arbre qui peuvent conduire à des feuilles constituant une solution finale optimisée, de manière à parcourir en priorité ces branches à plus forte probabilité de réussite, les valeurs des grandeurs recherchées étant considérées comme optimales lorsque des contraintes prédéterminées ne sont plus transgressées ou sont transgressées au minimum et une fonction objectif est minimisée, cette fonction objectif étant de la forme $g = α \times Régime + β \times Energie + γ \times Cible$

avec : α, β et γ sont des coefficients de pondération.These objects are achieved, in accordance with the invention, by means of a method of automatic optimization of a steady-state natural gas transmission network, the natural gas transmission network comprising both a set of passive structures such as pipelines or resistors, and a set of active structures comprising control valves, isolation valves, compressor stations with each at least one compressor, storage or supply devices, consumer devices , bypass elements of the compressor stations and bypass elements of the control valves, the passive structures and the active structures being interconnected by junctions, the optimization method comprising the determination of values for continuous variables such as the pressure and flow of natural gas at any point in the transmission system, and the determination of values for discrete variables es such as the start-up status of compressors, the state of opening of the compressor stations, the state of opening of the control valves, the state of the bypass elements of the compressor stations, the state of the bypass elements of the control valves, the orientation of the compressor stations and the orientation of the control valves,
characterized in that value ranges for continuous variables and sets of values for discrete variables are chosen as the initial state of optimization, in that value possibilities for the variables are explored by constructing as a tree with branches connected to nodes describing the combinations of values envisaged using a technique of separation of the variables, that is to say of cutting leading to the generation of new nodes in the tree, and evaluation, that is to say of determination with a strong probability of the branches of the tree which can lead to leaves constituting an optimized final solution, so as to go in priority to those branches with a higher probability of success, the values of the quantities sought are considered optimal when predetermined constraints are no longer transgressed or are transgressed to a minimum and an objective function is minimized, this objective function being of the form

boy Wut = α \times Diet + β \times Energy + γ \times Target

with: α, β and γ are weighting coefficients.

Régime représente un facteur de minimisation ou de maximisation de la pression en des points donnés du réseau tels que tout point aval d'un dispositif de stockage ou d'alimentation, tout point amont et tout point aval d'une station de compression ou d'une vanne de régulation, et tout point amont d'un dispositif de consommation,Regime represents a factor of minimization or maximization of the pressure at given points of the network such as any downstream point of a storage or supply device, any upstream point and any downstream point of a compressor station or a control valve, and any upstream point of a consumption device,

Energie représente un facteur de minimisation de la consommation d'énergie de compression,Energy represents a factor of minimization of the energy consumption of compression,

Cible représente un facteur de maximisation ou de minimisation du débit d'un tronçon du réseau situé entre deux jonctions ou de la pression d'une jonction particulière, et lesdites contraintes prédéterminées comprenant d'une part des contraintes d'égalité comprenant la loi de perte de charge dans les canalisations et la loi des noeuds régissant le calcul des réseaux, et d'autre part des contraintes d'inégalités comprenant des contraintes de débit minimal et maximal, des contraintes de pression minimale et maximale des ouvrages actifs ou passifs, des contraintes de puissance de compression des stations de compression.Target represents a factor for maximizing or minimizing the bit rate of a section of the network located between two junctions or the pressure of a particular junction, and said predetermined constraints comprising on the one hand equality constraints comprising the law of loss. in the pipes and the law of the nodes governing the calculation of the networks, and on the other hand inequality constraints including minimum and maximum flow constraints, minimum and maximum pressure constraints of the active or passive structures, constraints compression power compression stations.

Plus généralement, le problème de configuration optimale des ouvrages actifs est modélisé sous la forme d'un programme P₁ d'optimisation se présentant sous la forme suivante : $P_{1} {\begin{cases} \min_{(x s e)} f (x s) = g (x) + α \times {‖ s ‖}^{2} \\ C_{I} (x) + β . e \leq s_{I} \\ C_{E} (x) = s_{E} \\ x \in R^{n}, s_{I} \in R^{p}, s_{E} \in R^{q}, e \in {\{0 1\}}^{p} \end{cases}$

avec :

x est l'ensemble des variables des débits Q et des pressions P,
g(x) est la fonction objectif constituant le critère économique d'optimisation,
C_l(x) est l'ensemble des p contraintes d'inégalité linéaires et non linéaires sur les ouvrages actifs,
β est une matrice dont les coefficients sont nuls ou égaux aux valeurs maximales des contraintes,
cohérente mais le nombre de variables binaires est réellement : 3 x le nombre d'ouvrages actifs,
C_E(x) est l'ensemble des q contraintes d'égalité linéaires ou non linéaires,
s est une variable d'écart qui, lorsqu'elle est non nulle, représente la transgression d'une contrainte,
α est un coefficient représentant le degré d'autorisation de transgression de contraintes.

More generally, the problem of optimal configuration of active structures is modeled in the form of a program P ₁ optimization in the following form:

P_{1} {\begin{cases} \min_{(x s e)} f (x s) = boy Wut (x) + α \times {‖ s ‖}^{2} \\ {VS}_{I} (x) + β . e \leq s_{I} \\ {VS}_{E} (x) = s_{E} \\ x \in R^{not}, s_{I} \in R^{p}, s_{E} \in R^{q}, e \in {\{0 1\}}^{p} \end{cases}

with:

x is the set of variables of flow rates Q and pressures P,
g (x) is the objective function constituting the economic criterion of optimization,
C _l (x) is the set of p linear and nonlinear inequality constraints on active structures,
β is a matrix whose coefficients are equal to or less than the maximum values of the stresses,
consistent but the number of binary variables is really: 3 x the number of active works,
C _E (x) is the set of q linear or nonlinear equality constraints,
s is a variance variable which, when nonzero, represents the transgression of a constraint,
α is a coefficient representing the degree of authorization of constraint transgression.

Selon un mode particulier de réalisation, les variables sont représentées par des intervalles, la technique de séparation des variables est appliquée aux seules variables discrètes et des bornes de la fonction objectif sont calculées en utilisant l'arithmétique des intervalles.According to a particular embodiment, the variables are represented by intervals, the variable separation technique is applied to only the discrete variables and the limits of the objective function are calculated using the arithmetic of the intervals.

Selon un autre mode particulier de réalisation, les variables sont représentées par des intervalles, la technique de séparation des variables est appliquée à la fois aux variables discrètes et aux variables continues, la séparation comprenant le découpage de l'espace de définition des variables continues, l'exploration s'effectuant séparément sur des parties de l'ensemble réalisable et l'intervalle de variation de la fonction objectif étant évalué sur chacune de ces parties.According to another particular embodiment, the variables are represented by intervals, the technique of separation of the variables is applied to both the discrete variables and the continuous variables, the separation comprising the division of the definition space of the continuous variables, the exploration being performed separately on parts of the feasible set and the range of variation of the objective function being evaluated on each of these parts.

Dans ce cas, avantageusement, lors de l'exploration des possibilités de valeurs pour les variables avec une technique de séparation des variables et d'évaluation, on établit d'abord une liste de noeuds à explorer triés selon un critère de mérite M calculé pour chaque noeud, tant que la liste de noeuds à explorer n'est pas vide, pour chaque noeud courant, on évalue si ce noeud courant peut contenir une solution, si c'est le cas, on découpe l'intervalle correspondant à la variable considérée selon une loi de séparation pour établir une liste de noeuds fils, pour chaque noeud fils on évalue des bornes minimale et maximale de la fonction objectif et on évalue si le noeud fils peut améliorer la situation courante, si c'est le cas on effectue une propagation de la contrainte sur ses variables, si la propagation ne conduit pas à des intervalles vides, on évalue des bornes minimale et maximale de la fonction objectif et on vérifie qu'il n'y a pas d'impossibilité à ce que le noeud fils contienne au moins une solution faisable, on effectue un test pour déterminer s'il y a encore des valeurs discrètes non instanciées, c'est-à-dire des variables pour lesquelles aucune valeur précise et définitive n'a pu être décidée, on met à jour la meilleure solution courante s'il y a lieu et on calcule le mérite du noeud pour l'insérer dans la liste des feuilles, triée selon ce critère de mérite.In this case, advantageously, during the exploration of the possibilities of values for the variables with a technique of separation of the variables and evaluation, we first establish a list of nodes to explore sorted according to a merit criterion M calculated for each node, as long as the list of nodes to explore is not empty, for each current node, it is evaluated whether this current node can contain a solution, if it is the case, we cut the interval corresponding to the variable considered according to a separation law for establishing a list of child nodes, for each child node, minimum and maximum limits of the objective function are evaluated and it is judged whether the child node can improve the current situation, if it is the case a propagation of the constraint on its variables, if the propagation does not lead to empty intervals, one evaluates the minimum and maximum limits of the objective function and verifies that there is no impossibility that that the child node contains at least one feasible solution, a test is performed to determine whether there are still uninstantiated discrete values, i.e., variables for which no definite and definitive value could have been decided, we update the best current solution if necessary and calculate the merit of the node to insert it in the list of sheets, sorted according to this criterion of merit.

Le procédé selon l'invention peut en particulier mettre en oeuvre les caractéristiques avantageuses suivantes :The method according to the invention can in particular implement the following advantageous characteristics:

Le critère de mérite M est tel qu'un noeud est exploré en priorité lorsqu'il présente la plus faible borne minimale de la fonction objectif.The merit criterion M is such that a node is first explored when it has the lowest minimum bound of the objective function.

Lors des tests d'élimination des noeuds ne pouvant pas contenir l'optimum, on met en oeuvre l'une des méthodes consistant à utiliser la monotonie de la fonction objectif, à utiliser un test de contraintes transgressées ou à utiliser un test de valeur d'objectif moins bonne que la valeur courante.During the node elimination tests that can not contain the optimum, one of the methods is to use the monotony of the objective function, to use a test of transgressed constraints or to use a test of value. objective less good than the current value.

Lors de la séparation d'un noeud courant en noeuds fils, on divise le domaine de variation d'une ou plusieurs variables choisies selon des critères basés sur le diamètre d'intervalles rattachés aux variables.When separating a current node into child nodes, the variation domain of one or more variables chosen is divided according to criteria based on the interval diameter attached to the variables.

Le procédé comprend en outre un critère d'arrêt basé sur le temps d'exécution ou sur l'évaluation de certains diamètres d'intervalles.The method further comprises a stopping criterion based on the execution time or the evaluation of certain interval diameters.

En complément de la propagation des contraintes, on procède à une actualisation de la borne maximale de l'optimum de la fonction objectif en utilisant les conditions d'optimalité du problème d'optimisation dites de Fritz-John.In addition to the propagation of the constraints, the maximum bound of the optimum of the objective function is updated by using the optimality conditions of the so-called Fritz-John optimization problem.

Lorsqu'à un noeud du procédé de séparation et d'évaluation, toutes les variables discrètes ont été instanciées, on met en oeuvre en outre un processus d'optimisation non linéaire basé sur une méthode de points intérieurs.When at a node of the separation and evaluation method, all the discrete variables have been instantiated, a nonlinear optimization process based on an internal point method is also implemented.

De façon alternative, à chaque noeud du procédé de séparation et d'évaluation, on met en outre en oeuvre un processus d'optimisation non linéaire basé sur une méthode de points intérieurs.Alternatively, at each node of the separation and evaluation method, a nonlinear optimization process based on an internal point method is also implemented.

D'autres caractéristiques et avantages de l'invention ressortiront de la description suivante de modes particuliers de réalisation, donnés à titre d'exemples, en référence aux dessins annexés, sur lesquels :

la Figure 1 est un schéma-bloc montrant les modules principaux d'un système d'optimisation automatique d'un réseau de transport de gaz selon l'invention ;
la Figure 2 est une vue schématique d'un exemple de partie de réseau de transport de gaz ;
la Figure 3 est une vue schématique d'un exemple de configuration d'une station de compression se situant sur un point d'interconnexion d'un réseau de transport de gaz ;
la Figure 4 est une vue schématique montrant le processus d'exploration d'un arbre selon la technique de séparation des variables et d'évaluation ;
la Figure 5 est une vue schématique d'un exemple de partie de réseau à laquelle est appliqué le procédé d'optimisation selon l'invention,
la Figure 6 est un tableau donnant des exemples d'intervalles de pression d'initialisation pour différents noeuds de la partie de réseau de la figure 5 ;
la Figure 7 est un tableau donnant des exemples d'intervalles de débit d'initialisation pour différents arcs de la partie de réseau de la Figure 5 ;
la Figure 8 est un tableau donnant les résultats de tests effectués sur la partie de réseau de la figure 5 ;
la Figure 9 est un tableau donnant des résultats des intervalles de pression pour les différents noeuds de la partie du réseau de la Figure 5 dans les cas du tableau de la Figure 8 où la propagation n'est pas stoppée ;
la Figure 10 est un tableau donnant les résultats des intervalles de débit pour les différents arcs de la partie du réseau de la Figure 5 dans les cas du tableau de la Figure 8 où la propagation n'est pas stoppée ;
la Figure 11 est un organigramme illustrant un exemple de mise en oeuvre du procédé d'optimisation selon l'invention;
la Figure 12 est un diagramme montrant un arbre de calcul qui représente la propagation/rétro-propagation de contraintes ; et
la Figure 13 est une vue schématique d'un exemple de réseau de transport de gaz naturel auquel est applicable l'invention.

Other features and advantages of the invention will emerge from the following description of particular embodiments, given by way of example, with reference to the appended drawings, in which:

Figure 1 is a block diagram showing the main modules of an automatic optimization system of a gas transport network according to the invention;
Figure 2 is a schematic view of an exemplary gas transport network portion;
Figure 3 is a schematic view of an exemplary configuration of a compressor station located at an interconnection point of a gas transport network;
Figure 4 is a schematic view showing the process of exploring a tree according to the technique of separation of variables and evaluation;
FIG. 5 is a schematic view of an exemplary network portion to which the optimization method according to the invention is applied,
Fig. 6 is a table giving examples of initialization pressure intervals for different nodes of the network part of Fig. 5;
Fig. 7 is a table giving examples of initialization rate intervals for different arcs of the grating portion of Fig. 5;
Fig. 8 is a table giving the results of tests carried out on the network part of Fig. 5;
Figure 9 is a table giving results of the pressure intervals for the different nodes of the network part of Figure 5 in the cases of the table of Figure 8 where the propagation is not stopped;
Figure 10 is a table giving the results of the flow intervals for the different arcs of the network part of Figure 5 in the cases of the table of Figure 8 where the propagation is not stopped;
Figure 11 is a flowchart illustrating an exemplary implementation of the optimization method according to the invention;
Fig. 12 is a diagram showing a calculation tree which represents the propagation / backpropagation of constraints; and
Figure 13 is a schematic view of an example of natural gas transport network to which the invention is applicable.

La présente invention s'applique d'une manière générale à tous les réseaux de transport de gaz, notamment de gaz naturel, même si ces réseaux sont très étendus, à l'échelle d'un pays ou d'une région. De tels réseaux peuvent comprendre plusieurs milliers de canalisations, plusieurs centaines de vannes de régulation, plusieurs dizaines stations de compression, plusieurs centaines de ressources (points d'entrée du gaz sur le réseau) et plusieurs milliers de consommations (points de sortie du gaz sur le réseau).The present invention generally applies to all gas transport networks, including natural gas, even if these networks are very extensive, at the scale of a country or a region. Such networks may include several thousand pipes, several hundred control valves, several tens of compression stations, several hundred resources (gas entry points on the network) and several thousand consumptions (gas exit points on the network). the network).

Le procédé selon l'invention vise à déterminer automatiquement tous les degrés de liberté d'un réseau en régime permanent et ceci de manière optimale.The method according to the invention aims to automatically determine all the degrees of freedom of a steady-state network and this in an optimal manner.

Les valeurs sont optimales dans le sens où les contraintes ne sont pas transgressées et un critère économique est minimisé ou, si cela n'est pas possible, les contraintes sont transgressées au minimum.The values are optimal in the sense that the constraints are not transgressed and an economic criterion is minimized or, if this is not possible, the constraints are transgressed to a minimum.

Les degrés de liberté sont les pressions, les débits, les démarrages de compresseurs, les états ouvert/fermé, en ligne/en dérivation (by-pass) et les orientations directe ou inverse des ouvrages actifs.The degrees of freedom are pressures, flow rates, compressor starts, open / closed, in-line / bypass conditions, and direct or reverse orientations of active works.

Pour un réseau réel, il existe plusieurs centaines de variables à valeur entière (par exemple 1 pour ouvert et 0 pour fermé) en sus des quelques milliers de variables continues (pressions et débits).For a real network, there are several hundred variables with integer value (for example 1 for open and 0 for closed) in addition to a few thousand continuous variables (pressures and flows).

Le procédé selon l'invention rend possible le lancement de calcul en série, c'est-à-dire sans intervention humaine. Ce caractère autonome du calcul est d'un intérêt majeur dans un contexte de réseaux pouvant donner lieu à une multiplicité de scénarios d'acheminement.The method according to the invention makes it possible to launch the calculation in series, that is to say without human intervention. This autonomous character of the calculation is of major interest in a context of networks that can lead to a multiplicity of routing scenarios.

La Figure 1 est un schéma-bloc illustrant les principaux modules mis en oeuvre dans le cadre de la définition d'un réseau de transport de gaz.Figure 1 is a block diagram illustrating the main modules implemented in the context of the definition of a gas transport network.

Le module 5 constitue un modeleur qui est un ensemble permettant la modélisation du réseau. On entend par là sa description physique via ses ouvrages et sa structure (sous-réseaux connexes, blocs de pression,...). Ce modeleur inclut également de préférence des moyens de simulation (ou d'équilibrage) du réseau en débits et pressions.Module 5 is a modeler that is a set for modeling the network. By this is meant its physical description via its structures and structure (related sub-networks, pressure blocks, ...). This modeler also preferably includes means for simulating (or balancing) the network in flows and pressures.

Le module 8 constitue par sa part un coeur de calcul permettant d'assurer une optimisation du réseau.The module 8 is in turn a computational heart to ensure optimization of the network.

Le module d'optimisation 8 comprend essentiellement un solveur 6 dont les fonctions (notamment mise en oeuvre de la technique de séparation des variables et d'évaluation) seront explicitées plus loin et un solveur non linéaire convexe 7 qui peut agir en complément du solveur 6.The optimization module 8 essentially comprises a solver 6 whose functions (in particular implementation of the technique of separation of variables and evaluation) will be explained later and a convex nonlinear solver 7 which can act in addition to the solver 6 .

La Figure 2 montre de façon schématique une partie de réseau de transport de gaz comprenant différents points de prélèvements de gaz pour des consommations locales C. A chaque point de prélèvement est associée une contrainte de pression dépendant des besoins de consommation.Figure 2 schematically shows a portion of a gas transport network comprising different gas sampling points for local consumption C. At each sampling point is associated a pressure constraint depending on consumption needs.

La partie du réseau de transport comprend également des points d'approvisionnement en gaz pour fournir au réseau du gaz à partir de ressources locales R qui peuvent être par exemple des réserves de gaz stockées dans des cavités souterraines.The portion of the transmission network also includes gas supply points for supplying the gas network from local resources R which may be, for example, gas reserves stored in underground cavities.

La capacité du tronçon de réseau dépend à la fois du niveau des consommations C et des mouvements d'approvisionnement à partir des ressources R.The capacity of the network section depends both on the level of consumption C and on the supply movements from resources R.

Dans un réseau de transport de gaz, la pression du gaz décroît progressivement au cours du transit. Afin que le gaz soit acheminé avec un respect de la contrainte de pression admissible pour le consommateur, le niveau de pression doit être remonté régulièrement à l'aide de stations de compression réparties sur le réseau.In a gas transport network, the gas pressure gradually decreases during transit. In order for the gas to be conveyed with respect to the permissible pressure constraint for the consumer, the pressure level must be regularly raised using compression stations distributed over the network.

Chaque station de compression comprend au moins un compresseur et compte généralement de 2 à 12 compresseurs, la puissance totale des machines installées pouvant être comprise entre environ1 MW et 50 MW.Each compressor station comprises at least one compressor and generally has from 2 to 12 compressors, the total power of the machines installed being between about 1 MW and 50 MW.

La pression au refoulement des compresseurs ne doit pas dépasser la pression maximale de service (PMS) de la canalisation.The discharge pressure of the compressors must not exceed the maximum operating pressure (PMS) of the pipeline.

La Figure 3 illustre un exemple de configuration d'une station de compression qui se situe en même temps en un point d'interconnexion 1.0 du réseau. Une première canalisation 100 d'approvisionnement est raccordée au point d'interconnexion 1.0. Une deuxième canalisation d'approvisionnement sur laquelle est placée une vanne de régulation de pression 30 est également raccordée au point d'interconnexion 1.0. Un ou plusieurs compresseurs 40 sont disposés sur une troisième canalisation qui prend naissance au point d'interconnexion ou jonction 1.0.Figure 3 illustrates an example configuration of a compression station which is at the same time a point of interconnection 1.0 of the network. A first supply line 100 is connected to the interconnection point 1.0. A second supply pipe on which is placed a pressure regulating valve 30 is also connected to the interconnection point 1.0. One or more compressors 40 are disposed on a third pipe which originates at the point of interconnection or junction 1.0.

Selon un exemple de réalisation typique, on peut avoir une pression de 51 bar sur la première canalisation 100, une pression de 59 bar sur la deuxième canalisation en amont de la vanne de régulation 30, une pression de 51 bar sur la deuxième canalisation en aval de la vanne de régulation 30 et une pression de 67 bar sur la troisième canalisation en aval des compresseurs 40.According to a typical embodiment, it is possible to have a pressure of 51 bar on the first pipe 100, a pressure of 59 bar on the second pipe upstream of the regulating valve 30, a pressure of 51 bar on the second pipe downstream. of the control valve 30 and a pressure of 67 bar on the third pipe downstream of the compressors 40.

La présente invention vise à optimiser automatiquement les mouvements de gaz sur des réseaux complexes, le procédé offrant à la fois une grande robustesse et une grande précision.The present invention aims to automatically optimize gas movements on complex networks, the process offering both great robustness and high accuracy.

Dans la suite de la description, on considèrera que l'expression « ouvrage actif » englobe les vannes de régulation et les stations de compression ainsi que les vannes d'isolement, les ressources et les stockages.In the following description, it will be considered that the term "active structure" includes control valves and compressor stations as well as isolation valves, resources and storage.

L'expression « ouvrage passif » recouvre les canalisations et les résistances.The term "passive work" covers pipes and resistances.

Le procédé selon l'invention a pour but de rechercher les réglages adéquats pour les ouvrages actifs et d'établir une carte de débits et de pressions du réseau, pour optimiser un critère.The method according to the invention aims to find the appropriate settings for active structures and establish a map of flows and network pressures, to optimize a criterion.

Le critère est composé des trois termes différents :

le régime de pression : minimise ou maximise les pressions en aval des stockages et ressources, en amont et aval des stations de compression et des vannes de régulation et en amont des consommations,
l'énergie : minimise la consommation d'énergie de compression,
la cible : maximise ou minimise le débit d'un arc ou la pression d'un noeud particulier.

The criterion is composed of three different terms:

the pressure regime: minimizes or maximizes the downstream pressure of storage and resources, upstream and downstream of compressor stations and control valves and upstream of consumption,
energy: minimizes compression energy consumption,
the target: maximizes or minimizes the flow of an arc or the pressure of a particular node.

Dans le problème mathématique d'optimisation, ce critère est appelé fonction objectif. Dans cette fonction, chaque terme est pondéré par un coefficient (α, β et γ) qui lui donne plus ou moins d'importance : $g = α \times Régime + β \times Energie + γ \times cible$

In the mathematical problem of optimization, this criterion is called objective function. In this function, each term is weighted by a coefficient (α, β and γ) which gives it more or less importance:

boy Wut = α \times Diet + β \times Energy + γ \times target

Les degrés de liberté sont :

les pressions en chaque noeud,
les débits en chaque arc,

pour les variables continues, qui peuvent prendre toutes les valeurs comprises dans un intervalle.The degrees of freedom are:

the pressures in each node,
the flows in each arc,

for continuous variables, which can take any value within an interval.

Les degrés de liberté sont :

l'ouverture/la fermeture des ouvrages actifs,
le by-pass des stations de compression et des vannes de régulation,
l'orientation des stations de compression et des vannes de régulation,
le démarrage des compresseurs,

pour les paramètres discrets ou variables discrètes, qui ne peuvent prendre qu'un nombre fini de valeurs.The degrees of freedom are:

the opening / closing of active structures,
the bypass of compressor stations and control valves,
the orientation of compressor stations and control valves,
compressor startup,

for discrete parameters or discrete variables, which can only take a finite number of values.

Le but est de trouver les valeurs des variables qui minimisent le critère économique. La recherche des valeurs des variables est soumise à des contraintes de différents types :

des contraintes d'égalité : loi de perte de charge dans les canalisations, loi des noeuds. Ces contraintes sont intrinsèques au réseau, on ne peut donc les transgresser ;
des contraintes d'inégalités : contraintes de débit minimal et maximal, de pression minimale et maximale des ouvrages, contraintes de puissance de compression des stations, contraintes de vitesse minimale et maximale du gaz à chaque noeud, contraintes déprimogènes pour les vannes de régulation et pour les stations de compression, contraintes de pompage et de gavage des turbocompresseurs, contraintes de pressions de refoulement minimale et maximale des compresseurs, contraintes d'énergies journalières minimale et maximale des consommations, etc. Ces contraintes sont inhérentes aux ouvrages du réseau ou liées aux contraintes contractuelles du réseau (exemple : pression minimale pour un client); elles donnent des limites à ne pas dépasser, mais certaines peuvent être transgressées.

The goal is to find the values of the variables that minimize the economic criterion. The search for variable values is subject to constraints of different types:

constraints of equality: law of pressure loss in the pipes, law of the nodes. These constraints are intrinsic to the network, so we can not transgress them;
Inequality constraints: minimum and maximum flow rate constraints, minimum and maximum pressure of the structures, station compressive power constraints, minimum and maximum gas velocity constraints at each node, pressure constraints for the control valves and for compressor stations, pumping and force-feeding constraints for turbochargers, minimum and maximum compressive pressure constraints on compressors, minimum and maximum daily energy constraints, etc. These constraints are inherent in the network structures or related to the contractual constraints of the network (example: minimum pressure for a customer); they give limits not to be exceeded, but some can be transgressed.

Mathématiquement, ces contraintes sont de deux types : linéaires ou non linéaires.Mathematically, these constraints are of two types: linear or nonlinear.

Pour effectuer une modélisation d'un réseau de transport de gaz dans sa totalité, on peut considérer qu'à chaque état d'un ouvrage actif correspond une variable binaire e (qui prend la valeur 1 lorsque l'état est actif ou 0 dans le cas contraire, par exemple 1 pour ouvert et 0 pour fermé). On peut ainsi modéliser le choix entre chacun des états uniquement avec des contraintes linéaires. Le principe est illustré ci-dessous dans le cas d'une station de compression.To carry out a modeling of a gas transport network in its entirety, it can be considered that at each state of an active structure corresponds a binary variable e (which takes the value 1 when the state is active or 0 in the otherwise, for example 1 for open and 0 for closed). It is thus possible to model the choice between each state only with linear constraints. The principle is illustrated below in the case of a compressor station.

Exemple pour une station de compression :

Soit x = (Q,P_amont,P_aval) le triplet des variables continues de débits Q et de pressions P_amont et P_aval de la station de compression.
Soient e_f, e_b, e_d, e_i les 4 variables binaires associées aux 4 états alternatifs - fermé, bipassé, direct et inverse - ne pouvant se produire simultanément. Soient C_f (x), C_b(x), C_d(x), C_i(x), les 4 contraintes pour ces 4 états disjonctifs. Par exemple, pour l'état direct, C_d(x) est le vecteur des contraintes de débits minimal et maximal, de taux de compression minimal et maximal et de puissances minimale et maximale.
Soient C_{f max}, C_{b max}, C_{d max}, C_{i max} une estimation des valeurs maximales de ces contraintes quel que soit x. Dans l'exemple de l'état direct, C_{d max} est le vecteur des débits minimal et maximal, des taux de compression minimal et maximal et des puissances minimale et maximale.

Example for a compressor station:

Let x = (Q, P _upstream, P _downstream ) be the triplet of continuous variables of flow rates Q and pressures P _upstream and P _downstream of the compressor station.
Let e _f , _{b b} , e _d _{i be} the 4 binary variables associated with the 4 alternative states - closed, bypassed, direct and inverse - that can not occur simultaneously. Let C _f (x), C _b (x), C _d (x), C _i (x) be the 4 constraints for these 4 disjunctive states. For example, for the direct state, C _d (x) is the vector of minimum and maximum rate constraints, minimum and maximum compression ratio, and minimum and maximum powers.
Let C _max , C _max , C _{d max} , C _{max be} an estimation of the maximum values of these constraints whatever x. In the direct state example, C _{d max} is the vector of minimum and maximum flow rates, minimum and maximum compression ratios, and minimum and maximum powers.

Les contraintes linéaires vont donc s'écrire sous la forme :

C_f(x) ≤ (1 - e_f).C_{f max},
C_b(x) ≤ (1 - e_b).C_{b max},
C_d(x) ≤ (1 - e_d).C_{d max},
C_i(x) ≤ (1 - e_i).C_{i max},
e_f + e_b + e_d + e_i = 1 afin d'assurer le choix d'un et d'un seul des 4 états.

The linear constraints will thus be written in the form:

C _f (x) ≤ (1 - e _f ) .C _{f max} ,
C _b (x) ≤ (1 - e _b ) .C _{b max} ,
C _d (x) ≤ (1 - e _d ) .C _{d max} ,
C _i (x) ≤ (1 - e _i ) .C _{i max} ,
e _f + e _b + e _d + e _i = 1 to ensure the choice of one and only one of the four states.

En partant de ce principe, on peut également effectuer une modélisation en ne conservant que les 3 variables e_b, e_d, e_i, ce qui réduit la combinatoire.Starting from this principle, one can also carry out a modelization by keeping only the 3 variables e _b , e _d , e _i , which reduces the combinatorics.

Ces variables vont s'intégrer dans les contraintes de la manière suivante :

C_f(x) ≤ (e_b + e_d + e_i).C_{f max},
C_b(x) ≤ (1 - e_b).C_{b max},
C_d(x) ≤ (1 - e_d).C_{d max},
C_i(x) ≤ (1 - e_i).C_{i max},
e_b + e_d + e_i ≤ 1 afin d'assurer le choix entre un des 4 états, l'état fermé correspondant aux 3 variables nulles.

These variables will fit into the constraints as follows:

C _f (x) ≤ (e _b + e _d + e _i ) .C _{f max} ,
C _b (x) ≤ (1 - e _b ) .C _{b max} ,
C _d (x) ≤ (1 - e _d ) .C _{d max} ,
C _i (x) ≤ (1 - e _i ) .C _{i max} ,
e _b + e _d + e _i ≤ 1 to ensure the choice between one of the 4 states, the closed state corresponding to the 3 null variables.

Ainsi le problème de configuration optimale des ouvrages actifs est modélisé sous la forme d'un programme d'optimisation mixte (associant variables continues et variables binaires) non linéaire (car une partie des contraintes C_f(x), C_b(x), C_d(x), C_i(x) est non linéaire).Thus the problem of optimal configuration of active structures is modeled in the form of a non-linear mixed optimization program (associating continuous variables and binary variables) (because some of the constraints C _f (x), C _b (x), C _d (x), C _i (x) is nonlinear).

Le programme général peut donc s'écrire sous la forme suivante : $P_{0} {\begin{cases} \min_{(x e)} g (x) \\ C_{I} (x) + β . e \leq 0 \\ C_{E} (x) = 0 \\ x \in R^{n}, e \in {\{0 1\}}^{p} \end{cases}$

avec :

x, l'ensemble des variables des débits et de pressions (Q, P),
g(x), une fonction objectif a priori non linéaire. Il s'agit du critère économique (exemple : le coût de fonctionnement des ouvrages actifs tel que le gaz carburant consommé par les stations de compression),
C_i(x), l'ensemble des contraintes linéaires (contraintes de bornes) et non linéaires sur les ouvrages actifs ; ces contraintes sont des contraintes d'inégalité au nombre de p,
β, un vecteur dont les coefficients sont nuls ou égaux aux valeurs maximales des contraintes,
e, le vecteur des variables binaires, , de dimension p pour que l'équation le faisant intervenir soit cohérente mais le nombre de variables binaires est réellement : 3 x le nombre d'ouvrages actifs,
C_E(x), l'ensemble des contraintes d'égalité linéaires (exemple : loi des noeuds), et non linéaires (exemple: équations de perte de charge des canalisations). Elles sont au nombre de q.

The general program can be written in the following form:

P_{0} {\begin{cases} \min_{(x e)} boy Wut (x) \\ {VS}_{I} (x) + β . e \leq 0 \\ {VS}_{E} (x) = 0 \\ x \in R^{not}, e \in {\{0 1\}}^{p} \end{cases}

with:

x, the set of flow and pressure variables (Q, P),
g (x), a nonlinear objective function. This is the economic criterion (example: the operating cost of active structures such as the fuel gas consumed by compressor stations),
C _i (x), the set of linear constraints (boundary constraints) and nonlinear constraints on active structures; these constraints are inequality constraints in the number of p,
β, a vector whose coefficients are zero or equal to the maximum values of the stresses,
e, the vector of the binary variables,, of dimension p so that the equation involving it is coherent but the number of binary variables is really: 3 x the number of active works,
C _E (x), the set of linear equality constraints (example: node law), and non-linear constraints (example: pipe pressure loss equations). There are q of them.

Le procédé selon l'invention vise à fournir une réponse quel que soit l'état de saturation du réseau. C'est-à-dire qu'il est demandé au procédé d'autoriser, s'il ne peut faire autrement, de transgresser certaines contraintes afin de proposer un résultat, même en cas de saturation. L'autorisation de violation des contraintes est tempérée puisqu'on cherchera à la minimiser et qu'elle conduira de toute façon à un message de saturation. Prenant en compte cette demande, on modifie légèrement l'écriture du problème en introduisant les variables s qui, si elles sont non nulles, représentent la transgression des contraintes. $P_{1} {\begin{cases} \min_{(x s e)} f (x s) = g (x) + α \times {‖ s ‖}^{2} \\ C_{I} (x) + β . e \leq s_{I} \\ C_{E} (x) = s_{E} \\ x \in R^{n}, s_{I} \in R^{p}, s_{E} \in R^{q}, e \in {\{0 1\}}^{p} \end{cases}$

avec :

x est l'ensemble des variables des débits Q et des pressions P,
g(x) est la fonction objectif constituant le critère économique d'optimisation,
C_i(x) est l'ensemble des p contraintes d'inégalité linéaires et non linéaires sur les ouvrages actifs,
β est un vecteur dont les coefficients sont nuls ou égaux aux valeurs maximales des contraintes,
e est le vecteur des variables binaires, de dimension p pour que l'équation le faisant intervenir soit cohérente mais le nombre de variables binaires est réellement : 3 x le nombre d'ouvrages actifs,,
C_E(x) est l'ensemble des q contraintes d'égalité linéaires ou non linéaires,
s est une variable d'écart qui, lorsqu'elle est non nulle, représente la transgression d'une contrainte,
a est un coefficient représentant le degré d'autorisation de transgression de contraintes.

The method according to the invention aims to provide a response regardless of the state of saturation of the network. That is to say that the process is requested to allow, if it can not do otherwise, to violate certain constraints in order to propose a result, even in case of saturation. The authorization of violation of the constraints is tempered since one will seek to minimize it and that it will lead anyway to a message of saturation. Taking this request into account, the writing of the problem is slightly modified by introducing the variables s which, if they are non-zero, represent the transgression of the constraints.

P_{1} {\begin{cases} \min_{(x s e)} f (x s) = boy Wut (x) + α \times {‖ s ‖}^{2} \\ {VS}_{I} (x) + β . e \leq s_{I} \\ {VS}_{E} (x) = s_{E} \\ x \in R^{not}, s_{I} \in R^{p}, s_{E} \in R^{q}, e \in {\{0 1\}}^{p} \end{cases}

with:

x is the set of variables of flow rates Q and pressures P,
g (x) is the objective function constituting the economic criterion of optimization,
C _i (x) is the set of p linear and nonlinear inequality constraints on active structures,
β is a vector whose coefficients are zero or equal to the maximum values of the stresses,
e is the vector of the binary variables, of dimension p so that the equation involving it is coherent but the number of binary variables is really: 3 x the number of active works,
C _E (x) is the set of q linear or nonlinear equality constraints,
s is a variance variable which, when nonzero, represents the transgression of a constraint,
a is a coefficient representing the degree of authorization of constraint transgression.

Notons que, à variables binaires fixées, le programme P₁, qui n'est pas rigoureusement équivalent à P₀, a une solution proche de P₀ si le coefficient α est choisi suffisamment grand puisqu'alors les variables d'écart s_I et s_E sont recherchées très fortement proches de 0.Let us note that, with fixed binary variables, the program P ₁ , which is not rigorously equivalent to P ₀ , has a solution close to P ₀ if the coefficient α is chosen sufficiently large since then the variables of differences s _I and s _E are searched very strongly close to 0.

C'est un problème combinatoire de taille importante puisqu'il compte plusieurs centaines de variables à valeur entière en sus des quelques milliers de variables continues.This is a large combinatorial problem since it has several hundred full-valued variables in addition to a few thousand continuous variables.

Cette mixité du type des variables oblige à faire de l'optimisation combinatoire et continue. C'est pourquoi plusieurs méthodes mathématiques pouvant s'adapter à ces deux types d'optimisation sont de préférence combinées de façon hybride afin d'obtenir au final une solution exacte.This mixing of the type of the variables makes it necessary to make combinatorial and continuous optimization. This is why several mathematical methods that can be adapted to these two types of optimization are preferably combined in a hybrid manner in order to obtain an exact solution in the end.

Le procédé selon l'invention met en oeuvre en premier lieu une technique de séparation des variables et d'évaluation, dite « Branch & Bound » (notée par la suite B&B), c'est-à-dire « ramifier » et « borner ». Cette technique couvre une classe de méthodes d'optimisation capables de traiter des problèmes où interviennent des variables discrètes. Le caractère discret d'une variable s'oppose au caractère continu :

une variable continue peut prendre n'importe quelle valeur dans un intervalle donné. Dans le cadre du calcul de réseau, ce sera le cas des pressions exprimées en bar, par exemple : P ∈ [40,80],
une variable discrète ne peut prendre qu'un certain nombre de valeurs. Il s'agit souvent de variables binaires qui représentent par exemple le sens de fonctionnement d'une station de compression par exemple x = 0 (sens direct) ou x = 1 (sens indirect).

The method according to the invention implements in the first place a technique of separation of variables and evaluation, called "Branch &Bound" (denoted by the following B & B), that is to say "branch" and "bound"". This technique covers a class of optimization methods that can handle problems involving discrete variables. The discrete character of a variable is opposed to the continuous character:

a continuous variable can take any value in a given interval. As part of the network calculation, this will be the case for pressures expressed in bar, for example: P ∈ [40.80],
a discrete variable can only take a number of values. They are often binary variables which represent, for example, the operating direction of a compression station, for example x = 0 (forward direction) or x = 1 (indirect direction).

La méthode B&B est une méthode arborescente et consiste, au fur et à mesure qu'on construit l'arbre, à réduire le domaine de variation des variables. Cette méthode est couramment utilisée pour obtenir le minimum global d'un problème d'optimisation mettant en jeu des variables binaires.The B & B method is a tree method and consists, as the tree is constructed, in reducing the range of variation of the variables. This method is commonly used to obtain the overall minimum of an optimization problem involving binary variables.

Pour résoudre avec la méthode B&B un problème mixte qui est un problème traitant à la fois des variables discrètes et continues, plusieurs variantes peuvent être envisagées :

B&B₁: la méthode B&B ne sépare que sur les variables binaires. Les variables sont représentées par des intervalles. On pourra alors calculer les bornes de la fonction objectif en utilisant l'arithmétique des intervalles.
B&B₂: la méthode B&B sépare à la fois sur les variables binaires et continues ; ceci implique une représentation par intervalles. Dans ce cas, le principe de séparation (branch) va consister à découper l'espace de définition des variables continues au lieu de fixer les variables discrètes à l'une de leur valeur. Ainsi, on va explorer séparément des parties de l'ensemble réalisable et borner (bounding) l'intervalle de variation de la fonction objectif sur ces sous-parties.

To solve with the B & B method a mixed problem that is a problem dealing with both discrete and continuous variables, several variants can be envisaged:

B & B ₁ : the B & B method only separates on binary variables. The variables are represented by intervals. We can then calculate the limits of the objective function using the arithmetic of the intervals.
B & B ₂ : The B & B method separates both binary and continuous variables; this implies an interval representation. In this case, the separation principle (branch) will consist in cutting the definition space of the continuous variables instead of setting the discrete variables to one of their value. Thus, we will separately explore parts of the feasible set and bounding the range of variation of the objective function on these subparts.

La mise en place d'une méthode B&B de séparation des variables et d'évaluation nécessite donc le choix de stratégies concernant :

la sélection du noeud à examiner : en fonction de la date d'arrivée des noeuds dans la pile, de leur positionnement ou de la valeur d'une fonction de mérite calculée en chaque noeud candidat,
l'évaluation des bornes de la solution courante qui permet d'avancer dans la méthode B&B,
l'élimination des noeuds ne pouvant pas contenir l'optimum (test de contraintes transgressées, de valeur d'objectif moins bonne que la valeur courante, utilisation de la monotonie de la fonction objectif),
la séparation du noeud courant en noeuds fils (2 ou plus) en divisant le domaine de variation d'une ou plusieurs variables (choisie(s) selon des critères basés sur le diamètre d'intervalles rattachés à la (aux) variable(s), le diamètre ou la largeur d'un intervalle correspondant à la différence entre sa borne maximale et sa borne minimale),
le critère d'arrêt basé sur le temps d'exécution ou sur l'évaluation de certains diamètres.

The implementation of a B & B method of separation of variables and evaluation therefore requires the choice of strategies concerning:

the selection of the node to be examined: according to the arrival date of the nodes in the stack, their positioning or the value of a merit function calculated in each candidate node,
the evaluation of the limits of the current solution which allows to advance in the B & B method,
the elimination of the nodes that can not contain the optimum (test of transgressed constraints, of objective value less good than the current value, use of the monotony of the objective function),
separating the current node into child nodes (2 or more) by dividing the variation domain of one or more variables (chosen according to criteria based on the interval diameter attached to the variable (s) , the diameter or width of an interval corresponding to the difference between its maximum and minimum limits),
the stopping criterion based on the execution time or the evaluation of certain diameters.

Pour le problème de la configuration optimale des ouvrages actifs, les méthodes B&B consistent à fixer progressivement l'état des ouvrages actifs, et à évaluer à chaque étape, parmi ces combinaisons partielles, celles qui sont susceptibles de mener à la combinaison globale la plus favorable.For the problem of the optimal configuration of the active structures, the B & B methods consist in progressively fixing the state of the active structures, and to evaluate at each stage, among these partial combinations, those which are likely to lead to the most favorable overall combination. .

Un exemple va être décrit en référence à la Figure 4.An example will be described with reference to Figure 4.

On considère un réseau de gaz sur lequel on aurait plusieurs stations de compression. On cherche, par exemple, à minimiser le gaz carburant sur le réseau. Si l'on choisit au début de l'arbre de B&B, la station de compression n° 1 et que l'on teste la variable binaire associée à son état $(e_{d}^{1} = 1) .$

We consider a gas network on which we would have several compressor stations. One seeks, for example, to minimize the fuel gas on the network. If one chooses at the beginning of the B & B tree, the compression station n ° 1 and that one tests the binary variable associated with its state

(e_{d})

(_{1} = 1) .

$f_{\min}^{I}$

est la borne minimale de la fonction objectif calculée au noeud i, sachant l'ensemble des décisions qui ont déjà été prises.

f_{\min}^{I}

is the minimum bound of the objective function computed at node i, knowing all the decisions that have already been made.

$f_{\max}^{I}$

est la borne maximale de la fonction objectif associée à la meilleure combinaison d'états connus lors de l'exploration du noeud i.

f_{\max}^{I}

is the maximum bound of the objective function associated with the best combination of known states when exploring node i.

Si $f_{\min}^{1} > f_{\max}^{1} (avec f_{\max}^{1} = f_{\max}^{0}),$

alors il est sûr que la station 1 orientée en sens indirect

(e_{d}^{1} = 0)

ne conduit pas à la solution optimale.Yes

f_{\min}^{} > f_{\max}^{} (with f_{\max}^{} = f_{\max}^{}),

then it is sure that the station 1 oriented in indirect direction

(e_{d})

(_{1} = 0)

does not lead to the optimal solution.

En revanche, si $f_{\min}^{1} \leq f_{\max}^{1}$

l'exploration continue en fixant une autre variable binaire. Toutes les variables binaires sont ainsi fixées progressivement. Si on n'effectue aucune coupe dans une branche, on obtient une configuration réalisable, c'est-à-dire que l'on a fixé l'ensemble des variables binaires et que l'on respecte l'ensemble des contraintes.On the other hand, if

f_{\min}^{} \leq f_{\max}^{}

exploration continues by setting another binary variable. All binary variables are thus fixed progressively. If we do not make any cuts in a branch, we obtain a feasible configuration, that is to say that we fixed the set of binary variables and that we respect all the constraints.

Différentes techniques peuvent être associées à la technique de séparation des variables et d'évaluation.Different techniques can be associated with the technique of separation of variables and evaluation.

On peut en particulier utiliser la propagation de contraintes qui permet d'exploiter les informations de l'équation ou de l'inéquation pour diminuer les intervalles des variables de cette équation.One can in particular use the propagation of constraints which makes it possible to exploit the information of the equation or the inequation to decrease the intervals of the variables of this equation.

On ne considère que l'équation non linéaire C(x) et, de manière générale, on cherche à résoudre : $C (x) \in [a b] \subseteq IR où x \in X \subseteq {IR}^{n}$

avec :

IR est l'ensemble des intervalles,
X est un vecteur d'intervalles de dimension n.

We consider only the non linear equation C (x) and, in general, we try to solve:

VS (x) \in [at b] \subseteq IR or x \in X \subseteq {IR}^{not}

with:

IR is the set of intervals,
X is a vector of intervals of dimension n.

La propagation de contraintes peut être basée sur la construction d'un arbre de calcul qui représente C(x). Dans un premier temps, on calcule à partir des feuilles de l'arbre, qui sont les variables et les constantes, la valeur des noeuds intermédiaires et du noeud racine correspondant à la valeur de la contrainte (ce qui équivaut à appliquer les règles de l'arithmétique d'intervalle), puis on propage la valeur de l'intervalle de la contrainte, depuis la racine de l'arbre vers les feuilles pour réduire les espaces de définition des variables.Constraint propagation can be based on the construction of a computational tree that represents C (x) . First, we calculate from the leaves of the tree, which are the variables and the constants, the value of the intermediate nodes and the root node corresponding to the value of the constraint (which is equivalent to applying the rules of the interval arithmetic), then propagate the value of the constraint interval from the root of the tree to the leaves to reduce the definition spaces of the variables.

L'algorithme de propagation d'une contrainte sur ses variables est le suivant :

Etape 1, propagation : construction de l'arbre de calcul de la contrainte C, les feuilles sont les variables intervalles x _i ou des constantes réelles,
dans chaque noeud est stocké le résultat de l'opération partielle et unitaire qu'il représente, par exemple x _a + x _b,
le dernier calcul est effectué à la racine.
Etape 2, rétro-propagation : descente de l'arbre de la racine vers les feuilles.
A chaque noeud, on tente de réduire le résultat partiel calculé en 1.
Par exemple : x _a + x _b = [a,b] ⇒ x _a := ([a,b]-x _b) ∩ x _a et x _b := ([a,b]- x _a) ∩ x _b

The algorithm for propagating a constraint on its variables is as follows:

Step 1, propagation: construction of the computation tree of the constraint C, the leaves are the interval variables x _i or real constants,
in each node is stored the result of the partial and unit operation it represents, for example x _a + x _b ,
the last calculation is done at the root.
Step 2, back propagation: descent of the tree from the root to the leaves.
At each node, we try to reduce the partial result calculated in 1.
For example: x _a + x _b = [a, b] ⇒ x _a : = ([a, b] - x _b ) ∩ x _a and x _b : = ([a, b] - x _a ) ∩ x _b

La Figure 12 illustre la propagation/rétro-propagation des contraintes pour l'équation suivante donnée à titre d'exemple : $2 x_{3} x_{2} + x_{1} = 3 avec x_{1} = [1 3], pour i \in \{1 2 3\}$

Figure 12 illustrates the propagation / backpropagation of constraints for the following equation given as an example:

2 x_{3} x_{2} + x_{1} = 3 with x_{1} = [1 3], for i \in \{1 2 3\}

La première étape de l'algorithme est présentée dans la partie gauche de la Figure 12 : en partant des valeurs des variables et des constantes, on effectue chaque opération unitaire constituant l'expression jusqu'à obtenir la valeur du membre de gauche de l'expression en haut de l'arbre ; ce noeud est le noeud racine.The first step of the algorithm is presented in the left part of Figure 12: starting from the values of the variables and the constants, one carries out each unit operation constituting the expression until obtaining the value of the member of left of the expression at the top of the tree; this node is the root node.

La seconde étape de l'algorithme est explicitée par la partie droite de la Figure 12 : on veut le membre de gauche égal à une valeur précise, on redescend donc l'arbre à partir de la racine, grâce aux opérations inverses à celles utilisées dans la première partie, on cherche à réduire les intervalles de chaque noeud et particulièrement celui des variables. Dans l'exemple, la propagation a permis de réduire chaque intervalle des variables de [1,3] à [1,1], c'est-à-dire que les variables ont été instanciées à 1, uniquement grâce à la propagation.The second step of the algorithm is explained by the right part of Figure 12: we want the left-hand member equal to a precise value, we go back down the tree from the root, thanks to the operations reversed to those used in the first part, we try to reduce the intervals of each node and particularly that of the variables. In the example, the propagation allowed to reduce each interval of the variables from [1,3] to [1,1], that is to say that the variables were instantiated at 1, only thanks to the propagation.

L'algorithme de propagation sur l'ensemble des contraintes d'un problème s'effectue de la façon suivante :The propagation algorithm on all the constraints of a problem is carried out as follows:

1. Initialisation de la file des contraintes à propager 1. Initialization of the queue of constraints to propagate

Pour cela, on insère, sans doublon, dans une file triée selon un critère de mérite M toutes les contraintes.For that, one inserts, without doubloon, in a sorted file according to a criterion of merit M all the constraints.

2. Boucle sur la file des contraintes 2. Loop on the queue of constraints

          Tant que la file n'est pas vide {
    Extraction de la "meilleure" contrainte c (pour le critère M)
    Propagation de C
    Si la propagation a conduit à un intervalle vide pour au moins une
   variable {
     Sortie de la boucle : il n'y a pas de solution au problème
    }
    Sinon {
     Pour chaque variable modifiée par la propagation de C {
       Pour chaque contrainte où intervient cette variable {
         Si la contrainte n'est pas déjà résolue, ajout dans la file
           }
         }
       }
     }As long as the queue is not empty {  Extraction of the "best" constraint c (for the criterion M) Propagation of C If the propagation has led to an empty interval for at least one variable {Output of the loop: there is no solution to the problem} Otherwise { For each variable modified by the propagation of C {For each constraint in which this variable intervenes {If the constraint is not already resolved, add in the file}}}}

Selon un exemple de réalisation, seul « l'âge » de la contrainte intervient dans le critère de mérite M, i.e. la file est équivalente à une pile FIFO. Cependant, on peut utiliser un critère plus complexe. Par exemple, une variable très réduite par la propagation d'une contrainte pourrait conduire à insérer les contraintes dans lesquelles elle intervient dans la file avec un mérite élevé.According to an exemplary embodiment, only the "age" of the constraint intervenes in the merit criterion M, i.e. the queue is equivalent to a FIFO stack. However, a more complex criterion can be used. For example, a variable greatly reduced by the propagation of a constraint could lead to insert the constraints in which it intervenes in the queue with a high merit.

On notera qu'une contrainte est dite résolue lorsqu'elle est déjà vérifiée quelles que soient les valeurs que les variables prennent dans leurs intervalles (autrement dit si l'intervalle résultant de la propagation sur la contrainte ne contient que des valeurs acceptables).It will be noted that a constraint is said to be solved when it is already checked whatever the values that the variables take in their intervals (in other words if the interval resulting from the propagation on the constraint contains only acceptable values).

Pour une contrainte C de fonction d'inclusion C(X) = └C(X),C(X)┘, C est résolue si :

C est une contrainte d'égalité et C(X) = 0,
C est une contrainte d'inégalité positive et C(X) ≥ 0,
C est une contrainte d'inégalité négative et C(X) ≤ 0.

For a constraint C of inclusion function C (X) = └ C (X) , C (X) ┘ , C is resolved if:

C is an equality constraint and C (X) = 0,
C is a positive inequality constraint and C (X) ≥ 0,
C is a negative inequality constraint and C (X ) ≤ 0.

Lorsqu'une contrainte est résolue, sa propagation ne conduira plus à aucune réduction des intervalles de ses variables.When a constraint is resolved, its propagation will no longer lead to any reduction in the intervals of its variables.

La technique de propagation de contraintes peut être utilisée par exemple pour déterminer l'orientation des ouvrages actifs des réseaux de transport de gaz. Les ouvrages actifs peuvent être simplement considérés orientés en sens direct lorsque le débit est positif et en sens indirect lorsque le débit est négatif. On peut aussi effectuer une modélisation complète de la configuration des ouvrages actifs en faisant intervenir 3 ou 4 variables binaires comme indiqué plus haut. La mise en oeuvre de la technique de la propagation des contraintes peut être effectuée à l'aide d'une bibliothèque d'arithmétique des intervalles et de propagation des contraintes capable de traiter des variables discrètes.The constraint propagation technique can be used, for example, to determine the orientation of the active structures of the gas transport networks. Active structures can simply be considered oriented in a forward direction when the flow is positive and in the indirect direction when the flow is negative. It is also possible to perform a complete modeling of the configuration of the active structures by involving 3 or 4 binary variables as indicated above. The implementation of the constraint propagation technique can be performed using an interval arithmetic and constraint propagation library capable of processing discrete variables.

Les méthodes de propagation des contraintes peuvent d'une part servir à réduire la combinatoire dans des temps réduits, lors d'une première étape pouvant précéder un processus d'optimisation, exact ou approché, et d'autre part être intégrées aux méthodes B&B pour permettre en chaque noeud un meilleur calcul des bornes de la fonction objectif et éventuellement des coupes supplémentaires.Stress propagation methods can, on the one hand, serve to reduce the combinatorics in reduced times, in a first step that can precede an optimization process, exact or approximate, and on the other hand be integrated into the methods. B & B to allow in each node a better calculation of the terminals of the objective function and possibly additional cuts.

En particulier, dans ce dernier cas où la propagation des contraintes s'effectue au sein d'un noeud de l'arbre de recherche et est utilisée pour élaguer les noeuds qu'elle permet de déclarer infaisables, et pour diminuer le diamètre des intervalles des variables, on considère dans la file initiale des contraintes à propager, les contraintes où interviennent la ou les variables dont la séparation a conduit à la création du noeud en cours d'évaluation. Si ce noeud est la racine de l'arbre, alors on place toutes les contraintes dans la file.In particular, in the latter case where the propagation of the stresses is carried out within a node of the search tree and is used to prune the nodes it allows to declare infeasible, and to reduce the diameter of the intervals of variables, we consider in the initial queue constraints to propagate, the constraints in which intervene the variables whose separation led to the creation of the node being evaluated. If this node is the root of the tree, then we put all the constraints in the queue.

A titre d'exemple de mise en oeuvre d'une technique de propagation des contraintes, on se réfèrera aux Figures 5 à 10.As an example of implementing a stress propagation technique, reference is made to FIGS. 5 to 10.

Si l'on se reporte à la Figure 5, on voit un réseau simple de transport du gaz comprenant une ressource R, une consommation C, un premier compresseur CP1 et un deuxième compresseur CP2. Le réseau comprend des noeuds N₀ à N₄ (jonctions ou points d'interconnexion) et des arcs l à VII (canalisations ou tronçons comportant les compresseurs CP1, CP2, la ressource R et la consommation C).Referring to FIG. 5, there is shown a simple gas transport network comprising a resource R, a consumption C, a first compressor CP1 and a second compressor CP2. The network comprises nodes N ₀ to N ₄ (junctions or interconnection points) and arcs I to VII (pipes or sections comprising the compressors CP1, CP2, the resource R and the consumption C).

Le réseau définit cinq variables de pression aux noeuds N₀ à N₄ et sept variables de débit dans les arcs I à VII.The network defines five pressure variables at nodes N ₀ to N ₄ and seven flow variables in arcs I to VII.

La Figure 6 donne un exemple d'intervalles de pression d'initialisation (en bar) aux différents noeuds N₀ à N₄.Figure 6 gives an example of initialization pressure intervals (in bar) at the different nodes N ₀ to N ₄ .

La ressource A a une pression de consigne de 40 bar. C'est pourquoi son intervalle d'initialisation est un intervalle de largeur nulle.Resource A has a set pressure of 40 bar. This is why its initialization interval is an interval of zero width.

Le noeud N₄ de consommation a une pression minimale de livraison de 42 bar, d'où l'initialisation dans l'intervalle [40, 60].The N ₄ consumer node has a minimum delivery pressure of 42 bar, hence the initialization in the interval [40, 60].

La Figure 7 donne un exemple d'intervalles de débit d'initialisation (en m³/h) dans les arcs I à VII.Figure 7 gives an example of initialization rate intervals (in m ³ / h) in arcs I to VII.

La ressource R et la consommation C correspondant aux arcs I et VII ont des débits imposés à 800 000 m³/h. Leurs intervalles sont donc initialisés à des intervalles de largeur nulle.The resource R and the consumption C corresponding to arcs I and VII have flow rates imposed at 800,000 m ³ / h. Their intervals are therefore initialized at intervals of zero width.

Les arcs III et V sur lesquels se trouvent les compresseurs CP1 et CP2 respectivement présentent des intervalles de débit plus réduits que les arcs II, IV et VI correspondant à de simples canalisations.The arcs III and V on which the compressors CP1 and CP2 respectively are located have smaller flow intervals than arcs II, IV and VI corresponding to simple pipes.

Plusieurs tests sont effectués :

A. On teste tout d'abord toutes les combinaisons d'orientation des compresseurs CP1, CP2 (tests A1 à A4).
B. On laisse libre l'orientation du compresseur CP1 et on fixe celle du compresseur CP2 (tests B1 et B2).
C. On laisse libres les orientations des deux compresseurs CP1, CP2 (test C).

Several tests are carried out:

A. First of all, all the orientation combinations of the compressors CP1, CP2 (tests A1 to A4) are tested.
B. The orientation of the compressor CP1 is left free and that of the compressor CP2 (tests B1 and B2) is fixed.
C. The orientations of the two compressors CP1, CP2 (test C) are left free.

Les résultats de ces tests A1 à A4, B1, B2 et C sont présentés sur le tableau de la Figure 8.The results of these tests A1 to A4, B1, B2 and C are presented in the table of Figure 8.

Dans les trois cas où la propagation n'est pas stoppée (tests A1, B1 et C), on obtient les résultats identiques présentés sur les tableaux des Figures 9 et 10.In the three cases where the propagation is not stopped (tests A1, B1 and C), the identical results shown in the tables of FIGS. 9 and 10 are obtained.

La Figure 9 indique les intervalles de pression résultats (en bar) aux différents noeuds N₀ à N₄.Figure 9 shows the result pressure ranges (in bar) at the different nodes N ₀ to N ₄ .

La Figure 10 indique les intervalles de débit résultats (en m³/h) pour les différents arcs I à VII.Figure 10 shows the flow rate results (in m ³ / h) for the different arcs I to VII.

Dans ces exemples, on voit que les informations contenues dans les contraintes sont utilisées pour réduire les intervalles des variables et permettent également de fixer la valeur de certaines variables discrètes (ici l'orientation de chaque compresseur). En particulier, on voit que si l'orientation d'un ou des deux compresseurs est laissée libre, en appliquant la seule méthode de propagation des contraintes, on peut conclure que le compresseur libre doit s'orienter dans le sens direct.In these examples, we see that the information contained in the constraints are used to reduce the intervals of the variables and also make it possible to fix the value of certain discrete variables (here the orientation of each compressor). In particular, we see that if the orientation of one or both compressors is left free, applying the only method of propagation constraints, we can conclude that the free compressor must move in the forward direction.

La méthode de propagation des contraintes ainsi que la méthode de séparation des variables et d'évaluation (B&B) fait appel au calcul par intervalles dont les caractéristiques principales seront rappelées ci-dessous.The constraint propagation method as well as the method of separating variables and evaluation (B & B) uses interval calculus whose main characteristics will be recalled below.

En arithmétique par intervalles, on ne manipule plus des nombres, qui approchent plus ou moins fidèlement une valeur, mais des intervalles contenant cette valeur. Par exemple, on peut tenir compte d'une erreur de mesure en remplaçant une valeur mesurée x avec une incertitude ε par l'intervalle [x - ε,x + ε]. On peut également remplacer une valeur par sa plage de validité telle qu'une pression P d'une ressource représentée par un intervalle [4, 68] bar. Enfin, si l'on désire obtenir un résultat valide pour tout un ensemble de valeurs, on utilise un intervalle contenant ces valeurs. En effet, l'objectif de l'arithmétique par intervalles est de fournir des résultats qui contiennent à coup sûr la valeur ou l'ensemble cherché. On parle alors de résultats garantis, validés ou encore certifiés.In interval arithmetic, we no longer manipulate numbers, which approach a value more or less faithfully, but intervals containing this value. For example, a measurement error can be taken into account by replacing a measured value x with an uncertainty ε by the interval [x - ε, x + ε]. It is also possible to replace a value by its validity range such as a pressure P of a resource represented by an interval [4, 68] bar. Finally, if we want to obtain a valid result for a whole set of values, we use an interval containing these values. Indeed, the purpose of interval arithmetic is to provide results that surely contain the value or set sought. We then speak of guaranteed, validated or certified results.

Comme cela a été implicitement admis jusqu'à présent, les intervalles qui ne contiennent pas de « trou », sont des sous-ensembles fermés connexes de R. On notera IR l'ensemble des intervalles. On peut les généraliser en plusieurs dimensions : un vecteur intervalle x ∈ IR ⁿ est un vecteur dont les n composantes sont des intervalles et une matrice intervalle A ∈ IR ^mxn est une matrice dont les composantes sont des intervalles. Une représentation graphique d'un vecteur intervalle de IR, IR ² et IR ³ correspond respectivement à un segment de droite, à un rectangle et à un parallélépipède. Un vecteur intervalle est donc un hyper-parallélépipède. Par la suite, on utilisera indifféremment les termes de vecteur d'intervalles, de pavé, de boîte ou même d'intervalle.As has been implicitly accepted so far, the intervals that do not contain a "hole" are connected closed subsets of R. IR will be noted as the set of intervals. They can be generalized in several dimensions: an interval vector x ∈ IR ⁿ is a vector whose n components are intervals and an interval matrix A ∈ IR ^mxn is a matrix whose components are intervals. A graphical representation of an interval vector of IR, IR ² and IR ³ corresponds respectively to a line segment, a rectangle and a parallelepiped. An interval vector is therefore a hyper-parallelepiped. Subsequently, we will use the interval, block, box, or even interval vector terms interchangeably.

On note les objets intervalles par des caractères gras : x. On note x le minimum de x et x son maximum. On a alors x = [ x , x ], et on considère l'ordre partiel sur IR ⁿ : $X \leq Y \Leftrightarrow x_{i} \leq y_{i} pour i = 1.. n .$

We mark the interval objects with bold characters: x . We denote x the minimum of x and x its maximum. We then have x = [ x , x ], and we consider the partial order on IR ⁿ :

X \leq Y \Leftrightarrow x_{i} \leq {there}_{i} for i = 1 .. not .

On note w(x) la largeur de x, (avec w pour width) ou encore son diamètre : $w (x) = \overline{x} - \underset{̲}{x}$

We denote w ( x ) the width of x , (with w for width) or its diameter:

w (x) = \tilde{x} - \underset{}{x}

Le centre mid(x) et son rayon rad(x) sont définis par : $mid (x) = \frac{\overline{x} - \underset{̲}{x}}{2}$

rad (x) = \frac{\overline{x} - \underset{̲}{x}}{2} = \frac{w (x)}{2}

The center mid ( x ) and its radius rad ( x ) are defined by:

mid (x) = \frac{\tilde{x} - \underset{}{x}}{2}

rad (x) = \frac{\tilde{x} - \underset{}{x}}{2} = \frac{w (x) ()}{2}

Une fonction F : lR ⁿ → IR est une fonction d'inclusion de f sur X ∈ IR ⁿ. Si X ∈ X alors f(X) ∈ F(X).A function F : lR ⁿ → IR is an inclusion function of f on X ∈ IR ⁿ . If X ∈ X then f (X) ∈ F ( X ).

On désigne par l'adjectif « ponctuel » un objet numérique usuel (c'est-à-dire un nombre réel, ou un vecteur, une matrice de nombres réels) et on le confond avec l'intervalle de diamètre nul.The term "point" refers to a usual numerical object (that is, a real number, or a vector, a matrix of real numbers) and is confused with the zero-diameter interval.

Le résultat d'une opération ◇ entre deux intervalles x et y est le plus petit intervalle (au sens de l'inclusion) contenant tous les résultats de l'opération appliquée entre tous les éléments x de x et tous les éléments y de y, c'est-à-dire contenant l'ensemble : $\{x◊y; x \in x, y \in y\}$

The result of an operation ◇ between two intervals x and y is the smallest interval (in the sense of inclusion) containing all the results of the operation applied between all the elements x of x and all the elements y of y , that is, containing the set:

\{x◊y; x \in x, there \in there\}

De même, le résultat d'une fonction F(z) est le plus petit intervalle contenant l'ensemble : $\{f (Z); z \in z\}$

Similarly, the result of a function F (z) is the smallest interval containing the set:

\{f (Z); z \in z\}

Si l'on considère les opérateurs traditionnels +, -, x, ², / ou √, on peut définir les formules suivantes, plus utilisables en pratique que la définition théorique ci-dessus : $[\underset{̲}{x} \overline{x}] + [\underset{̲}{y} \overline{y}] = [\underset{̲}{x} + \underset{̲}{y}, \overline{x} + \overline{y}]$

[\underset{̲}{x} \overline{x}] - [\underset{̲}{y} \overline{y}] = [\underset{̲}{x} - \overline{y}, \overline{x} - \underset{̲}{y}]

[\underset{̲}{x} \overline{x}] \times [\underset{̲}{y} \overline{y}] = [\min (\underset{̲}{x} \times \underset{̲}{y}, \overline{x} \times \underset{̲}{y}, \underset{̲}{x} \times \overline{y}, \overline{x} \times \overline{y}), \max (\underset{̲}{x} \times \underset{̲}{y}, \overline{x} \times \underset{̲}{y}, \underset{̲}{x} \times \overline{y}, \overline{x} \times \overline{y})]

{[\underset{̲}{x} \overline{x}]}^{2} = {\begin{cases} [\min ({\underset{̲}{x}}^{2} {\overline{x}}^{2}), \max ({\underset{̲}{x}}^{2} {\overline{x}}^{2})] si 0 \notin [\underset{̲}{x} \overline{x}] \\ [\begin{matrix} 0, \max ({\underset{̲}{x}}^{2} {\overline{x}}^{2}) \end{matrix}] sinon \end{cases}

1 / [\underset{̲}{x} \overline{x}] = \min (1 / \underset{̲}{x}, 1 / \overline{x}), \max (1 / \underset{̲}{x}, 1 / \overline{x}) si 0 \notin [\underset{̲}{x} \overline{x}]

[\underset{̲}{x} \overline{x}] / [\underset{̲}{y} \overline{y}] = [\underset{̲}{x} \overline{x}] \times (1 / [\underset{̲}{y} \overline{y}]) si 0 \neq [\underset{̲}{y} \overline{y}]

\sqrt{[\underset{̲}{x} \overline{x}]} = [\sqrt{\underset{̲}{x}}, \sqrt{\overline{x}} si 0 \leq \underset{̲}{x}]

If we consider the traditional operators +, -, x, ² , / or √, we can define the following formulas, more usable in practice than the theoretical definition above:

[\underset{}{x} \tilde{x}] + [\underset{}{there} \tilde{there}] = [\underset{}{x} + \underset{}{there}, \tilde{x} + \tilde{there}]

[\underset{}{x} \tilde{x}] - [\underset{}{there} \tilde{there}] = [\underset{}{x} - \tilde{there}, \tilde{x} - \underset{}{there}]

[\underset{}{x} \tilde{x}] \times [\underset{}{there} \tilde{there}] = [\min (\underset{}{x} \times \underset{}{there}, \tilde{x} \times \underset{}{there}, \underset{}{x} \times \tilde{there}, \tilde{x} \times \tilde{there}), \max (\underset{}{x} \times \underset{}{there}, \tilde{x} \times \underset{}{there}, \underset{}{x} \times \tilde{there}, \tilde{x} \times \tilde{there})]

{[\underset{}{x} \tilde{x}]}^{2} = {\begin{cases} [\min ({\underset{}{x}}^{2} {\tilde{x}}^{2}), \max ({\underset{}{x}}^{2} {\tilde{x}}^{2})] if 0 \notin [\underset{}{x} \tilde{x}] \\ [\begin{matrix} 0, \max ({\underset{}{x}}^{2} {\tilde{x}}^{2}) \end{matrix}] if not \end{cases}

1 / [\underset{}{x} \tilde{x}] = \min (1 / \underset{}{x}, 1 / \tilde{x}), \max (1 / \underset{}{x}, 1 / \tilde{x}) if 0 \notin [\underset{}{x} \tilde{x}]

[\underset{}{x} \tilde{x}] / [\underset{}{there} \tilde{there}] = [\underset{}{x} \tilde{x}] \times (1 / [\underset{}{there} \tilde{there}]) if 0 \neq [\underset{}{there} \tilde{there}]

\sqrt{[\underset{}{x} \tilde{x}]} = [\sqrt{\underset{}{x}}, \sqrt{\tilde{x}} if 0 \leq \underset{}{x}]

Les propriétés algébriques traditionnelles (c'est-à-dire en arithmétique ponctuelle) telles que la réciprocité entre l'addition et la soustraction ou la distributivité de la multiplication par rapport à l'addition ne sont plus vérifiées :

la soustraction n'est plus la réciproque de l'addition. En effet : $x - x = \{x - y | x \in x, y \in x\} \supset \{x - x | x \in x\} = \{0\}$
la division n'est plus non plus la réciproque de la multiplication, par le même raisonnement que ci-dessus, on obtient : $x / x \supset \{1\}$
la multiplication d'un intervalle par lui-même n'est pas égale à l'élévation au carré. Prenons l'exemple où x = [-3,2] : $x \times x = [- 6, 9]$
$x^{2} = [0 9]$
la multiplication n'est pas distributive par rapport à l'addition. Prenons x = [-2,3], y = [1,4] et z = [-2,1] : $x \times (y + z) = [- 10, 15]$
$x \times y + x \times z = [14 16]$
La multiplication est en fait sous-distributive par rapport à l'addition, c'est-à-dire que : $x \times (y + z) \subset x \times y + x \times z$

Traditional algebraic properties (that is, point arithmetic) such as the reciprocity between addition and subtraction or multiplication distributivity over addition are no longer true:

subtraction is no longer the reciprocal of the addition. Indeed : $x - x = \{x - there | x \in x, there \in x\} \supset \{x - x | x \in x\} = \{0\}$
division is no longer the reciprocal of multiplication, by the same reasoning as above, we obtain: $x / x \supset \{1\}$
the multiplication of an interval by itself is not equal to squaring. Take the example where x = [-3,2]: $x \times x = [- 6, 9]$
$x^{2} = [0 9]$
the multiplication is not distributive compared to the addition. Take x = [-2,3], y = [1,4] and z = [-2,1]: $x \times (there + z) = [- 10, 15]$
$x \times there + x \times z = [14 16]$
The multiplication is in fact sub-distributive compared to the addition, that is to say that: $x \times (there + z) \subset x \times there + x \times z$

On peut aussi définir des fonctions élémentaires telles que le sinus, l'exponentielle, etc. prenant en argument des intervalles. Pour cela, on utilise la définition abstraite ci-dessus.We can also define elementary functions such as sine, exponential, and so on. taking in argument intervals. For this, we use the abstract definition above.

Si on s'intéresse à une fonction monotone, on déduit facilement les formules permettant de les calculer.If one is interested in a monotonous function, one can easily deduce the formulas for calculating them.

Par contre, on ne sait définir les fonctions élémentaires que sur des intervalles contenus dans leur domaine de définition : par exemple, le logarithme ne sera défini que pour des intervalles strictement positifs.On the other hand, it is not possible to define the elementary functions only on intervals contained in their domain of definition: for example, the logarithm will be defined only for strictly positive intervals.

L'arithmétique par intervalles permet de calculer avec des ensembles et d'obtenir des informations générales et précieuses pour l'optimisation globale d'une fonction.Interval arithmetic makes it possible to compute with sets and to obtain general and valuable information for the overall optimization of a function.

Pour éviter qu'il y ait une surestimation des résultats, il est préférable d'utiliser pour la fonction à prendre en compte une expression dans laquelle chaque variable n'apparaît qu'une fois.To avoid overestimating the results, it is preferable to use for the function to take into account an expression in which each variable appears only once.

Différentes méthodes de séparation des variables et d'évaluation (B&B) utilisant l'arithmétique des intervalles seront décrites ci-dessous.Different methods of separating variables and evaluation (B & B) using interval arithmetic will be described below.

Une méthode B&B peut être caractérisée en 5 étapes :

1. sélection : choix du noeud à examiner,
2. évaluation des bornes (bounding),
3. élimination : destruction des noeuds ne pouvant pas contenir l'optimum,
4. séparation : construction de 2 noeuds fils en divisant le domaine de variation d'une variable,
5. critère d'arrêt.

A B & B method can be characterized in 5 steps:

1. selection: choice of the node to examine,
2. Boundary evaluation,
3. elimination: destruction of nodes that can not contain the optimum,
4. separation: construction of 2 child nodes by dividing the variation domain of a variable,
5. stopping criterion.

Différentes solutions peuvent être choisies pour ces 5 étapes afin d'améliorer la qualité du procédé.Different solutions can be chosen for these 5 steps in order to improve the quality of the process.

Soit le problème d'optimisation min_x∈x f(X). Le vecteur d'intervalles de dimension n, X ∈ IR ⁿ, est la zone de recherche. La fonction f : Rⁿ→R est la fonction objectif.That is the optimization problem min _{x∈ x} f (X). The vector of intervals of dimension n, X ∈ IR ⁿ , is the search area. The function f : R ⁿ → R is the objective function.

On note f* le minimum global du problème, X* un point optimal tel que f(X*) = f*, et l'ensemble de ces points X* : $f^{*} = \min_{X \in X} f (X) et X^{*} = \{X \in X | f (X) = f^{*}\}$

We denote f * the global minimum of the problem, X * an optimal point such that f (X *) = f *, and all these points X *:

f^{*} = \min_{X \in X} f (X) and X^{*} = \{X \in X | f (X) = f^{*}\}

On note les objets intervalles par des caractères gras : x. On note x le minimum de x et x son maximum. On a alors x = [x , x ], et on considère l'ordre partiel sur IRⁿ : $X \leq Y \Leftrightarrow x_{i} \leq y_{i} pour i = 1.. n .$

We mark the interval objects with bold characters: x. We denote x the minimum of x and x its maximum. We then have x = [ x , x ], and we consider the partial order on IR ⁿ :

X \leq Y \Leftrightarrow x_{i} \leq {there}_{i} for i = 1 .. not .

We denote w ( x ) the width of x , (with w for width) or its diameter:

w (x) = \tilde{x} - \underset{}{x}

Le centre mid(x) et son rayon rad(x) sont définis par : $mid (x) = \frac{\overline{x} + \underset{̲}{x}}{2}$

rad (x) = \frac{\overline{x} - \underset{̲}{x}}{2} = \frac{w (x)}{2}

The center mid ( x ) and its radius rad ( x ) are defined by:

mid (x) = \frac{\tilde{x} + \underset{}{x}}{2}

rad (x) = \frac{\tilde{x} - \underset{}{x}}{2} = \frac{w (x) ()}{2}

Une fonction F : IRⁿ → IR est une fonction d'inclusion de ƒ sur X ∈ lR ⁿ. Si X ∈ X alors ƒ(X) ∈ F(X).An F: IR function ⁿ → IR is an inclusion function of ƒ on X ∈ lR ⁿ . If X ∈ X then ƒ (X) ∈ F (X) .

Voici différentes règles de sélection du noeud à examiner parmi la liste de noeuds en attente. Bien sûr, ces stratégies peuvent être combinées : par exemple la stratégie « Best-first » est souvent combinée avec la stratégie « Oldest-first » comme second critère s'il y a des ex-æquo.Here are different rules for selecting the node to examine from the list of pending nodes. Of course, these strategies can be combined: for example the "Best-first" strategy is often combined with the "Oldest-first" strategy as a second criterion if there are ties.

1. Plus ancien en premier (« Oldest-first »)1. Older first ("Oldest-first")

Cette stratégie consiste à examiner en premier le noeud qui a été créé le plus tôt.This strategy is to look first at the node that was created earlier.

2. Plus profond en premier (« Depth-First »)2. Deeper first ("Depth-First")

Cette stratégie consiste à examiner en premier le noeud qui se trouve au niveau le plus profond de l'arbre, i.e. le noeud qui a le plus d'ascendants.This strategy is to first examine the node at the deepest level of the tree, i.e. the node with the most ancestors.

3. Meilleur en premier (« Best-First ») [Règle de Moore-Skelboe]3. Best First ("Best-First") [Moore-Skelboe Rule]

Cette stratégie consiste à privilégier le noeud qui correspond au plus petit F(X), i.e. celui qui a la plus faible borne inférieure de l'optimum.This strategy consists of privileging the node that corresponds to the smallest F (X), i.e. the one with the lowest lower bound of the optimum.

4. Index de rejet (« Reject Index »)4. Reject Index (Reject Index) a. Optimum connuat. Optimum known

Définissons pour chaque noeud correspondant au vecteur d'intervalles X le paramètre : ${pf}^{*} (X) = \frac{f^{*} - \underset{̲}{F (X)}}{w (F (X))}$

Define for each node corresponding to the interval vector X the parameter:

{mp}^{*} (X) = \frac{f^{*} - \underset{}{F (X)}}{w (F (X))}

Remarquons que si w(F(X)) est nul, alors il n'y a pas lieu d'évaluer pƒ* puisque le noeud ne sera pas découpé.Note that if w (F (X)) is zero, then there is no need to evaluate pƒ * since the node will not be split.

Le noeud sélectionné est alors celui qui correspond à la plus forte valeur de pf*. Le calcul de ce paramètre nécessite cependant de connaître l'optimum à l'avance, ce qui n'est pas toujours le cas. C'est pourquoi des variantes du « Reject Index » s'appuyant sur des estimations de l'optimum ont été développées.The selected node is then the one that corresponds to the highest value of pf *. The calculation of this parameter, however, requires knowing the optimum in advance, which is not always the case. This is why variants of the "Reject Index" based on estimates of the optimum have been developed.

b. Optimum estiméb. Optimum estimated

La variante du paramètre pƒ* quand on ne connaît pas l'optimum à l'avance peut s'écrire : ${pf}^{*} (f_{k} X) = \frac{f_{k} - \underset{̲}{F (X)}}{w (F (X))}$

où k est l'indice de l'itération considérée. L'indice k correspond globalement au nombre de noeuds examinés et ƒ_k est une approximation de ƒ* à l'itération k.The variant of the parameter p ƒ * when we do not know the optimum in advance can be written:

{mp}^{*} (f_{k} X) = \frac{f_{k} - \underset{}{F (X)}}{w (F (X))}

where k is the index of the iteration considered. The index k corresponds globally to the number of nodes examined and ƒ _k is an approximation of ƒ * to the iteration k.

Remarquons que la règle « Best First » n'est donc jamais qu'un cas particulier de pƒ pour lequel ƒ_k = ƒ_k. En effet, si Y ₀ est l'intervalle du noeud présentant la plus faible borne inférieure de F (« meilleur noeud »), on a alors pƒ(Y ₀) = 0 et pƒ négatif pour tous les autres noeuds.Note that the "Best First" rule is therefore only a special case of p ƒ for which ƒ _k = ƒ _k . Indeed, if Y ₀ is the interval of the node with the lowest lower bound of F ("best node"), then p ƒ ( Y ₀ ) = 0 and p ƒ negative for all the other nodes.

D'autres possibilités pour f _k peuvent être : $f_{k} = \frac{\underset{̲}{F_{k}} + \overline{F_{k}}}{2}$

ou encore :

f_{k} = \underset{̲}{F_{k}}

Other possibilities for f _k can be:

f_{k} = \frac{\underset{}{F_{k}} + \tilde{F_{k}}}{2}

or :

f_{k} = \underset{}{F_{k}}

c. Avec contraintesvs. With constraints

Pour un problème contraint de la forme : ${\begin{cases} \min f (X) \\ C_{i} (X) \leq 0, i = 1.. p \\ X \in R^{n} \end{cases}$

For a constrained problem of the form:

{\begin{cases} \min f (X) \\ {VS}_{i} (X) \leq 0, i = 1 .. p \\ X \in R^{not} \end{cases}

les « Reject Index » définis ci-dessus ne tiennent pas compte du tout des contraintes et risquent de sélectionner des noeuds qui présentent de bonnes valeurs de pf mais conduiront à des noeuds infaisables.the Reject Indexes defined above do not take into account constraints at all and may select nodes that have good pf values but will lead to unfeasible nodes.

Certains auteurs proposent donc de définir un index de faisabilité construit de la façon suivante.Some authors propose to define a feasibility index built in the following way.

Pour une contrainte C_i et pour un noeud correspondant à un domaine de variation X, définissons : ${pu}_{C_{i}} (X) = \min (\frac{- \underset{̲}{C_{i} (X)}}{w (C_{i} (X))} 1)$

For a constraint C _i and for a node corresponding to a domain of variation X, let us define:

{could}_{{VS}_{i}} (X) = \min (\frac{- \underset{}{{VS}_{i} (X)}}{w ({VS}_{i} (X))} 1)

Dans le cas où w(C _i(X)) = 0, la faisabilité de la contrainte i peut être décidée directement, et pu _Ci(X) peut être fixé à 1 si X satisfait C_i, -1 sinon. Remarquons que si pu _ci(X) < 0, alors X, de manière certaine, ne satisfait pas C_i car C _i (X) > 0. Au contraire, si pu_Ci (X) = 1, alors C _i(X) ≤ 0 et donc X satisfait C_i de manière certaine. Dans tous les autres cas, l'état de transgression de C_i est indéterminé.In the case where w (C _i (X)) = 0, the feasibility constraint i can be determined directly, and could _Ci (X) may be set to 1 if X satisfies C _i, -1 otherwise. Note that if _one could (X) <0, then X, in some way, does not satisfy C _i C _i as (X)> 0. Conversely, if pu _Ci ( X ) = 1, then C _i ( X ) ≤ 0 and therefore X satisfies C _i with certainty. In all other cases, the state of transgression of C _i is undetermined.

Pour les X qui ne sont pas « certainement infaisables », c'est-à-dire pour lesquels ∀ i= 1..p, pu _ci(X) ≥ 0, définissons un index de faisabilité global pour l'ensemble des p contraintes : $pu (X) = \prod_{i = 1}^{p} {pu}_{C_{i}} (X)$

For X that are not "certainly infeasible", that is to say for which ∀ i = 1..p, could _this (X) ≥ 0, define a comprehensive feasibility index for the set of p :

could (X) = Π_{i = 1}^{p} {could}_{{VS}_{i}} (X)

Ainsi construit, cet index global possède 2 propriétés :

pu(X) = 1 ⇔ X est « certainement faisable », i.e. faisable de manière certaine
pu(X) ∈ [0,1] ⇔ X est indéterminé.

Thus constructed, this global index has 2 properties:

pu ( X ) = 1 ⇔ X is "certainly feasible", ie feasible with certainty
pu ( X ) ∈ [0,1] ⇔ X is undetermined.

Cela permet alors de définir un index de rejet modifié intégrant l'index de faisabilité : $pupf (f_{k} X) = pu (X) \times pf (f_{k} X)$

This then makes it possible to define a modified rejection index integrating the feasibility index:

pupf (f_{k} X) = could (X) \times mp (f_{k} X)

Si pu(X) = 1, i.e. si X est « certainement faisable », alors on se ramène au simple « Reject Index ». En revanche, si X est indéterminé, ce nouvel index prend en compte le degré de faisabilité de X. Cela permet de définir une nouvelle règle de sélection de noeud : on sélectionne le noeud avec la plus grande valeur de pupf.If pu (X) = 1, ie if X is "certainly feasible", then we can go back to the simple "Reject Index". On the other hand, if X is undetermined, this new index takes into account the degree of feasibility of X. This allows you to define a new node selection rule: you select the node with the largest pupf value.

Un dernier critère permet d'hybrider le critère pupf avec le critère « best-first » classique basé sur la valeur de F(X) : ${pupfb}^{*} (f_{k} X) = {\begin{cases} \frac{\underset{̲}{F (X)}}{{pupf}^{*} (f_{k} X)} si {pupf}^{*} (f_{k} X) \neq 0 \\ M si {pupf}^{*} (f_{k} X) = 0 \end{cases}$

avec M une valeur très grande fixée préalablement.A last criterion makes it possible to hybridize the criterion pupf with the criterion "best-first" classic based on the value of F (X):

{pupfb}^{*} (f_{k} X) = {\begin{cases} \frac{\underset{}{F (X)}}{{pupf}^{*} (f_{k} X)} if {pupf}^{*} (f_{k} X) \neq 0 \\ M if {pupf}^{*} (f_{k} X) = 0 \end{cases}

with M a very large value fixed beforehand.

En effet si pupf(f _k,X) = 0, alors soit pf(f _k,X) = 0, ce qui signifie - dans le cas où f _k =f- que l'on est sûr de ne pas améliorer f; soit pu(f _k,X) = 0 ce qui signifie qu'il existe au moins une contrainte telle que C _i (X) = 0. De tels X ne semblent pas très prometteurs. C'est pourquoi on fixe M à une valeur très grande.Indeed if pupf ( f _k , X ) = 0, then either pf ( f _k , X ) = 0, which means - in the case where f _k = f- that we are sure not to improve f ; let pu ( f _k , X ) = 0 which means that there exists at least one constraint such that C _i ( X ) = 0. Such X does not seem very promising. That is why we fix M to a very large value.

On considérera maintenant l'étape d'évaluation.We will now consider the evaluation step.

Dans cette étape, il s'agit d'évaluer les bornes de la fonction objectif, mais aussi celles des contraintes s'il y en a. Pour les méthodes de B&B utilisant l'arithmétique d'intervalles, les fonctions d'inclusion sont généralement obtenues par extension « naturelle » des fonctions usuelles.In this step, it is a question of evaluating the limits of the objective function, but also those of the constraints if there are some. For B & B methods using interval arithmetic, the inclusion functions are usually obtained by "natural" extension of the usual functions.

Exemple :Example:

Si f : x → x² - e^x et x = [-5,2], alors F : x → x ² - e^x est une fonction d'inctusion de f sur x avec : $x^{2} = {[\underset{̲}{x} \overline{x}]}^{2} = {\begin{cases} [\min ({\underset{̲}{x}}^{2} {\overline{x}}^{2}), \max ({\underset{̲}{x}}^{2} {\overline{x}}^{2})] si 0 \notin [\underset{̲}{x} \overline{x}] \\ [\begin{matrix} 0, \max ({\underset{̲}{x}}^{2} {\overline{x}}^{2}) \end{matrix}] sinon \end{cases}$

et

e^{x} = e^{[\underset{̲}{x} \overline{x}]} = [e^{\underset{̲}{x}} e^{\overline{x}}]

If f : x → x ² - e ^x and x = [-5,2], then F : x → x ² - e ^x is an inctusion function of f on x with:

x^{2} = {[\underset{}{x} \tilde{x}]}^{2} = {\begin{cases} [\min ({\underset{}{x}}^{2} {\tilde{x}}^{2}), \max ({\underset{}{x}}^{2} {\tilde{x}}^{2})] if 0 \notin [\underset{}{x} \tilde{x}] \\ [\begin{matrix} 0, \max ({\underset{}{x}}^{2} {\tilde{x}}^{2}) \end{matrix}] if not \end{cases}

and

e^{x} = e^{[\underset{}{x} \tilde{x}]} = [e^{\underset{}{x}} e^{\tilde{x}}]

Pour l'étape d'élimination, plusieurs méthodes sont possibles.For the elimination step, several methods are possible.

1. Test de faisabilité1. Feasibility test

Si le problème est un problème soumis à p contraintes d'inégalité C_i : ${\begin{cases} \min f (X) \\ C_{i} (X) \leq 0, i = 1.. p \\ X \in R^{n} \end{cases}$

If the problem is a problem subject to p inequality constraints C _i :

{\begin{cases} \min f (X) \\ {VS}_{i} (X) \leq 0, i = 1 .. p \\ X \in R^{not} \end{cases}

Soit C _i une fonction d'inclusion de la contrainte C_i. A chaque examen d'un noeud correspondant au domaine de variation de X, les p contraintes C _i(X) sont évaluées. Si 3 i ∈ {1,p}/[-∞,0] ∩ C _i(X) = ∅, alors il est sûr que le noeud ne peut contenir de solution faisable. Il peut donc être élagué.Let C _{i be} an inclusion function of the constraint C _i . At each examination of a node corresponding to the domain of variation of X, the p constraints C _i ( X ) are evaluated. If 3 i ∈ {1, p} / [- ∞, 0] ∩ C _i ( X ) = ∅, then it is certain that the node can not contain a feasible solution. It can therefore be pruned.

2. Test de coupure (« Cut Off Test »)2. Cut Off Test

C'est le critère le plus simple et le plus connu d'élimination : il s'agit de rejeter tous les noeuds pour lesquels f^* ≤ f < F(X), où f est la borne supérieure courante de l'optimum.This is the simplest and best known criterion of elimination: it is to reject all the nodes for which f ^* ≤ f < F ( X ) , where f is the current upper bound of the optimum.

3. Test de point milieu (« Middle Point Test »)3. Middle Point Test

Certaines publications ne font pas de distinction entre le « Cut Off Test » et le « Middle Point Test » (MPT). Le MPT ne serait en fait qu'un moyen supplémentaire de calculer une borne supérieure de f*. Le « Cut Off Test » consiste à prendre initialement comme borne supérieure F(X) puis de l'actualiser à chaque division d'intervalle. Pour un problème contraint, l'actualisation n'est possible que lorsqu'on sait que X contient au moins un point faisable. Dans le MPT, on prend f(mid(X)) qui est aussi une borne supérieure de l'optimum. S'il s'agit d'un problème contraint, il est toutefois nécessaire de s'assurer que mid(X) est un point faisable.Some publications do not distinguish between the "Cut Off Test" and the "Middle Point Test" (MPT). The MPT would be just another way of calculating an upper bound of f *. The "Cut Off Test" consists of initially taking as the upper limit F (X) then update it at each interval division. For a constrained problem, the update is possible only when one knows that X contains at least one feasible point. In the MPT, we take f (mid ( X )) which is also an upper bound of the optimum. If it is a constrained problem, it is however necessary to make sure that mid (X) is a feasible point.

4. Test de monotonie4. Monotonicity test

Pour un problème non contraint, si la fonction objectif est strictement monotone par rapport à la composante x _i d'un vecteur d'intervalle X, alors l'optimum ne peut se trouver à l'intérieur de x _i. Pour déterminer si f est strictement monotone par rapport aux composantes de X, on évalue les n composantes de la fonction d'inclusion du gradient de f sur X. Si pour i, l'intervalle résultant ne contient pas la valeur 0, alors f est strictement monotone par rapport à x_i.For a non-constrained problem, if the objective function is strictly monotonic with respect to the component x _i of an interval vector X , then the optimum can not lie within x _i . To determine if f is strictly monotonic with respect to components of X , we evaluate the n components of the inclusion function of the gradient of f on X. If for i, the resulting interval does not contain the value 0, then f is strictly monotonic with respect to x _i .

Dans ce cas, on peut réduire la composante x _i à un réel : x _i se réduit à x _i si la i^ème composante de la fonction d'inclusion du gradient est un intervalle qui a une borne supérieure strictement négative, et x _i se réduit à x _i si la i^ème composante de la fonction d'inclusion du gradient est un intervalle qui a une borne inférieure strictement positive.In this case, we can reduce the component x _i to a real: x _i is reduced to x _i if the i ^th component of the inclusion function of the gradient is an interval that has a strictly negative upper bound, and x _i is reduced to x _i if the i ^th component of the inclusion function of the gradient is an interval that has a strictly positive lower bound.

Pour l'étape de séparation, plusieurs méthodes sont également envisageables :For the separation step, several methods are also possible:

1. Bissection sur une variable1. Bissection on a variable

Dans toutes les règles suivantes, on sélectionne la variable j qui maximise une fonction de mérite D. On sépare donc sur la variable j telle que j=arg(max_i=1.nD(i)).In all the following rules, one selects the variable j which maximizes a function of merit D. One thus separates on the variable j such that j = arg (max _{i = 1.n} D (i)).

a. Plus grand diamètreat. Larger diameter

La fonction de mérite est ici simplement le diamètre de la variable : D(i)=ω(x_i ). La difficulté d'utiliser cette fonction de mérite est liée à la nécessité de s'affranchir des facteurs d'échelles. Par exemple, si l'on traite un problème de calcul de réseau, il faudra bien échelonner les variables pour pouvoir comparer les diamètres des pressions avec ceux des variables binaires.The merit function is here simply the diameter of the variable: D (i) = ω ( x _i ). The difficulty of using this merit function is related to the need to overcome scale factors. For example, if you are dealing with a network computation problem, you will have to scale the variables in order to compare the pressure diameters with those of the binary variables.

Pour pouvoir contourner cet obstacle, une règle proche de celle-ci et ne faisant pas non plus intervenir d'information sur les dérivées peut être définie : $D (i) = {\begin{matrix} w (x_{i}) si 0 \in x_{i} \\ \frac{w (x_{i})}{mig (x_{i})} si 0 \notin x_{i} \end{matrix}$

avec mig(X) = min_x∈Xi|X|. On pourrait utiliser la magnitude : mag(X) = max_x∈Xi|X|.To be able to circumvent this obstacle, a rule close to this one and not making use of information on the derivatives can be defined:

D (i) = {\begin{matrix} w (x_{i}) if 0 \in x_{i} \\ \frac{w (x_{i})}{mig (x_{i})} if 0 \notin x_{i} \end{matrix}

with mig ( X ) = min _{x∈ X i} | X | We could use the magnitude: mag ( X ) = max _{x∈ X i} | X |

Cette variante permet ainsi de normaliser le diamètre des intervalles considérés.This variant thus makes it possible to standardize the diameter of the intervals considered.

b. Règle de Hansenb. Hansen Rule

Ici, $D (i) = w (x_{i}) \times (▿ F_{i} (X))$

où ∇F _i est la i^ème composante de la fonction d'inclusion du gradient de f. L'idée est de séparer sur la variable qui a le plus d'impact sur f. Right here,

D (i) = w (x_{i}) \times (▿ F_{i} (X))

where ∇ F _i is the ^ith component of the inclusion function of the gradient of f. The idea is to separate on the variable that has the most impact on f.

c. Règle de Ratzvs. Rule of Ratz

Ici, $D (i) = w [(x_{i} - mid (x_{i})) \times ▿ F_{i} (X)]$

Right here,

D (i) = w [(x_{i} - mid (x_{i})) \times ▿ F_{i} (X)]

L'idée sous-jacente est ici de réduire le diamètre de w(F(X)) qui, après calcul, se ramène à la somme sur toutes les directions du terme D(i).The underlying idea here is to reduce the diameter of w (F (X)) which, after calculation, is reduced to the sum over all the directions of the term D (i).

d. Règle de Ratz bisd. Ruler of Ratz bis

L'idée sous-jacente est la même, mais on va jusqu'au second ordre : $D (i) = w [(x_{i} - mid (x_{i})) \times (▿ f_{i} (mid (x_{i})) + \frac{1}{2} \sum_{k = 1}^{n} H_{ik} (x_{i} - mid (x_{i})))]$

où H_ik est l'élément de coordonnées (i,k) de la matrice des dérivées secondes (hessien) de f.The underlying idea is the same, but we go to the second order:

D (i) = w [(x_{i} - mid (x_{i})) \times (▿ f_{i} (mid (x_{i})) + \frac{1}{2} Σ_{k = 1}^{not} H_{ik} (x_{i} - mid (x_{i})))]

where H _ik is the coordinate element (i, k) of the second derivative matrix (hessian) of f .

Pour des méthodes qui calculent de toute façon le gradient et le hessien, par différentiation automatique, cette règle n'est pas beaucoup plus coûteuse que les autres.For methods that calculate the gradient and the hessian anyway, by automatic differentiation, this rule is not much more expensive than the others.

2. Multi-section2. Multi-section a. Multi-section statiqueat. Static multi-section

Jusqu'ici, nous avons considéré qu'à partir d'un noeud, 2 noeuds fils étaient créés grâce à la bissection du pavé X ∈ IR ⁿ dans une direction unique. Cependant, il peut être pertinent de retenir plusieurs directions de séparation. Par exemple, on peut couper l'intervalle de variation de chaque variable en 2, 2ⁿ noeuds fils sont alors créés. On peut aussi couper l'intervalle d'une direction en 3 parties, créant ainsi 3 noeuds fils, ou encore les intervalles de 2 variables en 3, créant 3² fils, etc.So far, we have considered that from a node, 2 child nodes were created by the bisection of the block X ∈ IR ⁿ in a single direction. However, it may be relevant to remember several directions of separation. For example, we can cut the range of variation of each variable in 2, 2 ⁿ son nodes are then created. We can also cut the interval of a direction in 3 parts, thus creating 3 son nodes, or the intervals of 2 variables in 3, creating 3 ² son, etc.

b. Multi-section adaptativeb. Adaptive multi-section

Notons (a) la règle de plus grand diamètre présentée en 1.a, (b) la règle qui sépare les intervalles de toutes les variables en 2, (c) la règle qui sépare les intervalles de toutes les variables en 3.Let us note (a) the rule of largest diameter presented in 1.a, (b) the rule which separates the intervals of all the variables in 2, (c) the rule which separates the intervals of all the variables in 3.

Une règle hybride (adaptative) utilisera 3 paramètres P₁, P₂ et pf pour déterminer quelle règle utiliser.A hybrid rule (adaptive) will use 3 parameters P ₁ , P ₂ and pf to determine which rule to use.

Les paramètres p₁ et p₂ sont deux seuils qui seront à ajuster. pf est le « Reject index » défini plus haut, et est une fonction du noeud considéré.The parameters p ₁ and p ₂ are two thresholds that will have to be adjusted. pf is the "Reject index" defined above, and is a function of the considered node.

Les noeuds qui auront un « Reject Index » pf < p₁ seront séparés suivant la règle (a), ceux tels que p₁ <pf< p₂ seront séparés suivant la règle (b) et ceux tels que pf > p₂ seront séparés selon la règle (c).The nodes which will have a "Reject Index" pf < p ₁ will be separated according to the rule (a), those such that p ₁ < pf <p ₂ will be separated according to the rule (b) and those such that pf> p ₂ will be separated according to rule (c).

Une telle règle peut tout à fait être définie à partir de variantes de pf, telles que pupf défini plus haut par exemple.Such a rule can quite well be defined from variants of pf , such as pupf defined above for example.

Différents critères d'arrêt peuvent être utilisés.Different stopping criteria can be used.

1. Diamètre de la zone de recherche 1. Diameter of the search area

Un critère d'arrêt peut être l'examen d'un noeud N tel que w(X) ≤ ε où X est l'intervalle de variations des variables pour N. Bien sûr, cela suppose un bon échelonnement des variables.A stopping criterion can be the examination of a node N such that w ( X ) ≤ ε where X is the interval of variations of the variables for N. Of course, this supposes a good staggering of the variables.

2. Diamètre de la fonction objectif 2. Diameter of the objective function

Un critère d'arrêt peut être l'examen d'un noeud N tel que w(F(X)) ≤ ε où X est l'intervalle de variations des variables pour N.A stopping criterion can be the examination of a node N such that w ( F ( X )) ≤ ε where X is the interval of variations of the variables for N.

3. Temps maximal d'exécution 3. Maximum execution time

Un critère d'arrêt complémentaire peut être un temps maximal d'exécution au-delà duquel on arrête l'algorithme, quels que soient les résultats obtenus. Un critère d'arrêt de ce type est nécessaire en complément éventuellement d'un autre pour éviter des explorations trop longues.A complementary stop criterion can be a maximum execution time beyond which the algorithm is stopped, whatever the results obtained. A stopping criterion of this type is necessary in addition to possibly another to avoid too long explorations.

On décrira maintenant en référence à la Figure 11 un exemple d'organigramme illustrant la méthode de B&B (séparation des variables et évaluation) et de propagation des contraintes appliquée dans un solveur pour rechercher une solution optimale et exacte dans le cadre de la configuration d'un réseau de transport de gaz.An exemplary flowchart illustrating the B & B (variable separation and evaluation) and constraint propagation method applied in a solver for an optimal and accurate solution in the context of the configuration of the present invention will now be described with reference to FIG. a gas transportation network.

Pour mettre en oeuvre cette technique, une bibliothèque d'intervalles est mise en place pour permettre la gestion des variables exprimées sous forme de nombres ou d'intervalles.To implement this technique, an interval library is set up to allow management of variables expressed as numbers or ranges.

Par ailleurs des procédures de différentiation automatique basées sur des arbres de calcul permettent de calculer les valeurs des dérivées premières et secondes à partir d'une expression mathématique.In addition, automatic differentiation procedures based on calculation trees make it possible to calculate the values of the first and second derivatives from a mathematical expression.

Des moyens sont également mis en oeuvre pour effectuer le calcul de développements de Taylor aux ordres 1 et 2.Means are also implemented to calculate Taylor's developments at orders 1 and 2.

Dans l'organigramme de la Figure 11, les étapes 201, 202 et 203 correspondent à des étapes globales du procédé de B&B, tandis que les étapes 204, 206, 208, 211, 212, 214 sont appliquées à chaque pas du procédé de B&B. Les références 205, 207, 209, 210 correspondent à des tests aboutissant à une réponse oui ou non qui permet le choix de la procédure à suivre.In the flowchart of FIG. 11, steps 201, 202 and 203 correspond to global steps of the B & B process, while steps 204, 206, 208, 211, 212, 214 are applied to each step of the B & B process. . The references 205, 207, 209, 210 correspond to tests leading to a yes or no answer that allows the choice of the procedure to be followed.

De façon plus particulière, l'étape 201 correspond au choix de la meilleure feuille de l'arbre à explorer. L'étape 202 consiste en une séparation en noeuds fils. L'étape 203 comprend une série d'opérations effectuées pour chaque noeud fils.In particular, step 201 corresponds to the choice of the best leaf of the tree to explore. Step 202 consists of a separation into child nodes. Step 203 includes a series of operations performed for each child node.

Ainsi, à l'étape 203, il y a d'abord passage à une étape 204 de calcul des bornes, puis un test d'élagage 205 est ensuite effectué. Si la réponse est oui, il y a retour à l'étape 203 pour traiter un autre noeud fils. Si la réponse au test 205 est non, il y a passage à une étape 206 de propagation/rétro-propagation telle que celle proposée par exemple par F. Messine.Thus, in step 203, there is first a step 204 for calculating the terminals, and then a pruning test 205 is then performed. If the answer is yes, there is a return to step 203 to process another child node. If the answer to the test 205 is no, there is a step 206 propagation / back propagation as proposed for example by F. Messina.

Après l'étape 206, on effectue un nouveau test d'élagage 207. Si la réponse est oui, il y a retour à l'étape 203, si en revanche la réponse est non, on peut passer directement à un autre test 210, mais selon un mode de réalisation préférentiel, on procède d'abord à l'étape 208 à une résolution du système d'optimalité de Fritz-John, qui sera décrite plus en détail plus loin. En sortie de l'étape 208, un nouveau test d'élagage 209 permet de revenir à l'étape 203 si la réponse est oui et de passer au test 210 si la réponse est non (absence d'élagage).After step 206, a new pruning test 207 is carried out. If the answer is yes, there is a return to step 203, if on the other hand the answer is no, it is possible to go directly to another test 210. but according to a preferred embodiment, step 208 is first carried out to a resolution of the Fritz-John optimality system, which will be described in more detail below. At the exit of step 208, a new pruning test 209 makes it possible to return to step 203 if the answer is yes and to proceed to test 210 if the answer is no (no pruning).

Le test 210 permet d'examiner si les variables discrètes sont toutes instanciées ou non.The test 210 makes it possible to examine whether the discrete variables are all instantiated or not.

Si toutes les variables discrètes ne sont pas toutes instanciées, on passe à une étape 211 de mise à jour éventuelle de la meilleure solution, puis à une étape 212 de calcul du mérite du noeud pour insertion dans la file des feuilles et on revient à l'étape 203 de calcul pour un autre noeud fils.If all the discrete variables are not all instantiated, proceed to a step 211 of possible update of the best solution, then to a step 212 of calculating the merit of the node for insertion in the file of the sheets and returns to the computation step 203 for another child node.

Si le test 210 permet de déterminer que toutes les variables discrètes sont instanciées, on peut passer à une étape 214 de mise à jour éventuelle de la meilleure solution et on revient à l'étape 203 de calcul pour un autre noeud fils, sans calcul de mérite ni sous-arbre.If the test 210 makes it possible to determine that all the discrete variables are instantiated, it is possible to proceed to a step 214 of possibly updating the best solution and we return to the calculation step 203 for another child node, without calculation of merit or subtree.

A titre de variante, si le test 210 permet de déterminer que toutes les variables discrètes sont instanciées, on peut d'abord passer à une étape 213 de mise en oeuvre d'un solveur non linéaire qui permet d'effectuer une optimisation non linéaire basée par exemple sur une méthode des points intérieurs.Alternatively, if the test 210 makes it possible to determine that all the discrete variables are instantiated, one can first go on to a step 213 of implementing a nonlinear solver which makes it possible to perform a nonlinear optimization based on for example on a method of interior points.

Après l'étape 213, il y a passage à l'étape 214 précédemment décrite.After step 213, there is passage to step 214 previously described.

L'exemple de la Figure 11, sans les étapes 208, 209 et 213 est à nouveau explicité ci-dessous.The example of Figure 11, without steps 208, 209 and 213 is again explained below.

On part d'une liste triée de noeuds à explorer (étape 201). Le tri est effectué selon un mérite calculé pour chaque noeud. On peut par exemple effectuer une exploration selon la méthode du meilleur en premier (« Best First ») qui a été mentionnée plus haut. Dans ce cas, un noeud est exploré en priorité lorsqu'il présente la plus faible borne min de la fonction objectif.We start from a sorted list of nodes to explore (step 201). Sorting is done according to a calculated merit for each node. For example, an exploration can be done using the best first method mentioned above. In this case, a node is first explored when it has the lowest min terminal of the objective function.

Plusieurs fois au cours du procédé, on effectue un test d'élagage (étapes 205, 207). Si le noeud ne peut pas améliorer la solution courante, il ne sera pas exploré plus avant.Several times during the process, a pruning test is performed (steps 205, 207). If the node can not improve the current solution, it will not be explored further.

Le principe du procédé B&B est de partager un noeud en noeuds fils (étape 202). A titre d'exemple, on choisit la loi de séparation suivante : on sépare l'intervalle de la variable du noeud courant qui a le plus grand diamètre (la plus grande différence entre la borne supérieure et la borne inférieure de son intervalle) en deux intervalles. On range alors ces deux nouveaux noeuds dans une liste de noeuds fils du noeud courant. Puis pour chaque noeud fils (étape 203), on évalue la fonction objectif, c'est-à-dire qu'on évalue les bornes de la fonction objectif à partir des intervalles des variables de ce noeud (étape 204).The principle of the B & B method is to share a node in child nodes (step 202). For example, we choose the following separation law: we separate the interval of the variable from the current node that has the largest diameter (the largest difference between the upper bound and the lower bound of its interval) in two intervals. We then put these two new nodes in a list of child nodes of the current node. Then for each child node (step 203), the objective function is evaluated, that is to say that evaluates the boundaries of the objective function from the intervals of the variables of that node (step 204).

L'algorithme résultant peut être par exemple le suivant

The resulting algorithm can be for example the following

A titre de variante, un noeud pourrait être séparé en plus de deux noeuds fils (multi-section, par exemple quadri-section).Alternatively, a node could be separated into more than two child nodes (multi-section, for example quadri-section).

On indiquera ci-dessous quelques compléments concernant l'étape 208 de résolution du système d'optimalité de Fritz-John qui peut apporter une réponse au problème d'actualisation de la borne max de l'optimum en permettant de statuer sur la faisabilité d'un noeud.We will indicate below some complements concerning the step 208 of solving the Fritz-John optimality system which can provide an answer to the problem of updating the max bound of the optimum by making it possible to decide on the feasibility of a knot.

Considérons le problème d'optimisation suivant : ${\begin{matrix} \begin{array}{l} \min f (X) \\ C_{I} (X) \leq 0, i = 1.. p \\ C_{E} (X) \leq 0, i = 1.. q \\ X \in R^{n} \end{array} \end{matrix}$

Consider the following optimization problem:

{\begin{matrix} \begin{array}{l} \min f (X) \\ {VS}_{I} (X) \leq 0, i = 1 .. p \\ {VS}_{E} (X) \leq 0, i = 1 .. q \\ X \in R^{not} \end{array} \end{matrix}

L'approche la plus naturelle pour résoudre ce problème d'optimisation est de considérer le système d'équations issu des conditions d'optimalité de Karush-Kuhn-Tucker (KKT). Cependant, ces conditions d'optimalité présentent l'inconvénient de produire un système d'équations dégénéré dès lors que certaines contraintes sont linéairement dépendantes en la solution. Pour obtenir une approche plus robuste, on utilise les conditions d'optimalité de Fritz-John présentées ci-dessous.The most natural approach to solving this optimization problem is to consider the system of equations derived from the Karush-Kuhn-Tucker (KKT) optimality conditions. However, these optimality conditions have the disadvantage of producing a system of degenerate equations as soon as certain constraints are linearly dependent in the solution. To obtain a more robust approach, the Fritz-John optimality conditions presented below are used.

Les conditions de Fritz-John énoncent qu'ii existe λ₀,..., λ_p et µ₁,... µ_q qui vérifient le système d'optimalité suivant : ${\begin{matrix} λ_{0} ▿ f (X) + \sum_{i = 1}^{p} λ_{i} ▿ C_{I}^{i} (X) + \sum_{j = 1}^{p} μ_{i} ▿ C_{E}^{j} (X) & = 0 \\ λ_{i} C_{I}^{i} (X) & = 0, i & = 1.. p \\ C_{E}^{j} (X) & = 0, j & = 1.. q \\ λ_{i} & \geq 0, i & = 1.. p \end{matrix}$

The conditions of Fritz-John state that there exists λ ₀ , ..., λ _p and μ ₁ , ... μ _q which verify the following optimality system:

{\begin{matrix} λ_{0} ▿ f (X) + Σ_{i = 1}^{p} λ_{i} ▿ {VS}_{I}^{i} (X) + Σ_{j = 1}^{p} μ_{i} ▿ {VS}_{E}^{j} (X) & = 0 \\ λ_{i} {VS}_{I}^{i} (X) & = 0, i & = 1 .. p \\ {VS}_{E}^{j} (X) & = 0, j & = 1 .. q \\ λ_{i} & \geq 0, i & = 1 .. p \end{matrix}

Remarquons que les multiplicateurs µ_i peuvent être positifs ou négatifs tandis que les multiplicateurs λ_i sont exclusivement positifs.Note that the multipliers μ _i can be positive or negative while the multipliers λ _i are exclusively positive.

Une première différence entre les conditions de KKT et celles de Fritz-John réside dans ce que ces dernières introduisent le multiplicateur de Lagrange λ₀ ≠ 1.A first difference between the conditions of KKT and those of Fritz-John is that they introduce the Lagrange multiplier λ ₀ ≠ 1.

Une deuxième différence concernant toujours les multiplicateurs de Lagrange est que, pour les conditions de Fritz-John, les multiplicateurs λ_i et µ_j et peuvent être initialisés, respectivement, avec les intervalles [0,1] et [-1,1] tandis que, pour les conditions KKT, les multiplicateurs λ_i et µ_j sont initialisés, respectivement, avec les intervalles [0,+∝] et [-∝,+∝].A second difference always concerning the Lagrange multipliers is that, for the conditions of Fritz-John, the multipliers λ _i and μ _j and can be initialized, respectively, with the intervals [0,1] and [-1,1] while that for the conditions KKT, the multipliers λ _i and μ _j are initialized, respectively, with the intervals [0, + α] and [-α, + α].

Les conditions d'optimalité de Fritz-John n'incluent pas, à l'origine, de condition de normalisation. Dans ce cas, on peut remarquer qu'il y a (n + p + q + 1) variables et (n + p + q) équations, donc plus de variables que d'équations. Aussi, peut-on considérer la condition de normalisation suivante : $\begin{matrix} λ_{0} + \dots + λ_{p} + e_{1} μ_{1} + \dots + e_{q} μ_{q} = 1 où e_{j} = [1, 1 + ε_{0}], j = 1.. q \end{matrix}$

où ε₀ est le plus petit nombre tel que, selon la précision machine, 1 + ε₀ soit strictement supérieur à 1.
ou :

\begin{matrix} λ_{0} + \dots + λ_{p} + μ_{1}^{2} + \dots + μ_{q}^{2} = 1 \end{matrix}

The optimality conditions of Fritz-John do not initially include a normalization condition. In this case, we can notice that there are (n + p + q + 1) variables and (n + p + q) equations, so more variables than equations. Also, can we consider the following normalization condition:

\begin{matrix} λ_{0} + ... + λ_{p} + e_{1} μ_{1} + ... + e_{q} μ_{q} = 1 or e_{j} = [1, 1 + ε_{0}], j = 1 .. q \end{matrix}

where ε ₀ is the smallest number such that, according to the machine accuracy, 1 + ε ₀ is strictly greater than 1.
or :

\begin{matrix} λ_{0} + ... + λ_{p} + μ_{1}^{2} + ... + μ_{q}^{2} = 1 \end{matrix}

Dans le cas d'un problème d'optimisation en intervalles : ${\begin{matrix} \begin{array}{l} \min F (X) \\ C_{I} (X) \leq 0, i = 1.. p \\ C_{E} (X) \leq 0, i = 1.. q \\ X \in {IR}^{n} \end{array} \end{matrix}$

In the case of an interval optimization problem:

{\begin{matrix} \begin{array}{l} \min F (X) \\ {VS}_{I} (X) \leq 0, i = 1 .. p \\ {VS}_{E} (X) \leq 0, i = 1 .. q \\ X \in {IR}^{not} \end{array} \end{matrix}

C'est un Programme Sous Contraintes en Intervalles (ICSP).It is a Program Under Constraint in Intervals (ICSP).

Notons alors : $R_{1} (Λ M) = λ_{0} + \dots + λ_{p} + e_{1} μ_{1} + \dots + e_{q} μ_{q} - 1$

et

R_{2} (Λ M) = λ_{0} + \dots + λ_{p} + μ_{1}^{2} + \dots + μ_{q}^{2} - 1

où

Λ = {(λ_{0} \dots λ_{p})}^{T} et M = {(μ_{0} \dots μ_{q})}^{T}

Let's note then:

R_{1} (Λ M) = λ_{0} + ... + λ_{p} + e_{1} μ_{1} + ... + e_{q} μ_{q} - 1

and

R_{2} (Λ M) = λ_{0} + ... + λ_{p} + μ_{1}^{2} + ... + μ_{q}^{2} - 1

or

Λ = {(λ_{0} ... λ_{p})}^{T} and M = {(μ_{0} ... μ_{q})}^{T}

(CN1) s'écrit alors: $R_{1} (Λ M) = 0$

et (CN2) :

R_{2} (Λ M) = 0

(CN1) is written then:

R_{1} (Λ M) = 0

and (CN2):

R_{2} (Λ M) = 0

Pour résoudre le système des conditions d'optimalité de Fritz-John, on pose : $t = {(X Λ M)}^{T}$

et :

Φ (t) = (\begin{matrix} R_{k} (t) \\ λ_{0} ▿ f (X) + \sum_{i = 1}^{p} λ_{i} ▿ C_{I}^{i} (X) + \sum_{j = 1}^{q} μ_{i} ▿ C_{E}^{j} (X) \\ λ_{1} C_{I}^{i} (X) \\ ⋮ \\ λ_{p} C_{I}^{p} (X) \\ C_{E}^{1} (X) \\ ⋮ \\ C_{E}^{q} (X) \end{matrix}) où k = 1 ou 2

To solve the system of optimality conditions of Fritz-John, we pose:

t = {(X Λ M)}^{T}

and

Φ (t) = (\begin{matrix} R_{k} (t) \\ λ \\ _{0} ▿ f (X) + Σ_{i = 1}^{p} λ_{i} ▿ {VS}_{I}^{i} (X) + Σ_{j = 1}^{q} μ_{i} ▿ {VS}_{E}^{j} (X) \\ λ \\ _{1} {VS}_{I}^{i} (X) \\ ⋮ \\ λ_{p} {VS}_{I}^{p} (X) \\ {VS}_{E}^{1} (X) \\ ⋮ \\ {VS}_{E}^{q} (X) \end{matrix}) where k = 1 or 2

Notons t ^l une boîte de dimension N, où N = n + p + q + 1, contenant t. Soit J la jacobienne de Φ. Pour i, j = 1..N : $J_{ij} (t t^{I}) = \frac{\partial}{\partial t_{j}} Φ_{i} (T_{1} \dots T_{i} t_{j + 1} \dots t_{N})$

Let t ^{l be} a box of dimension N, where N = n + p + q + 1, containing t. Let J be the Jacobian of Φ. For i, j = 1..N:

J_{ij} (t t^{I}) = \frac{\partial}{\partial t_{j}} Φ_{i} (T_{1} ... T_{i} t_{j + 1} ... t_{NOT})

Les j premiers arguments de J _ij(t,t') sont des intervalles, les suivants sont des réels. En utilisant la normalisation linéaire (CN1), la jacobienne de Φ ne fera intervenir les multiplicateurs de Lagrange que sous forme de réels et non d'intervalles. Ainsi pour résoudre Φ(t) = 0, il n'est nul besoin d'initialiser un intervalle pour les multiplicateurs.The first j arguments of J _ij ( t , t ') are intervals, the following are reals. By using linear normalization (CN1), the Jacobian of Φ will only use Lagrange multipliers as real numbers, not intervals. So to solve Φ (t) = 0, there is no need to initialize an interval for the multipliers.

Utiliser (CN2) implique que les multiplicateurs de Lagrange apparaissent dans la jacobienne en intervalles et augmente les risques d'obtenir une matrice singulière. Une méthode de Newton peut alors soit échouer, soit être peu efficace. Dans ce cas, il faut envisager de découper les intervalles. Mais partager les intervalles des multiplicateurs implique, a priori, énormément de calculs supplémentaires.Using (CN2) implies that the Lagrange multipliers appear in the Jacobian at intervals and increases the chances of obtaining a singular matrix. A method of Newton can then either fail or be inefficient. In this case, it is necessary consider cutting the intervals. But to share the intervals of the multipliers implies, a priori, a lot of additional calculations.

D'où la préconisation d'utiliser (CN1) et l'ordre des variables de t comme indiqué ci-dessus. D'autant que (CN1) présente un caractère linéaire favorable.Hence the recommendation to use (CN1) and the order of the variables of t as indicated above. Especially since (CN1) has a favorable linear character.

En utilisant (CN1), certaines méthodes de Newton ne nécessitent pas d'initialiser un intervalle pour les multiplicateurs de Lagrange. Il peut pourtant être intéressant d'en disposer dans certains cas. En particulier, on peut avoir besoin d'une estimation des valeurs des multiplicateurs, ce qui est le cas dans le problème de calcul de réseau. Une telle estimation pour un multiplicateur peut être obtenue en retenant le milieu de son intervalle ; il faut donc disposer d'un encadrement. Pour le déterminer, on peut utiliser la méthode suivante :By using (CN1), some of Newton's methods do not need to initialize an interval for Lagrange multipliers. It may be interesting to have some in some cases. In particular, it may be necessary to estimate the values of the multipliers, which is the case in the network calculation problem. Such an estimate for a multiplier can be obtained by retaining the middle of its range; therefore, you need a framework. To determine it, we can use the following method:

On pose : $A (X) = [\begin{matrix} 1 & 1 & \dots & 1 & e_{1} & \dots & e_{q} \\ ▿ f (X) & ▿ C_{I}^{1} (X) & \dots & ▿ C_{I}^{p} (X) & ▿ C_{E}^{1} (X) & \dots & ▿ C_{E}^{q} (X) \end{matrix}]$

We pose:

AT (X) = [\begin{matrix} 1 & 1 & \dots & 1 & e_{1} & \dots & e_{q} \\ ▿ f (X) & ▿ {VS}_{I}^{1} (X) & \dots & ▿ {VS}_{I}^{p} (X) & ▿ {VS}_{E}^{1} (X) & \dots & ▿ {VS}_{E}^{q} (X) \end{matrix}]

Si on résout : $A (X) (\begin{matrix} Λ \\ M \end{matrix}) = (\begin{matrix} 1 \\ 0 \\ ⋮ \\ 0 \end{matrix})$

on obtiendra l'encadrement voulu pour les multiplicateurs de Lagrange.If we solve:

AT (X) (\begin{matrix} Λ \\ M \end{matrix}) = (\begin{matrix} 1 \\ 0 \\ ⋮ \\ 0 \end{matrix})

we will obtain the required framework for the Lagrange multipliers.

L'utilisation des conditions d'optimalités de Fritz-John au sein du solveur peut être de deux utilités. La première est qu'elles peuvent réduire davantage l'espace des solutions en complément ou en remplacement de la propagation des contraintes à partir d'un certain niveau de l'arbre de la méthode B&B. La seconde vient du fait que la résolution des conditions d'optimalité de Fritz-John est un opérateur de Newton. Peut alors s'appliquer le théorème de Moore-Nickel qui énonce que si un opérateur de Newton permet de réduire un intervalle de définition d'une variable au moins, alors l'espace des solutions courant contient forcément un optimum. Ainsi la résolution de ces conditions d'optimalité peut aussi être un critère de mise à jour de la borne max de l'optimum de la fonction objectif.The use of the Fritz-John optimality conditions within the solver can be of two uses. The first is that they can further reduce the solution space in addition to or as a replacement for constraint propagation from a certain level of the B & B tree. The second comes from the fact that the resolution of the optimality conditions of Fritz-John is a Newton operator. The Moore-Nickel theorem can then be applied, stating that if a Newton operator makes it possible to reduce a definition interval of at least one variable, then the space of the current solutions necessarily contains an optimum. Thus the resolution of these optimality conditions can also be a criterion for updating the max limit of the optimum of the objective function.

La résolution du système linéaire ci-dessus (SL) peut être effectuée, par exemple, par la méthode itérative de Gauss-Seidel (ou de propagation des contraintes) ou par la méthode LU.The resolution of the linear system above (SL) can be performed, for example, by the iterative Gauss-Seidel (or constraint propagation) method or by the LU method.

On se donne un système linéaire tel que celui posé par la linéarisation des conditions d'optimalité d'un problème d'optimisation, de la forme : $\begin{matrix} A . X + B = 0 \end{matrix}$

We give ourselves a linear system such as that posed by the linearization of the optimality conditions of an optimization problem, of the form:

\begin{matrix} AT . X + B = 0 \end{matrix}

A est une matrice m × n de réels ou d'intervalles, X est le vecteur des variables de dimension n, B est un vecteur de dimension m de réels ou d'intervalles. A is an m × n matrix of reals or intervals, X is the vector of variables of dimension n, B is a vector of dimension m of reals or intervals.

La méthode Gauss-Seidel est une méthode itérative faisant suite à une amélioration de la méthode de Jacobi.The Gauss-Seidel method is an iterative method following an improvement of the Jacobi method.

Une méthode itérative pour résoudre un système linéaire tel que (SL) consiste à construire une suite de vecteurs X _k qui converge vers la solution X*. Dans la pratique, les méthodes itératives sont rarement utilisées pour résoudre les systèmes linéaires de petites dimensions car, dans ce cas, elles sont généralement plus coûteuses que les méthodes directes. Toutefois, ces méthodes s'avèrent efficaces (en termes de coût) dans les cas où le système linéaire (SL) est de grande dimension et contient un grand nombre de coefficients nuls. On dit alors que la matrice A est creuse ; c'est le cas lors d'un calcul de réseau.An iterative method for solving a linear system such as (SL) consists of constructing a sequence of vectors X _k which converges to the solution X *. In practice, iterative methods are rarely used to solve linear systems of small dimensions because, in this case, they are usually more expensive than direct methods. However, these methods are efficient (in terms of cost) in cases where the linear system (SL) is large and contains a large number of zero coefficients. It is said then that the matrix A is hollow; this is the case during a network calculation.

La méthode itérative de Jacobi consiste à résoudre l'i^ème équation en fonction de X _i pour obtenir : $X_{i} = \frac{B_{i}}{A_{ii}} - \sum_{\begin{matrix} j = 1 \\ j \neq i \end{matrix}}^{n} \frac{A_{ij} \times X_{j}}{A_{ii}}$

The iterative method of Jacobi consists in solving the i ^th equation as a function of X _i to obtain:

X_{i} = \frac{B_{i}}{{AT}_{ii}} - Σ_{\begin{matrix} j = 1 \\ j \neq i \end{matrix}}^{not} \frac{{AT}_{ij} \times X_{j}}{{AT}_{ii}}

On construit le terme X ^k à partir des composants de X ^k-1 : $X_{i}^{k} = \frac{B_{i}}{A_{ii}} - \sum_{\begin{matrix} j = 1 \\ j \neq i \end{matrix}}^{n} \frac{A_{ij} \times X_{j}^{k - 1}}{A_{ii}}$

The term X ^{k is constructed} from the components of X ^k-1 :

X_{i}^{k} = \frac{B_{i}}{{AT}_{ii}} - Σ_{\begin{matrix} j = 1 \\ j \neq i \end{matrix}}^{not} \frac{{AT}_{ij} \times X_{j}^{k - 1}}{{AT}_{ii}}

Or, lors du calcul de X ^k, les composantes $X_{j}^{k}$

pour j < i sont connues. La méthode Gauss-Seidel substitue

X_{j}^{k}

par

X_{j}^{k - 1}

pour j < i.However, when calculating X ^k , the components

X_{j}^{k}

for j <i are known. The Gauss-Seidel method substitutes

X_{j}^{k}

by

X_{j}^{k - 1}

for j <i.

Dans le problème du calcul de réseau, les éléments de A, X et B sont des intervalles. L'algorithme est donc le suivant :

In the network computation problem, the elements of A, X and B are intervals. The algorithm is as follows:

La méthode LU décompose la matrice A du système (SL) selon le produit suivant : $A = L . U$

où L est une matrice triangulaire inférieure à diagonale unité :

L = (\begin{matrix} 1 & 0 & \dots & 0 \\ L_{21} & 1 & ⋱ & ⋮ \\ ⋮ & ⋱ & ⋱ & 0 \\ L_{n 1} & \dots & L_{nn - 1} & 1 \end{matrix})

et U est une matrice triangulaire supérieure :

U = (\begin{matrix} U_{11} & \dots & U_{1 n} \\ ⋮ & ⋱ & ⋮ \\ 0 & \dots & U_{nn} \end{matrix})

The LU method breaks down system A matrix (SL) according to the following product:

AT = The . U

where L is a triangular matrix lower than diagonal unit:

The = (\begin{matrix} 1 & 0 & \dots & 0 \\ {The}_{21} & 1 & ⋱ & ⋮ \\ ⋮ & ⋱ & ⋱ & 0 \\ {The}_{not 1} & \dots & {The}_{nn - 1} & 1 \end{matrix})

and U is an upper triangular matrix:

U = (\begin{matrix} U_{11} & \dots & U_{1 not} \\ ⋮ & ⋱ & ⋮ \\ 0 & \dots & U_{nn} \end{matrix})

Le système devient donc : $\begin{matrix} L . U . X = B \end{matrix}$

que l'on peut décomposer en deux systèmes :

{\begin{matrix} L . Y = B \\ U . X = Y \end{matrix}

The system becomes:

\begin{matrix} The . U . X = B \end{matrix}

that can be broken down into two systems:

{\begin{matrix} The . Y = B \\ U . X = Y \end{matrix}

La résolution de (SL1) puis de (SL2) est grandement facilitée par la forme triangulaire de L et U.The resolution of (SL1) then of (SL2) is greatly facilitated by the triangular shape of L and U.

La Figure 13 montre un exemple de réseau auquel est applicable le procédé d'optimisation automatique selon l'invention.Figure 13 shows an example of a network to which the automatic optimization method according to the invention is applicable.

Ce réseau comprend un ensemble de points d'interconnexion (jonctions ou noeuds) 1.1 à 1.13 qui permettent de relier entre elles des canalisations passives 101 à 112 ou des tronçons de canalisation comportant des ouvrages actifs tels que des vannes de régulation 31, 32, une station de compression 41, une vanne d'isolement 51, des consommations 61 à 65 ou des ressources 21, 22.This network comprises a set of interconnection points (junctions or nodes) 1.1 to 1.13 which make it possible to connect between them passive pipes 101 to 112 or pipe sections comprising active structures such as control valves 31, 32, a compression station 41, an isolation valve 51, consumptions 61 to 65 or resources 21, 22.

Des conduits de by-pass 31A, 32A, 41A sont associés aux vannes de régulation 31, 32 et à la station de compression 41.By-pass ducts 31A, 32A, 41A are associated with the control valves 31, 32 and with the compression station 41.

Claims

Procédé d'optimisation automatique d'un réseau de transport de gaz naturel en régime permanent, le réseau de transport de gaz naturel comprenant à la fois un ensemble d'ouvrages passifs tels que des canalisations (101 à 112) ou des résistances, et un ensemble d'ouvrages actifs comprenant des vannes de régulation (31, 32), des vannes d'isolement (51), des stations de compression (41) avec chacune au moins un compresseur, des dispositifs de stockage ou d'alimentation (21, 22), des dispositifs de consommation (61 à 65), des éléments (41A) de dérivation des stations de compression (41) et des éléments (31A, 32A) de dérivation des vannes de régulation (31, 32), les ouvrages passifs et les ouvrages actifs étant reliés entre eux par des jonctions (1.1 à 1.13), le procédé d'optimisation comprenant la détermination de valeurs pour des variables continues telles que la pression et le débit du gaz naturel en tout point du réseau de transport, et la détermination de valeurs pour des variables discrètes telles que l'état de démarrage des compresseurs, l'état d'ouverture des stations de compression, l'état d'ouverture des vannes de régulation, l'état des éléments de dérivation des stations de compression, l'état des éléments de dérivation des vannes de régulation, l'orientation des stations de compression et l'orientation des vannes de régulation,
caractérisé en ce que l'on choisit comme état initial de l'optimisation des intervalles de valeurs pour les variables continues et des ensembles de valeurs pour les variables discrètes, en ce que l'on explore les possibilités de valeurs pour les variables en construisant au fur et à mesure un arbre avec des branches reliées à des noeuds décrivant les combinaisons de valeurs envisagées en utilisant une technique de séparation des variables, c'est-à-dire de découpage conduisant à la génération de nouveaux noeuds dans l'arbre, et d'évaluation, c'est-à-dire de détermination avec une probabilité forte des branches de l'arbre qui peuvent conduire à des feuilles constituant une solution finale optimisée, de manière à parcourir en priorité ces branches à plus forte probabilité de réussite, les valeurs des grandeurs recherchées étant considérées comme optimales lorsque des contraintes prédéterminées ne sont plus transgressées ou sont transgressées au minimum et une fonction objectif est minimisée, cette fonction objectif étant de la forme

g = α \times Régime + β \times Energie + γ \times Cible

avec : α, β et γ sont des coefficients de pondération

Régime représente un facteur de minimisation ou de maximisation de la pression en des points donnés du réseau tels que tout point aval d'un dispositif de stockage ou d'alimentation, tout point amont et tout point aval d'une station de compression ou d'une vanne de régulation, et tout point amont d'un dispositif de consommation,

Energie représente un facteur de minimisation de la consommation d'énergie de compression,

Cible représente un facteur de maximisation ou de minimisation du débit d'un tronçon du réseau situé entre deux jonctions ou de la pression d'une jonction particulière, et lesdites contraintes prédéterminées comprenant d'une part des contraintes d'égalité comprenant la loi de perte de charge dans les canalisations et la loi des noeuds régissant le calcul des réseaux, et d'autre part des contraintes d'inégalités comprenant des contraintes de débit minimal et maximal, des contraintes de pression minimale et maximale des ouvrages actifs ou passifs, des contraintes de puissance de compression des stations de compression.

A method of automatically optimizing a steady-state natural gas transmission network, the natural gas transmission network comprising both a set of passive structures such as pipes (101 to 112) or resistors, and a set of active structures comprising control valves (31, 32), isolation valves (51), compressor stations (41) each with at least one compressor, storage or supply devices (21, 22), consumption devices (61 to 65), elements (41A) for bypassing the compressor stations (41) and elements (31A, 32A) for bypassing the control valves (31, 32), the passive structures and the active structures being interconnected by junctions (1.1 to 1.13), the optimization method comprising the determination of values for continuous variables such as the pressure and the flow rate of the natural gas at any point of the transport network, and the determination of value rs for discrete variables such as the start state of the compressors, the opening state of the compressor stations, the opening state of the control valves, the state of the bypass elements of the compressor stations, state of the bypass elements of the control valves, the orientation of the compressor stations and the orientation of the control valves,
characterized in that value ranges for continuous variables and sets of values for discrete variables are chosen as the initial state of optimization, in that value possibilities for the variables are explored by constructing as a tree with branches connected to nodes describing the combinations of values envisaged using a technique of separation of the variables, that is to say of cutting leading to the generation of new nodes in the tree, and evaluation, that is to say of determination with a strong probability of the branches of the tree which can lead to leaves constituting an optimized final solution, so as to go in priority to those branches with a higher probability of success, the values of the quantities sought are considered optimal when predetermined constraints are no longer transgressed or are transgressed to a minimum and an objective function is minimized, this objective function being of the form

boy Wut = α \times Diet + β \times Energy + γ \times Target

with: α, β and γ are weighting coefficients

Regime represents a factor of minimization or maximization of the pressure at given points of the network such as any downstream point of a storage or supply device, any upstream point and any downstream point of a compressor station or a control valve, and any upstream point of a consumption device,

Energy represents a factor of minimization of the energy consumption of compression,

Target represents a factor for maximizing or minimizing the bit rate of a section of the network located between two junctions or the pressure of a particular junction, and said predetermined constraints comprising on the one hand equality constraints comprising the law of loss. in the pipes and the law of the nodes governing the calculation of the networks, and on the other hand inequality constraints including minimum and maximum flow constraints, minimum and maximum pressure constraints of the active or passive structures, constraints compression power compression stations.

Procédé selon la revendication 1, caractérisé en ce que le problème de configuration optimale des ouvrages actifs est modélisé sous la forme d'un programme P₁ d'optimisation se présentant sous la forme suivante :

P_{1} {\begin{cases} \min_{(x s e)} f (x s) = g (x) + α \times {‖ s ‖}^{2} \\ C_{I} (x) + β . e \leq s_{I} \\ C_{E} (x) = s_{E} \\ x \in R^{n}, s_{I} \in R^{p}, s_{E} \in R^{q}, e \in {\{0 1\}}^{p} \end{cases}

avec : x est l'ensemble des variables des débits Q et des pressions P,

g(x) est la fonction objectif constituant le critère économique d'optimisation,

C₁(x) est l'ensemble des p contraintes d'inégalité linéaires et non linéaires sur les ouvrages actifs,

β est un vecteur dont les coefficients sont nuls ou égaux aux valeurs maximales des contraintes,

e est le vecteur des variables binaires,

C_E(x) est l'ensemble des q contraintes d'égalité linéaires ou non linéaires,

s est une variable d'écart qui, lorsqu'elle est non nulle, représente la transgression d'une contrainte,

α est un coefficient représentant le degré d'autorisation de transgression de contraintes.

Method according to Claim 1, characterized in that the problem of optimum configuration of the active structures is modeled in the form of an optimization program P ₁ in the following form:

P_{1} {\begin{cases} \min_{(x s e)} f (x s) = boy Wut (x) + α \times {‖ s ‖}^{2} \\ {VS}_{I} (x) + β . e \leq s_{I} \\ {VS}_{E} (x) = s_{E} \\ x \in R^{not}, s_{I} \in R^{p}, s_{E} \in R^{q}, e \in {\{0 1\}}^{p} \end{cases}

with: x is the set of variables of flow rates Q and pressures P,

g (x) is the objective function constituting the economic criterion of optimization,

C ₁ (x) is the set of p linear and nonlinear inequality constraints on active structures,

β is a vector whose coefficients are zero or equal to the maximum values of the stresses,

e is the vector of the binary variables,

C _E (x) is the set of q linear or nonlinear equality constraints,

s is a variance variable which, when nonzero, represents the transgression of a constraint,

α is a coefficient representing the degree of authorization of constraint transgression.

Procédé selon la revendication 1 ou la revendication 2, caractérisé en ce que les variables sont représentées par des intervalles, en ce que la technique de séparation des variables est appliquée aux seules variables discrètes et en ce que des bornes de la fonction objectif sont calculées en utilisant l'arithmétique des intervalles.Method according to Claim 1 or Claim 2, characterized in that the variables are represented by intervals, in that the variable separation technique is applied to the discrete variables only and that limits of the objective function are calculated in using the arithmetic of the intervals.

Procédé selon la revendication 1 ou la revendication 2, caractérisé en ce que les variables sont représentées par des intervalles, en ce que la technique de séparation des variables est appliquée à la fois aux variables discrètes et aux variables continues, la séparation comprenant le découpage de l'espace de définition des variables continues, l'exploration s'effectuant séparément sur des parties de l'ensemble réalisable et l'intervalle de variation de la fonction objectif étant évalué sur chacune de ces parties.Method according to Claim 1 or Claim 2, characterized in that the variables are represented by intervals, in that the variable separation technique is applied to both the discrete and continuous variables, the separation comprising the division of the variables. the space of definition of the continuous variables, the exploration being carried out separately on parts of the realizable set and the interval of variation of the objective function being evaluated on each of these parts.

Procédé selon la revendication 4, caractérisé en ce que lors de l'exploration des possibilités de valeurs pour les variables avec une technique de séparation des variables et d'évaluation, on établit d'abord une liste de noeuds à explorer triés selon un critère de mérite M calculé pour chaque noeud, tant que la liste de noeuds à explorer n'est pas vide, pour chaque noeud courant, on évalue si ce noeud courant peut contenir une solution, si c'est le cas, on découpe l'intervalle correspondant à la variable considérée selon une loi de séparation pour établir une liste de noeuds fils, pour chaque noeud fils on évalue des bornes minimale et maximale de la fonction objectif et on évalue si le noeud fils peut améliorer la situation courante, si c'est le cas on effectue une propagation de la contrainte sur ses variables, si la propagation ne conduit pas à des intervalles vides, on évalue des bornes minimale et maximale de la fonction objectif et on vérifie qu'il n'y a pas d'impossibilité à ce que le noeud fils contienne au moins une solution faisable, on effectue un test pour déterminer s'il y a encore des valeurs discrètes non instanciées, c'est-à-dire des variables pour lesquelles aucune valeur précise et définitive n'a pu être décidée, on met à jour la meilleure solution courante s'il y a lieu et on calcule le mérite du noeud pour l'insérer dans la liste des feuilles, tirée selon ce critère de mérite.Method according to Claim 4, characterized in that, when exploring the possible values for the variables with a variable separation and evaluation technique, a list of nodes to be searched is first sorted according to a criterion of merit M calculated for each node, as long as the list of nodes to explore is not not empty, for each current node, it is evaluated whether this current node can contain a solution, if it is the case, we cut the interval corresponding to the variable considered according to a separation law to establish a list of child nodes, for each child node evaluates the minimum and maximum limits of the objective function and it is evaluated whether the child node can improve the current situation, if it is the case one carries out a propagation of the constraint on its variables, if the propagation does not lead at empty intervals, minimum and maximum bounds of the objective function are evaluated and it is verified that it is not impossible for the child node to contain at least one feasible solution, a test is made to determine there are still uninstantiated discrete values, ie variables for which no definite and definitive value could be decided, the best current solution is updated if necessary and calcu the merit of the node to insert it in the list of sheets, drawn according to this criterion of merit.

Procédé selon la revendication 5, caractérisé en ce que le critère de mérite M est tel qu'un noeud est exploré en priorité lorsqu'il présente la plus faible borne minimale de la fonction objectif.Method according to claim 5, characterized in that the merit criterion M is such that a node is first explored when it has the lowest minimum bound of the objective function.

Procédé selon la revendication 5 ou la revendication 6, caractérisé en ce que lors des tests d'élimination des noeuds ne pouvant pas contenir l'optimum, on met en oeuvre l'une des méthodes consistant à utiliser la monotonie de la fonction objectif, à utiliser un test de contraintes transgressées ou à utiliser un test de valeur d'objectif moins bonne que la valeur courante.Method according to Claim 5 or Claim 6, characterized in that, during the node elimination tests which can not contain the optimum, one of the methods of using the monotonicity of the objective function is used. use a test of transgressed constraints or to use an objective value test that is not as good as the current value.

Procédé selon l'une quelconque des revendications 5 à 7, caractérisé en ce que lors de la séparation d'un noeud courant en noeuds fils, on divise le domaine de variation d'une ou plusieurs variables choisies selon des critères basés sur le diamètre d'intervalles rattachés aux variables.Process according to any one of Claims 5 to 7, characterized in that, during the separation of a current node into child nodes, the variation domain of one or more variables chosen is divided according to criteria based on the diameter of intervals attached to the variables.

Procédé selon l'une quelconque des revendications 5 à 8, caractérisé en ce qu'il comprend, en outre, un critère d'arrêt basé sur le temps d'exécution ou sur l'évaluation de certains diamètres d'intervalles.Method according to any one of claims 5 to 8, characterized in that it further comprises a stopping criterion based on the execution time or the evaluation of certain interval diameters.

Procédé selon l'une quelconque des revendications 5 à 10, caractérisé en ce qu'en complément de la propagation des contraintes, on procède à une actualisation de la borne maximale de l'optimum de la fonction objectif en utilisant les conditions d'optimalité du problème d'optimisation dites de Fritz-John.Method according to one of Claims 5 to 10, characterized in that, in addition to the propagation of the constraints, the maximum bound of the optimum of the objective function is updated by using the optimality conditions of the so-called Fritz-John optimization problem.

Procédé selon i'une quelconque des revendications 5 à 10, caractérisé en ce que lorsqu'à un noeud du procédé de séparation et d'évaluation, toutes les variables discrètes ont été instanciées, on met en oeuvre en outre un processus d'optimisation non linéaire basé sur une méthode de points intérieurs.Process according to any one of Claims 5 to 10, characterized in that when at a node of the separation and evaluation process all the discrete variables have been instantiated, a further optimization process is carried out. linear based on a method of inner points.

Procédé selon l'une quelconque des revendications 5 à 9, caractérisé en ce qu'à chaque noeud du procédé de séparation et d'évaluation, on met en outre en oeuvre un processus d'optimisation non linéaire basé sur une méthode de points intérieurs.Method according to one of Claims 5 to 9, characterized in that, at each node of the separation and evaluation method, a non-linear optimization process based on an internal point method is furthermore carried out.