Keywords: Neurodynamic Chance-Constrained Geometric Optimization Geometric optimization, Chance constraints, Recurrent neural networks, Lyapunov Theory, Ordinary differential equations

In many real-world scenarios, decision-makers face uncertainties that can affect the outcomes of their decisions. These uncertainties arise from various sources, such as variability in demand, fluctuating market conditions, or incomplete information about system parameters. Traditional deterministic optimization approaches assume that all parameters are known with certainty, which may not accurately reflect the reality of the problem. Chance-constrained optimization provides a more realistic and robust approach by explicitly accounting for the uncertainty in decision-making. Geometric programming is often misunderstood as a technique exclusively designed for posynomial problems. However, it is a versatile mathematical theory with significant value in addressing a broad range of separable problems. In fact, its true strength lies in its ability to effectively tackle seemingly inseparable problems by leveraging their linear algebraic structure. This general applicability of geometric programming makes it a valuable tool for studying and solving various optimization problems, extending its practical usefulness beyond its initial perception. Recurrent neural networks (RNNs) offer a biologically inspired computational framework with great optimization potential. By emulating the interconnected structure of neurons in the brain, RNNs excel in modeling complex and dynamic systems. This capability allows them to capture temporal dependencies and feedback loops,

Titre: Optimisation Géométrique Neurodynamique avec des Contraintes en Probabilité Mots clés: Optimisation géométrique, Contraintes en probabilité, Réseaux de neurones récurrents, Théorie de Lyapunov, Equations différentielles ordinaires.

Résumé: Dans de nombreux problèmes réels, les décideurs sont confrontés à des incertitudes qui peuvent affecter les résultats de leurs décisions. Ces incertitudes découlent de diverses sources, telles que la variabilité de la demande, les conditions fluctuantes du marché ou des informations incomplètes sur les paramètres du système. Les approches traditionnelles d'optimisation déterministe supposent que tous les paramètres sont connus avec certitude, ce qui peut ne pas refléter avec précision la réalité du problème. L'optimisation sous contraintes de probabilité offre une approche plus réaliste et robuste en tenant compte explicitement de l'incertitude dans la prise de décision. La programmation géométrique est souvent mal comprise comme une technique exclusivement conçue pour les problèmes posynômes. Cependant, c'est une théorie mathématique polyvalente qui a une valeur significative pour résoudre un large éventail de problèmes. En fait, sa véritable force réside dans sa capacité à résoudre efficacement des problèmes en apparence inséparables en exploitant leur structure algébrique linéaire. Cette applicabilité générale de la programmation géométrique en fait un outil précieux pour étudier et résoudre divers problèmes d'optimisation, étendant ainsi son utilité pratique au-delà de sa perception initiale. Les réseaux de neurones récurrents offrent un cadre de calcul inspiré de la biologie avec un grand potentiel d'optimisation. En imitant la structure interconnectée des neurones du cerveau, les réseaux de neurones récurrents excellent dans la modélisation de systèmes complexes et dynamiques.

Cette capacité leur permet de capturer les dépendances tem-porelles et les boucles de rétroaction, ce qui les rend bien adaptés aux scénarios d'optimisation impliquant des prises de décision séquentielles ou des processus itératifs. De plus, l'un des principaux avantages des approches neurodynamiques est leur faisabilité de mise en oeuvre matérielle. L'objectif principal de cette thèse est de développer des algorithmes neurodynamiques efficaces et performants pour résoudre des problèmes d'optimisation géométrique avec des contraintes de probabilité.

La thèse commence par les programmes géométriques avec des contraintes de probabilité impliquant des variables aléatoires indépendantes.

De plus, un type spécifique de programmes géométriques appelés programmes rectangulaires est également examiné en détail. L'objectif est de comprendre les caractéristiques et les complexités associées à cette sous-classe de programmes géométriques. Ensuite, la thèse explore l'application de la théorie des copules pour aborder les programmes géométriques avec des contraintes de probabilité impliquant des variables aléatoires dépendantes. La théorie des copules fournit un cadre mathématique pour modéliser et analyser la structure de dépendance entre les variables aléatoires, améliorant ainsi la compréhension et l'optimisation de ces problèmes. Enfin, la thèse examine l'optimisation géométrique robuste, qui prend en compte les distributions incertaines des variables aléatoires. Cette approche vise à développer des algorithmes d'optimisation résistant à l'incertitude dans les distributions de probabilité sous-jacentes, garantissant des solutions plus fiables et stables.

making them well-suited for optimization scenarios that involve sequential decision-making or iterative processes. Moreover, one of the key advantages of neurodynamic approaches is their hardware implementation feasibility. The primary objective of this thesis is to develop neurodynamic algorithms that are efficient and effective in solving chance-constrained geometric optimization problems. The thesis begins by focusing on chance-constrained geometric programs involving independent random variables. In addition, a specific type of geometric programs known as rectangular programs is also examined in detail. The objective is to understand the characteristics and complexities associated with this subclass of geometric programs. Subsequently, the thesis explores applying copula theory to address chance-constrained geometric programs with dependent random variables. Copula theory provides a mathematical framework for modeling and analyzing the dependence structure between random variables, thereby enhancing the understanding and optimization of these problems. Lastly, the thesis investigates distributionally robust geometric optimization, which considers uncertain distributions of random variables. This approach focuses on developing optimization algorithms that are robust against uncertainty in the underlying probability distributions, ensuring more reliable and stable solutions.

Synthèse

Dans son exploration des techniques visant à minimiser les coûts dans les problèmes de conception en ingénierie, Zener a initié les premières recherches dans ce qui est maintenant connu sous le nom de programmation géométrique. Le travail de Zener, ainsi que les articles ultérieurs de Duffin et Peterson ont posé les bases fondamentales de ce domaine. Le terme "programmation géométrique" a été adopté en raison du rôle significatif joué par l'inégalité arithmético-géométrique dans son développement initial. Au départ, la programmation géométrique se concentrait principalement sur la minimisation de fonctions posynomiales tout en respectant des contraintes d'inégalité. Ainsi, on aurait pu l'appeler programmation posynomiale plutôt que programmation géométrique. La programmation géométrique offre un moyen d'exprimer et d'analyser de nombreux problèmes d'optimisation importants dans un format séparable, même lorsqu'ils sont généralement considérés comme inséparables. L'aspect fondamental de cette approche réside dans l'exploitation des propriétés linéaires inhérentes au problème en cours.

Traditionnellement, la programmation géométrique a été appliquée à des problèmes avec des valeurs de paramètres connues et précises. Cependant, dans la vie réelle, les valeurs observées de ces paramètres sont souvent imprécises. Ces données incertaines peuvent prendre différentes formes, notamment des plages bornées, des intervalles, des ensembles flous ou des variables aléatoires. Afin de gérer ces incertitudes, des extensions et des modifications de l'approche traditionnelle de la programmation géométrique ont été développées. Ces variantes intègrent des techniques de l'analyse par intervalles, de la programmation stochastique ou de l'optimisation robuste. En tenant compte de la nature incertaine des paramètres, ces approches fournissent des solutions plus robustes et fiables aux problèmes de programmation géométrique de la vie réelle. Par conséquent, avec l'incorporation de données incertaines ou imprécises, la programmation géométrique devient un outil polyvalent et adaptable pour résoudre des problèmes d'optimisation dans des applications pratiques.

Dans la programmation géométrique stochastique, il est admis que certains coefficients et/ou exposants dans le problème de programmation géométrique ne sont pas connus précisément, et leur connaissance incomplète est représentée à l'aide de modèles probabilistes ou aléatoires. L'incertitude associée à ces coefficients peut provenir de diverses sources telles que des erreurs de mesure, de la variabilité des données ou du caractère intrinsèquement aléatoire du système modélisé. Au lieu de supposer des valeurs déterministes pour ces paramètres incertains, la programmation géométrique stochastique intègre leur caractère aléatoire dans la formulation du problème.

Tout au long de la vie, les réseaux cérébraux sont constamment en effervescence avec des ondes d'activité. D'innombrables signaux s'harmonisent et oscillent, formant l'architecture fonctionnelle sous-jacente qui façonne chaque aspect de l'existence. Les émotions, les niveaux de stress, les préférences musicales, les aspirations et les rêves trouvent leur origine dans le comportement dynamique de ces réseaux, complété par des systèmes de mémoire adaptative qui cherchent à mieux s'aligner sur l'environnement en constante évolution de chacun. L'ampleur considérable et l'intégration complexe des réseaux dans le cerveau humain posent des défis redoutables pour comprendre si des formes traditionnelles de calcul s'y produisent et comment. Les preuves indiquent que le traitement de l'information se déroule à divers niveaux, englobant des protéines modulatrices au sein des cellules individuelles, des microcircuits corticaux et des réseaux fonctionnels étendus couvrant l'ensemble du cerveau. Cependant, la compréhension expérimentale du fonctionnement du cerveau progresse lentement. Heureusement, des ingénieurs ingénieux ont élaboré des algorithmes qui simulent partiellement certains aspects de ces réseaux. Bien qu'aucun modèle unique ne puisse pleinement encapsuler la complexité pure et le comportement du cerveau humain, ces outils offrent aux chercheurs une précieuse fenêtre sur la manière dont l'information pourrait être calculée et finalement représentée dans l'activité de réseaux distribués.

Une catégorie importante de problèmes logiques découlant de scénarios du monde réel peut être formulée comme des problèmes d'optimisation, qui peuvent être compris qualitativement comme une quête de la solution optimale. Ces problèmes sont courants dans des domaines tels que l'ingénierie et le commerce, ainsi que dans des défis perceptuels qui nécessitent une résolution rapide par les systèmes nerveux. Le concept fondamental derrière l'approche des réseaux neuronaux pour l'optimisation est de créer une fonction d'énergie non négative et d'établir un système dynamique qui imite un réseau neuronal artificiel. Ce système dynamique prend généralement la forme d'équations différentielles ordinaires du premier ordre. L'objectif est que le système dynamique converge vers un état statique ou un point d'équilibre, qui correspond à la solution du problème d'optimisation sous-jacent, à partir d'un point initial. Un avantage notable de l'utilisation de réseaux neuronaux pour résoudre des problèmes d'optimisation est leur implémentabilité matérielle. En d'autres termes, ces réseaux neuronaux peuvent être mis en oeuvre à l'aide de circuits intégrés. De plus, les réseaux neuronaux conçus pour une mise en oeuvre par circuit offrent des capacités de traitement en temps réel, renforçant ainsi leur applicabilité pratique. sions de la programmation géométrique sont discutées, qui assouplissent la contrainte posynomiale et élargissent la gamme de fonctions pouvant être optimisées. Le chapitre se conclut par une introduction à l'optimisation sous contraintes de probabilité, également connue sous le nom de programmation probabiliste ou stochastique, qui traite des incertitudes et des risques dans la prise de décision. Il donne un aperçu des concepts de programmation sous contraintes de probabilité, des résultats établis et des approches computationnelles pour résoudre les problèmes d'optimisation sous contraintes de probabilité. Les applications du monde réel de la programmation sous contraintes de probabilité sont également examinées, mettant en évidence son utilité dans la gestion des incertitudes et la prise de décisions robustes.

Le chapitre 3 fournit un contexte théorique sur les réseaux neuronaux dynamiques et l'optimisation biconvexe. Le chapitre commence par mettre en avant la complexité et l'adaptabilité des réseaux cérébraux, ce qui inspire l'utilisation de réseaux neuronaux pour l'optimisation. Il met en évidence le potentiel des réseaux neuronaux pour imiter des systèmes dynamiques et converger vers des solutions optimales. Le chapitre explore l'histoire des modèles de réseaux neuronaux pour l'optimisation, en commençant par le modèle de Hopfield et Tank pour le Problème du Voyageur de Commerce. Il discute des extensions et des avancées réalisées dans les approches de réseaux neuronaux, notamment l'incorporation de paramètres de pénalité et le développement de modèles pour les problèmes d'optimisation contraints et non lisses. Plusieurs modèles de réseaux neuronaux spécifiques sont présentés, chacun conçu pour aborder différents types de problèmes d'optimisation. Ces modèles comprennent des réseaux neuronaux récurrents basés sur des pénalités pour les problèmes d'optimisation contraints, des approches neurodynamiques pour les problèmes d'optimisation convexe avec des contraintes générales, et des réseaux récurrents à une couche pour les problèmes d'optimisation convexe non lisses. Les avantages de ces modèles de réseaux neuronaux sont mis en avant, tels que leur efficacité computationnelle, leur simplicité, et leur capacité à traiter une large gamme de problèmes d'optimisation. Le chapitre présente également les applications des réseaux neuronaux récurrents dans des problèmes du monde réel, tels que l'optimisation de portefeuille dynamique, l'allocation de force de freinage électro-hydraulique, le pronostic de la santé des éoliennes, et le fonctionnement optimal des microgrids électriques. Pour offrir une compréhension plus approfondie des réseaux neuronaux récurrents, deux modèles spécifiques sont présentés dans les sous-sections suivantes. Ces modèles démontrent l'application des réseaux neuronaux récurrents dans la résolution de problèmes d'optimisation non linéaires convexes. Le chapitre passe ensuite à l'optimisation biconvexe. Nous mettons en évidence plusieurs avancées et applications de l'optimisation biconvexe au cours des dernières années. Le chapitre souligne la nature complexe de la résolution des problèmes de programmation biconvexe, nécessitant l'exploration de nouvelles théories et méthodes de solution. Dans cette section, deux résultats théoriques significatifs liés à l'optimum partiel de la programmation biconvexe sont présentés, en utilisant la fonction de pénalité de l'objectif. Le premier résultat établit l'équivalence entre la condition partielle de Karush-Kuhn-Tucker (KKT) et la propriété d'exactitude par-tielle de la fonction de pénalité de l'objectif dans la programmation biconvexe. Cette découverte fournit un aperçu précieux de la relation entre la condition KKT et le comportement de la fonction de pénalité de l'objectif, permettant une compréhension plus approfondie du processus d'optimisation. Le deuxième résultat démontre l'équivalence entre la condition de stabilité partielle et la propriété d'exactitude partielle de la fonction de pénalité de l'objectif dans la programmation biconvexe. Ce résultat met en lumière l'importance de la condition de stabilité par rapport à la fonction de pénalité de l'objectif, offrant des informations précieuses sur le comportement de convergence des algorithmes visant à résoudre l'optimum partiel de la programmation biconvexe. Ces résultats théoriques offrent des garanties cruciales pour la convergence des algorithmes conçus pour traiter l'optimum partiel des problèmes de programmation biconvexe. En établissant l'équivalence entre les conditions clés et les propriétés de la fonction de pénalité de l'objectif, cette thèse fournit une base solide pour le développement d'algorithmes efficaces et performants capables de converger de manière fiable vers des solutions optimales pour la programmation biconvexe.

Le chapitre 4 se concentre sur les problèmes d'optimisation géométrique conjoints sous contraintes de probabilité où les coefficients suivent une distribution normale et les vecteurs de lignes de la matrice sont indépendants. Le chapitre commence par introduire les programmes géométriques et leurs formulations équivalentes déterministes. Ensuite, le chapitre explore le concept d'approximations convexes, qui sont des techniques utilisées pour transformer des problèmes d'optimisation non convexes en problèmes convexes. En approximant le problème d'origine par une formulation convexe, il devient possible d'appliquer des algorithmes d'optimisation convexes existants pour trouver des solutions optimales. Cette approche améliore la traitabilité des problèmes d'optimisation géométrique conjoints sous contraintes de probabilité. Contrairement aux méthodes existantes qui reposent sur des approximations convexes, nous proposons une approche de réseau neuronal dynamique récurrent pour résoudre le programme géométrique stochastique. L'efficacité de cette approche est évaluée à travers des expériences numériques, notamment un problème d'optimisation de forme tridimensionnelle et un problème de transport multidimensionnel. Ces expériences démontrent la capacité de l'approche du réseau neuronal dynamique à traiter efficacement les problèmes d'optimisation géométrique conjoints sous contraintes de probabilité. De plus, le chapitre introduit une approche neurodynamique spécialement conçue pour résoudre des programmes rectangulaires, qui sont un cas particulier des programmes géométriques. Nous transformons le problème stochastique en un problème déterministe et utilisons une transformation logarithmique combinée avec l'inégalité arithmético-géométrique pour le convertir en un problème biconvexe. L'approche proposée utilise des réseaux neuronaux pour mettre à jour de manière itérative les variables du problème d'optimisation et trouver des solutions optimales. En appliquant l'approche neurodynamique, le chapitre démontre la capacité à résoudre efficacement des programmes rectangulaires. Enfin, le chapitre présente une étude de cas sur la maximisation du Rapport Signal sur Bruit d'Interférence (SINR) pour les systèmes Massive Multiple Input Multiple Output (MIMO).

Cette étude de cas met en évidence l'application pratique de l'optimisation géométrique conjointe sous contraintes de probabilité dans les systèmes de communication sans fil. L'objectif est d'optimiser l'allocation de ressources dans un système MIMO pour maximiser le SINR, améliorant ainsi les performances du système.

Le chapitre 5 se concentre sur les problèmes d'optimisation conjointe sous contraintes de probabilité dépendantes. Il introduit la théorie des copules et les vecteurs aléatoires elliptiquement symétriques comme fondements théoriques pour la modélisation et la résolution de ces problèmes. Le chapitre présente une approche de réseau neuronal dynamique pour résoudre les programmes linéaires et les programmes géométriques, en fournissant des formulations équivalentes déterministes pour les deux. Des expériences numériques sont menées pour valider l'efficacité et l'efficience de l'approche proposée dans la résolution de problèmes d'optimisation conjointe sous contraintes de probabilité dépendantes dans des applications pratiques. Le chapitre commence par introduire la théorie des copules et les vecteurs aléatoires elliptiquement symétriques comme base théorique pour aborder de tels problèmes d'optimisation. Ces concepts fournissent la fondation pour la modélisation et la résolution de problèmes impliquant des variables aléatoires dépendantes. Le chapitre explore ensuite les programmes linéaires et leurs formulations équivalentes déterministes dans le contexte des problèmes d'optimisation conjointe sous contraintes de probabilité dépendantes. En formulant le problème de manière déterministe, il devient plus adaptable aux techniques de résolution. Nous proposons ensuite une approche de réseau neuronal dynamique pour résoudre efficacement ces programmes linéaires. Cette approche offre une méthode novatrice pour traiter les problèmes d'optimisation conjointe sous contraintes de probabilité dépendantes sans recourir aux algorithmes d'optimisation conventionnels. Pour valider l'efficacité de la méthode proposée, le chapitre présente des expériences numériques. Ces expériences impliquent la résolution de différents problèmes d'optimisation conjointe sous contraintes de probabilité dépendantes en utilisant l'approche de réseau neuronal dynamique proposée. En plus des programmes linéaires, le chapitre couvre également les programmes géométriques. Le chapitre présente les formulations équivalentes déterministes pour les programmes géométriques et étend l'approche du réseau neuronal dynamique pour résoudre ces problèmes. Cette extension permet la résolution efficace de problèmes d'optimisation géométrique conjointe sous contraintes de probabilité dépendantes. Tout comme dans le chapitre précédent, des expériences numériques sont menées pour évaluer la performance de la méthode proposée pour les programmes géométriques. Ces expériences visent à évaluer l'exactitude et l'efficacité de la méthode dans la résolution de problèmes d'optimisation géométrique conjointe sous contraintes de probabilité dépendantes dans des scénarios pratiques.

Le chapitre 6 examine les problèmes d'optimisation géométrique robustes à la distribution. Il introduit le sujet et se concentre sur les programmes linéaires, en proposant une approche neurodynamique duplex pour résoudre ces problèmes. Le chapitre fournit une analyse de la convergence et mène des expériences numériques pour évaluer les performances de la méthode proposée. De plus, le chapitre explore les programmes géométriques robustes à la distribution avec des ensembles d'incertitude définis par les deux premiers moments d'ordre et présente une approche dynamique de réseau neuronal récurrent pour les résoudre. Des expériences numériques sont menées pour valider l'efficacité de cette approche. Le chapitre commence par une introduction au sujet, en donnant un aperçu des défis et des motivations derrière l'optimisation robuste à la distribution. Le chapitre se concentre ensuite sur les programmes linéaires et présente une approche neurodynamique duplex pour résoudre les problèmes d'optimisation géométrique robustes à la distribution dans ce contexte. L'approche proposée utilise des architectures de réseaux neuronaux et des systèmes dynamiques pour mettre à jour de manière itérative les variables de décision et atteindre la convergence. Le chapitre fournit une analyse détaillée de la convergence, démontrant la stabilité et les propriétés de convergence de l'approche neurodynamique duplex. Pour évaluer les performances de la méthode proposée, le chapitre mène des expériences numériques. Ces expériences résolvent divers programmes linéaires robustes à la distribution en utilisant l'approche neurodynamique duplex. Les résultats des expériences sont analysés pour évaluer l'exactitude, l'efficacité et l'efficacité de la méthode proposée dans des scénarios pratiques. En plus des programmes linéaires, le chapitre explore les programmes géométriques robustes à la distribution avec des ensembles d'incertitude définis par les deux premiers moments d'ordre. Le chapitre présente une approche dynamique de réseau neuronal récurrent spécialement conçue pour résoudre ces programmes géométriques robustes à la distribution. Cette approche tire parti de l'architecture de réseau neuronal récurrent pour capturer la dynamique temporelle du processus d'optimisation. Des expériences numériques sont menées pour valider l'efficacité de l'approche proposée pour les programmes géométriques robustes à la distribution. Ces expériences résolvent divers problèmes d'optimisation avec des ensembles d'incertitude définis par les deux premiers moments d'ordre en utilisant l'approche dynamique de réseau neuronal récurrent. Les résultats des expériences fournissent des informations sur les performances et l'efficacité de la méthode proposée pour résoudre les problèmes d'optimisation géométrique robustes à la distribution. Le chapitre 7 sert de résumé complet des principales découvertes et contributions présentées tout au long de la thèse. Il commence par revisiter les objectifs de recherche et met en évidence l'importance des approches neurodynamiques proposées pour résoudre les problèmes d'optimisation géométrique sous contraintes de probabilité. Le chapitre résume ensuite les principales contributions de la recherche, en mettant l'accent sur les méthodologies et techniques novatrices développées pour résoudre différents types de problèmes d'optimisation géométrique. Il discute des avantages des approches neurodynamiques proposées, tels que leur capacité à gérer la non-convexité, à traiter des contraintes probabilistes conjointes et à fournir des solutions robustes. De plus, le chapitre explore les implications et les applications potentielles des approches neurodynamiques proposées. Il discute de la manière dont ces approches peuvent être appliquées dans des scénarios réels dans différents domaines, tels que le transport, l'optimisation de forme, le traitement du signal, et bien plus encore. Les avantages potentiels de l'intégration de ces approches dans les processus de prise de décision sont mis en évidence, no-tamment l'amélioration de l'exactitude, de l'efficacité et de la robustesse des solutions d'optimisation. Dans la section de conclusion du chapitre, des orientations futures possibles pour la recherche dans le domaine de l'optimisation géométrique sont discutées. Ces orientations peuvent inclure une investigation plus approfondie de types spécifiques de programmes géométriques, l'exploration de l'intégration de contraintes ou objectifs supplémentaires, l'extension des approches neurodynamiques pour traiter des problèmes à plus grande échelle, ou l'adaptation des méthodologies à des domaines d'application spécifiques.

Acknowledgements "Sometimes, others help us by simply being in our lives." -Dr. Ahmed Khaled Tawfik

At the forefront of my appreciation stands Professor Abdel Lisser, whose mentorship has exceeded the conventional bounds of a supervisor. His support, guidance, and wholehearted assistance across all facets of my research have left an indelible mark on completing this work. His mentorship is beyond measure, and I consider myself very fortunate to have been supervised during those three years by Professor Lisser.

I would also like to extend my thanks to the members of the jury, for their time, expertise, contributions and thoughtful feedback. I would also like to acknowledge my colleagues, Dawen, Shangyuan and Nam, whose insightful dialogues and thoughtprovoking discussions have been a constant source of inspiration. The collaborative atmosphere we shared has pushed me to engage in critical thinking and has added depth to my research expedition.

To my parents, Mom and Dad, I owe an immeasurable debt of gratitude for their unconditional love, support, and the countless sacrifices they have made on my behalf shaping not only my academic pursuits but the very essence of my being.

I extend my warmest thanks to my wonderful mother-in-law, whose tender morning messages and steadfast support have been a source of immeasurable comfort and positivity throughout my days. Her thoughtful gestures and continuous encouragement have woven a tapestry of warmth and brightness into my daily life.

I must take a moment to express my profound gratitude to my beloved sisters, Chaimae and Oumayma, for the presence they have maintained and the unflagging support they have provided throughout this remarkable endeavor. Their enduring companionship has been a beacon of strength, a steady source of encouragement, and a reminder of the profound bonds that family shares during both the joys and challenges of life's journey.

Lastly, I want to convey my deepest gratitude to my husband, Badreddine. His understanding, unconditional support and boundless love have been my steadfast pillars of strength. I am profoundly thankful for being my confidant, my cheerleader, my refuge in times of both celebration and uncertainty and simply for his enduring presence in my life.

To all those mentioned and to the countless others who have contributed to my academic and personal growth, your roles have been immeasurable. This thesis is a testament to the support and love that have surrounded me, and I am forever grateful for each of you.

with love, Siham

-Introduction

This thesis aims to explore the application of neurodynamic approaches in the field of chance-constrained geometric optimization. The thesis focuses on developing efficient algorithms and solution methods for solving complex optimization problems with geometric constraints under uncertain conditions. The thesis consists of seven chapters, including a comprehensive literature review, theoretical background, and various case studies. This introduction provides an overview of the research plan and the key elements discussed in each chapter.

Chapter 2 provides a literature review on geometric optimization. The chapter discusses solution approaches and algorithms for solving geometric programs, including dual and primal methods. Various algorithms and techniques have been reviewed, such as condensation, interior-point methods, and primal-dual algorithms. The chapter also explores the application of geometric programming in different domains, including mechanical engineering, chemical engineering, power control, and finance. Furthermore, extensions to geometric programming are discussed, which relax the posynomial constraint and expand the range of functions that can be optimized. The chapter concludes with an introduction to chance-constrained optimization, also known as probabilistic or stochastic programming, which addresses uncertainties and risks in decision-making. It provides an overview of chance-constrained programming concepts, established results, and computational approaches for solving chance-constrained optimization problems. Real-world applications of chance-constrained programming are also surveyed, highlighting its usefulness in handling uncertainties and making robust decisions.

Chapter 3 provides a theoretical background on dynamical neural networks and biconvex optimization. The chapter begins by emphasizing the complexity and adaptability of the brain's networks, which inspire the use of neural networks for optimization. It highlights the potential of neural networks to mimic dynamic systems and converge to optimal solutions. The chapter explores the history of neural network models for optimization, starting with Hopfield and Tank's model for the Traveling Salesman Problem. It discusses the extensions and advancements made in neural network approaches, including the incorporation of penalty parameters and the development of models for constrained and nonsmooth optimization problems. Several specific neural network models are presented, each designed to address different types of optimization problems. These models include penalty-based recurrent neural networks for constrained optimization problems, neurodynamic approaches for convex optimization problems with general constraints, and one-layer recurrent networks for non-smooth convex optimization problems. The advantages of these neural network models are highlighted, such as their computational efficiency, simplicity, and ability to handle a broad range of optimization problems. The chapter also showcases the applications of recurrent neural networks in real-world problems, such as dynamic portfolio optimization, electro-hydraulic braking force allocation, wind turbine health prognosis, and optimal operation of electrical microgrids. To provide a deeper understanding of recurrent neural networks, two specific models are presented in the subsequent subsections. These models demonstrate the application of recurrent neural networks in solving convex nonlinear optimization problems. The chapter then transitions to biconvex optimization. We highlight several advancements and applications of biconvex optimization in recent years. The chapter emphasizes the challenging nature of solving biconvex programming problems, necessitating the exploration of new theories and solution methods. In this section, two significant theoretical results related to the partial optimum of biconvex programming are presented, utilizing the objective penalty function. The first result establishes the equivalence between the partial Karush-Kuhn-Tucker (KKT) condition and the partial exactness property of the objective penalty function in biconvex programming. This finding provides a valuable insight into the relationship between the KKT condition and the behavior of the objective penalty function, enabling a deeper understanding of the optimization process. The second result demonstrates the equivalence between the partial stability condition and the partial exactness property of the objective penalty function in biconvex programming. This result sheds light on the significance of the stability condition in relation to the objective penalty function, offering valuable insights into the convergence behavior of algorithms aiming to solve the partial optimum of biconvex programming. These theoretical results offer crucial guarantees for the convergence of algorithms designed to tackle the partial optimum of biconvex programming problems. By establishing the equivalence between key conditions and properties of the objective penalty function, this thesis provides a solid foundation for developing efficient and effective algorithms that can reliably converge towards optimal solutions for biconvex programming.

Chapter 4 focuses on joint chance-constrained geometric optimization problems where the coefficients are normally distributed, and the matrix row vectors are independent. The chapter starts by introducing geometric programs and their deterministic equivalent formulations. The chapter then explores the concept of convex approximations, which are techniques used to transform non-convex optimization problems into convex ones. By approximating the original problem with a convex formulation, it becomes possible to apply existing convex optimization algorithms to find optimal solutions. This approach enhances the tractability of joint chance-constrained geometric optimization problems. Unlike existing methods that rely on convex approximations, we propose a recurrent dynamical neural network approach to solve the stochastic geometric program. The effectiveness of this approach is evaluated through numerical experiments, including a threedimensional shape optimization problem and a multidimensional transportation problem. These experiments demonstrate the ability of the dynamical neural network approach to handle joint chance-constrained geometric optimization problems efficiently. Furthermore, the chapter introduces a neurodynamic approach specifically designed to solve rectangular programs, which are a special case of geometric programs. We transform the stochastic problem into a deterministic one and utilize a logarithmic transformation combined with the arithmetic-geometric mean inequality to convert it into a biconvex problem. The proposed approach utilizes neural networks to iteratively update the vari-ables of the optimization problem and find optimal solutions. By applying the neurodynamic approach, the chapter demonstrates the capability to solve rectangular programs effectively. Finally, the chapter presents a case study on maximizing the Signal to Interference Noise Ratio (SINR) for Massive Multiple Input Multiple Output (MIMO) systems. This case study highlights the practical application of joint chance-constrained geometric optimization in wireless communication systems. The goal is to optimize the allocation of resources in a MIMO system to maximize the SINR, thereby improving the system's performance.

Chapter 5 focuses on dependent joint chance-constrained optimization problems. It introduces copula theory and elliptically symmetric random vectors as the theoretical foundations for modeling and solving these problems. The chapter presents a dynamical neural network approach to solve linear programs and geometric programs, providing deterministic equivalent formulations for both. Numerical experiments are conducted to validate the effectiveness and efficiency of the proposed approach in solving dependent joint chance-constrained optimization problems in practical applications. The chapter begins by introducing copula theory and elliptically symmetric random vectors as the theoretical underpinnings for addressing such optimization problems. These concepts provide the foundation for modeling and solving problems involving dependent random variables. The chapter proceeds to explore linear programs and their deterministic equivalent formulations in the context of dependent joint chance-constrained optimization problems. By formulating the problem deterministically, it becomes more amenable to solution techniques. We then propose a dynamical neural network approach to solve these linear programs efficiently. This approach offers a novel method to handle dependent joint chance-constrained optimization problems without resorting to conventional optimization algorithms. To validate the effectiveness of the proposed method, the chapter presents numerical experiments. These experiments involve solving various dependent joint chance-constrained optimization problems using the proposed dynamical neural network approach. In addition to linear programs, the chapter also covers geometric programs. The chapter presents the deterministic equivalent formulations for geometric programs and extends the dynamical neural network approach to address these problems. This extension allows for the efficient solution of dependent joint chanceconstrained geometric optimization problems. Similar to the previous Chapter, numerical experiments are conducted to evaluate the performance of the proposed approach for geometric programs. These experiments aim to assess the accuracy and efficiency of the method in solving dependent joint chance-constrained geometric optimization problems in practical scenarios.

Chapter 6 investigates distributionally robust geometric optimization problems. It introduces the topic and focuses on linear programs, proposing a neurodynamic duplex approach to solve these problems. The chapter provides a convergence analysis and conducts numerical experiments to evaluate the performance of the proposed method. Furthermore, the chapter explores distributionally robust geometric programs with uncertainty sets defined by the first two order moments and presents a dynamical recurrent neural network approach to solve them. Numerical experiments are conducted to validate the effectiveness of this approach. The chapter starts with an introduction to the topic, providing an overview of the challenges and motivations behind distributionally robust optimization. The chapter then focuses on linear programs and presents a neurodynamic duplex approach to solve distributionally robust geometric optimization problems in this context. The proposed approach utilizes neural network architectures and dynamic systems to update the decision variables and achieve convergence iteratively. The chapter provides a detailed convergence analysis, demonstrating the stability and convergence properties of the neurodynamic duplex approach. To assess the performance of the proposed method, the chapter conducts numerical experiments. These experiments solve various distributionally robust linear programs using the neurodynamic duplex approach. The results of the experiments are analyzed to evaluate the accuracy, efficiency, and effectiveness of the proposed method in practical scenarios. In addition to linear programs, the chapter explores distributionally robust geometric programs with uncertainty sets defined by the first two order moments. The chapter presents a dynamical recurrent neural network approach specifically designed to solve these distributionally robust geometric programs. This approach leverages the recurrent neural network architecture to capture the temporal dynamics of the optimization process. Numerical experiments are conducted to validate the effectiveness of the proposed approach for distributionally robust geometric programs. These experiments solve various optimization problems with uncertainty sets defined by the first two order moments using the dynamical recurrent neural network approach. The results of the experiments provide insights into the performance and efficacy of the proposed method in addressing distributionally robust geometric optimization problems.

Chapter 7 serves as a comprehensive summary of the main findings and contributions presented throughout the thesis. It begins by revisiting the research objectives and highlighting the significance of the proposed neurodynamic approaches for solving chanceconstrained geometric optimization problems. The chapter then summarizes the key contributions of the research, emphasizing the novel methodologies and techniques developed to address various types of geometric optimization problems. It discusses the advantages of the proposed neurodynamic approaches, such as their ability to handle nonconvexity, address joint probabilistic constraints, and provide robust solutions. Furthermore, the chapter explores the implications and potential applications of the proposed neurodynamic approaches. It discusses how these approaches can be applied in realworld scenarios across different domains, such as transportation, shape optimization, signal processing, and more. The potential benefits of incorporating these approaches into decision-making processes are highlighted, including improved accuracy, efficiency, and robustness of optimization solutions. In the concluding section of the chapter, possible future directions for research in the field of geometric optimization are discussed. These directions may include further investigation into specific types of geometric programs, exploring the integration of additional constraints or objectives, extending the neurodynamic approaches to handle larger-scale problems, or adapting the methodologies for specific application domains.

-Literrature review

. Geometric Optimization

In his investigation of techniques to minimize costs in engineering design problems, Zener [START_REF] Zener | A mathematical aid in optimizing engineering designs[END_REF][START_REF] Zener | A further mathematical aid in optimizing engineering designs[END_REF] sparked the initial research in what is now known as geometric programming (GP). Zener's work, along with subsequent papers by Duffin and Peterson [START_REF] Duffin | Dual programs and minimum cost[END_REF][START_REF] Duffin | Cost minimization problems treated by geometric means[END_REF][START_REF] Duffin | Duality theory for geometric programming[END_REF], as well as the book by Duffin, Peterson, and Zener [START_REF] Perlumutter | Geometric programming[END_REF], laid the foundational groundwork for this field. The term "geometric programming" was adopted due to the significant role played by the arithmetic-geometric mean inequality in its early development. Initially, the focus of geometric programming was primarily on minimizing posynomial functions while adhering to inequality constraints. Hence, it could have been named posynomial programming instead of geometric programming. Geometric programming offers a means to express and analyze numerous significant optimization problems in a separable format, even when they are typically considered inseparable. The fundamental aspect of this approach lies in leveraging the linear properties inherent in the problem at hand.

. Geometric Optimization: Definition and Concepts

Let x 1 , x 2 , ..., x n denote n real positive variable, then a real-valued function f in the form

f (x) = cx α 1 1 x α 2 2 ...x αn n
where c > 0 and α i , i = 1, ..., n are real constants is called monomial function. The constant c within the monomial is commonly known as the coefficient, while the constants α 1 , ..., α n are referred to as the exponents, signifying the respective powers of the variables in the monomial. For instance, the monomial 4.3x 5 1 x -0.0005 2 of variables x 1 and x 2 where the coefficient is 4.3 and the exponents are 5 and -0.0005.

A function f : R n -→ R n is called a posynomial, if it takes the form of a non-negative sum of monomials. The term "posynomial" is coined by combining the words "positive" and "polynomial," implying a mathematical expression that possesses both the properties of positivity and being a polynomial.

f (x) = n i=1 c i M j=1 x α ij j 3.2x 1 x 0.4 2 + 1.5x 2 1 x -2 2 and √ 3x 1 + x 1 x -2 2 x -2
3 are posynmials. Geometric programming is used to minimize functions that are formulated as posynomials while satisfying constraints of the same posynomial form. A general form of a geometric program can be then stated as follows

min x∈R M ++ i∈I 0 c i M j=1 x α ij j , subject to i∈I k c i M j=1 x α ij j ≤ 1, k = 1, ...., K, (2.1)
where x is a strictly positive M -dimensional vector, the exponents α ij are arbitrary real numbers and the coefficients c i are positive.

Monomials and posynomials are typically non-convex functions. However, it is possible to transform them into convex functions by applying a simple variable transformation. This transformation allows us to convert the original non-convex functions into equivalent convex forms, which can then be effectively optimized using convex optimization techniques. Program (2.1) has then an equivalent convex program given by taking t i = log(x i ). min

t i∈I 0 c i exp   M j=1 α ij t j   , subject to i∈I k c i exp   M j=1 α ij t j   ≤ 1, k = 1, ...., K, (2.2)
The equivalence between geometric programs and convex programs has significant implications. One of the most crucial consequences is that the principles and techniques of convex programming can be directly applied to geometric programs, particularly when expressed in their equivalent logarithmic form. This enables the utilization of general theories and algorithms from convex programming to solve geometric programs effectively. Nevertheless, in certain real-world scenarios, it may be impractical to express the problem in its conventional geometric format, rendering it non-convex. Consequently, solving such problems becomes notably challenging, even with approximate methods. In such instances, it becomes highly advantageous to employ straightforward approaches capable of generating a favorable suboptimal solution, if not the optimal solution. Geometric programming is named as such because it is founded on the principles of the arithmetic-geometric mean inequality. This inequality plays a fundamental role in formulating and solving geometric programming problems. By leveraging the arithmeticgeometric inequality, geometric programming offers a powerful framework for optimizing and solving problems involving posynomials and monomials, ultimately leading to efficient solutions in various domains. In [START_REF] Perlumutter | Geometric programming[END_REF], Duffin, Peterson, and Zener employ an extension of the weighted arithmetic-geometric mean inequality to derive lower bounds on the minimum value of a geometric program (2.1). The inequality states that for every x > 0 and y > 0 vectors in R n , we have

n i=1 x λ i ≥ n i=1 x i y i y i λ λ , (2.3)
where n i=1 y i = λ.

We can apply the the arithmetic-geometric mean inequality to the objective function of (2.1), we obtain

i∈I 0 c i M j=1 x α ij j ≥ i∈I 0 c i M j=1 x α ij j δ i δ i , (2.4) with i∈I 0 δ i = 1.
The constraints can be expressed as

1 ≥   i∈I k c i M j=1 x α ij j   λ k ≥ i∈I k c i M j=1 x α ij j δ i δ i λ λ k k , k = 1, ..., K (2.5)
Those inequalities are crucial keys for the duality in geometric programming, more details can be found in [START_REF] Ecker | Geometric programming: Methods, computations and applications[END_REF]. The primal problem (2.1) is often characterized by its complexity, while the dual version of the problem is comparatively simpler to solve. In practice, it is common to focus on solving the dual problem, as it provides valuable insights. The dual problem of (2.1) consists on maximizing the following dual function

K k=0 i∈I k c i δ i δ i K k=1 λ λ k k (2.6)
subject to the following normality and orthogonality conditions

i∈I 0 δ i = 1 (2.7) K k=0 i∈I k a ij λ i = 0, j = 1, ..., M. (2.8)
where i∈I k δ i = λ k .

. Solution Approaches and Algorithms

The solution of geometric programs using numerical methods is a highly interesting topic, and researchers have developed various specialized algorithms for this purpose. A fundamental question arises regarding whether to solve the primal or the dual program. In the case of regular geometric programs, the dual program has a concave objective function and linear constraints, while the primal program typically contains nonlinear constraints. Consequently, it may seem reasonable to solve regular geometric programs through the dual approach. However, in practice, this approach is not always recommended because if certain dual variables need to be zero in the optimal solution, the objective function's gradient becomes unbounded, leading to significant numerical challenges. Moreover, since programs more general than regular geometric programs either lack a dual representation or do not provide apparent computational advantages through the dual, several algorithms based on the primal approach have been developed. Many of these algorithms employ condensation, a technique that combines a sum of positive terms into a single term, enabling the approximation of a nonconvex program with a convex one. There exist several computational methods to solve geometric programs, those methods can be divided into two families : dual methods that solve the dual program, and subsequently, the duality relations are utilized to derive a solution for the primal geometric program and primal methods.

The first work that proposed an algorithm for posynomial programming was the dual method of Frank [START_REF] Frank | An algorithm for geometric programming[END_REF]. Rijckaert and Martens [START_REF] Rijckaert | A condensation method for generalized geometric programming[END_REF] propose a dual condensation algorithm to come up with a solution to generalized geometric programs. Bricker and Rajgopal [START_REF] Bricker | Yet another geometric programming dual algorithm[END_REF] describe an algorithm for the geometric programming dual problem using an adaptation of the generalized LP algorithm. Kortanek et al. [START_REF] Kortanek | An infeasible interior-point algorithm for solving primal and dual geometric programs[END_REF] present an infeasible interior-point algorithm for solving primal and dual geometric programs. By employing the Zadeh's extension principle and transforming the primal problem of fuzzy geometric programming into its dual form, Saraj [START_REF] Saraj | Solving a posynomial geometric programming problem with fully fuzzy approach[END_REF] facilitates the conversion of the dual form into a pair of mathematical programs. To achieve this, he applies the α-cut on the objective function and the r-cut on the constraints in the dual form of geometric programming. This process allows us to obtain acceptable (α, r) optimal values, which represent solutions that strike a balance between precision (α) and satisfaction of constraints (r) in the fuzzy geometric programming problem. Rosenberg [START_REF] Rosenberg | On solving a primal geometric program by partial dual optimization[END_REF] introduces a novel technique that aims to generate an initial estimate of a primal solution for a posynomial geometric program using a dual feasible point that satisfies a mild condition. It is important to note that while the dual feasible point may not be optimal, this technique provides a means to derive a preliminary primal solution. By leveraging this approach, it becomes possible to obtain a starting point for further optimization iterations, which can subsequently lead to the discovery of improved primal solutions for the posynomial geometric program.

Even if during the initial phase greater emphasis was placed on the development of dual methods, Numerous techniques have been developed to solve the primal posynomial program directly. Dinkel et al. [START_REF] Dinkel | Computational aspects of cuttingplane algorithms for geometric programming problems[END_REF] apply four cutting plane algorithms to solve posynomial geometric programs. Chun-Feng et al [START_REF] Chun-Feng | Global optimization for sum of geometric fractional functions[END_REF] present an efficient branch and bound algorithm for globally solving the sum of geometric fractional functions under geometric constraints. Melissaratos and Souvaine [START_REF] Melissaratos | On solving geometric optimization problems using shortest paths[END_REF] employ shortest path theory to deal with a class of geometric optimization problems. Sassen et al. [START_REF] Sassen | Geometric optimization using nonlinear rotation-invariant coordinates[END_REF] use demonstrate the benefits of using nonlinear rotation-invariant coordinates on a collection of geometric optimization problems. Ojha and Ota [START_REF] Ojha | Multi-objective geometric programming problem with karushkuhntucker condition using-constraint method[END_REF] propose a method that addresses the convexity issues associated with the weighted sum technique and known as the ϵ-constraint method. This approach involves minimizing a primary objective while expressing the remaining objectives as inequality constraints. By treating the additional objectives as constraints, the -constraint method enables the exploration of the trade-offs between the primary objective and the secondary objectives in a more flexible manner. This approach allows for a more nuanced optimization process that considers both the primary objective and the satisfaction of the specified constraints.

The degree of difficulty of GP is defined as the difference between the total number of monomial terms in the objective and constraints, and one plus the number of variables. When the degree of difficulty is zero, solving a GP is equivalent to solving a system of linear equations. Remarkably, modern numerical methods demonstrate good performance regardless of the degree of difficulty. These methods can effectively handle GP problems with varying levels of complexity, making them suitable for a wide range of applications. In this survey, we take an in-depth look at two algorithms that make use of modern convex optimization techniques, namely the interior-point method [START_REF] Bürgisser | Interior-point methods for unconstrained geometric programming and scaling problems[END_REF] and an infeasible Primal-Dual algorithm [START_REF] Kortanek | An infeasible interior-point algorithm for solving primal and dual geometric programs[END_REF].

• Barrier-based interior-point method. The standard barrier-based interior-point method, commonly used for convex optimization, can be directly applied to GP problems. This approach offers a worst-case polynomial-time complexity and exhibits highly efficient performance that gracefully scales with the problem size in practical scenarios. By leveraging the interior-point method, solving GP problems becomes more tractable and efficient, enabling the optimization of larger-scale problems within reasonable computational timeframes. Let us consider the following GP problem formulated in convex form

min t f 0 (x) = i∈I 0 c i exp   M j=1 α ij x j   , subject to f k (x) = i∈I k c i exp   M j=1 α ij x j   ≤ 1, k = 1, ..., K, . (2.9)
The fundamental concept behind the barrier method is to solve a sequence of unconstrained problems by incorporating the constraints into a modified objective function. This new objective function is a weighted combination of the original objective function and a barrier function ϕ that represents the constraints. As the weight parameter t assigned to the original objective function increases, the unconstrained problem progressively approximates the original problem more closely. This iterative process allows for the optimization of the problem while gradually considering the constraints and finding solutions that satisfy them. We can describe the barrier method algorithm as follows Given a strictly feasible point x, and given initial values t 0 > 0, µ > 1, and an error tolerance ϵ > 0, the following steps are repeated:

-Compute x * (t) by minimizing the function tf 0 (x) + ϕ(x), starting at x. This is an unconstrained, smooth, convex minimization problem that can be solved using iterative methods such as the gradient descent method or Newton's method.

-Set x to the newly obtained solution x * (t).

-If the barrier parameter K t is less than or equal to the specified error tolerance ϵ, exit the loop.

-Update t, i.e., t = µt.

• The infeasible Primal-Dual algorithm. The infeasible primal-dual method is an optimization approach that generates subfeasible solutions, meaning solutions that do not satisfy all the constraints exactly. However, the method ensures that the primal and dual objective function values of these subfeasible solutions converge to the respective primal and dual optimal values.

The method involves solving a sequence of optimization problems, starting from an initial infeasible solution. The objective is to improve both the primal and dual objective function values iteratively while gradually satisfying the constraints.

By iteratively adjusting the solution and updating the Lagrange multipliers, the method aims to converge towards a solution that is feasible and has optimal objective function values for both the primal and dual problems.

While the generated solutions may not be strictly feasible, the infeasible primal-dual method provides a way to obtain increasingly better approximations to the primal and dual optimal values as the iterations progress. The details of the algorithm are given in [START_REF] Kortanek | An infeasible interior-point algorithm for solving primal and dual geometric programs[END_REF].

. Stochastic geometric programming

Traditionally, geometric programming has been applied to problems with known and precise parameter values. However, in real-life applications, the observed values of these parameters are often imprecise or vague. Such uncertain data can take different forms, including bounded ranges, intervals, fuzzy sets, or random variables.

In order to handle these uncertainties, extensions and modifications to the traditional geometric programming approach have been developed. These variations incorporate techniques from interval analysis, fuzzy optimization, stochastic programming, or robust optimization. By considering the uncertain nature of the parameters, these approaches provide more robust and reliable solutions to real-life geometric programming problems.

Therefore, with the incorporation of uncertain or imprecise data, geometric programming becomes a versatile and adaptable tool for solving optimization problems in practical applications. It enables decision-makers to make informed decisions and find optimal solutions even in the presence of uncertainty.

In stochastic geometric programming, it is acknowledged that certain coefficientsc i and/or exponents a ij in the geometric programming problem are not known precisely, and their incomplete knowledge is represented using probabilistic or random models. The uncertainty associated with these coefficients can arise from various sources such as measurement errors, variability in data, or inherent randomness in the system being modeled. Instead of assuming deterministic values for these uncertain parameters, stochastic geometric programming incorporates their randomness into the problem formulation.

The origins of stochastic geometric programming can be traced back to the work of Avriela and Wilde [START_REF] Avriel | Stochastic geometric programming[END_REF], where the exponents a ij in the geometric programming problem are deterministic, while the coefficients c i are considered as positive random variables.

In the case of individual probabilistic constraints, the constraints of (2.1) can be replaced by

P   i∈I k c i M j=1 x α ij j ≤ 1   ≥ 1 -ϵ k , k = 1, ...., K (2.10)
where 1 -ϵ k are the tolerance levels.

Jagannathan [START_REF] Jagannathan | A stochastic geometric programming problem with multiplicative recourse[END_REF] considers the coefficients and the exponents as random variables with known joint distribution functions. This allowed to model the uncertainty associated with these parameters and analyze the properties of the problem under this stochastic setting. The approach consists of a multiplicative recourse stochastic model. In this model, the recourse variables are introduced to rectify or adjust the solution in a proportional manner to account for the possible violations of the constraints due to the randomness in the parameters. Dupačová [START_REF] Dupačová | Stochastic geometric programming with an application[END_REF] presents stochastic geometric programming in the context of metal cutting optimization. Liu [START_REF] Liu | Posynomial geometric programming with parametric uncertainty[END_REF] introduces a procedure that allows for the derivation of lower and upper bounds for the objective function in a posynomial geometric programming problem when there is uncertainty in the cost and constraint parameters. To address this uncertainty, a transformation of the imprecise geometric program into a collection of conventional geometric programs was proposed. Each program within this collection corresponds to a specific combination of parameter values from the provided ranges. By solving these conventional geometric programs, the objective value can be computed. Liu [START_REF] Liu | Geometric programming with fuzzy parameters in engineering optimization[END_REF] presents a procedure to determine the fuzzy objective value in a fuzzy posynomial geometric programming problem. The fuzziness arises from the exponents of decision variables in the objective function, the cost and constraint coefficients, and the right-hand sides, all of which are represented as fuzzy numbers. The approach is based on Zadeh's extension principle, which allows for the transformation of the fuzzy geometric programming problem into a pair of two-level mathematical programs. By utilizing a duality algorithm and a simple algorithm, this pair of two-level mathematical programs is further transformed into a pair of conventional geometric programs. Liu et al. [START_REF] Liu | Stochastic geometric optimization with joint probabilistic constraints[END_REF] focus is on geometric programs with joint probabilistic constraints. Specifically, when the stochastic parameters follow independent normal distributions, an approximation approach is proposed. This approach involves approximating the original problem using piecewise linear functions and transforming it into a convex geometric program. Additionally, a sequential convex optimization algorithm is developed to find an upper bound. Khanjani et al. [START_REF] Khanjani | Fuzzy chance-constrained geometric programming: The possibility, necessity and credibility approaches[END_REF] formulate three variants of chance-constrained GP based on different perspectives: possibility, necessity, and credibility. Each variant represents a different way of handling uncertainty and capturing the level of confidence in achieving the desired outcomes.

. Applications of geometric programming

Since its appearance, geometric optimization has found many successful applications in various fields such as mechanical engineering [START_REF] Janushevskis | Shape Optimization of Mechanical Components for Measurement Systems[END_REF][START_REF] Fedorov | Geometric optimization of a self-adaptive robotic leg[END_REF][START_REF] Jha | A Discrete Data Base Multiple Objective Optimization of Milling Operation Through Geometric Programming[END_REF][START_REF] Lin | Finding multiple optimal solutions of signomial discrete programming problems with free variables[END_REF][START_REF] Tsai | Global optimization for signomial discrete programming problems in engineering design[END_REF], chemical engineering [START_REF] Clasen | The solution of the chemical equilibrium programming problem with generalized benders decomposition[END_REF][START_REF] Wall | Or practice-solving complex chemical equilibria using a geometric-programming based technique[END_REF], heat exchanger system optimization [START_REF] Avriel | An extension of geometric programming with applications in engineering optimization[END_REF], power control [START_REF] Li | Power control in cooperative cognitive radio networks by geometric programming[END_REF][START_REF] Stupar | Posynomial models of inductors for optimization of power electronic systems by geometric programming[END_REF], communication network systems [START_REF] Chiang | Geometric programming duals of channel capacity and rate distortion[END_REF][START_REF] Guan | Optimal Resource Allocation for video streaming over cognitive radio network via Geometric Programming[END_REF] and Finance [START_REF] Kabanov | A geometric approach to portfolio optimization in models with transaction costs[END_REF].

To showcase the versatility of geometric programming and highlight the formulation of specific models, we will delve into three detailed applications. These applications serve as examples that demonstrate the wide range of domains where geometric programming finds utility.

Design of a Two-bar Truss

To achieve minimum weight, we need to design the two-bar truss, as depicted in Figure 2.1, for a vertical load of 2P . The truss members have a tubular section with a mean diameter d and wall thickness t. The maximum allowable stress in each member is denoted as σ s . Using geometric programming, we aim to determine the values of h (height) and d (diameter) considering the following data • P is the vertical load (33000 lb)

• t is the wall thickness (0.1 in)

• b is the span (30 in)

• σ s is the maximum permissible stress (60000 psi)

• ρ is the density (0.3 lb/in³)

The objective function for this problem is defined as (2.12)

f (d, h) = 0.188d * 900 + h 2
The stress constraint ensures that the maximum stress is within the permissible limit of 60000. The goal of the geometric programming problem is to find the optimal values of d and h that minimize the objective function while satisfying the stress constraint. The square root functions in equations (2.11)

δ 01 -δ 11 = 0.
The solutions of (2.13) are given by δ * 01 = 1, δ * 11 = 1,δ * 21 = 0.5 and δ * 22 = 0.5. The maximum value of f is then given by f * = 0.188

Helical springs design

To determine the minimum volume design of the cone clutch depicted in Figure 2.2 while ensuring the transmission of a specified minimum torque, we can employ the outer and inner radii of the cone, denoted as R 1 and R 2 , as design variables. The objective function Figure 2.2: Schematic of a cone clutch [START_REF] Kulkarni | Modified probability collectives approach for solving engineering problems[END_REF] can then be formulated as follows

f (R 1 , R 2 ) = R 3 1 -R 3 2 (2.14)
subject to the following constraints

R 2 1 + R 1 R 2 + R 2 2 R 1 + R 2 ≥ 5, (2.15)
R 1 R 2 ≥ 2.
(2.16)

The problem was addressed using an iterative approach known as complementary geometric programming [START_REF] Avriel | Complementary geometric programming[END_REF]. Consequently, the final solution is determined asR * 1 = 4.2874, R * 2 = 2.1437 and f * = 68.916.

Production Model Under Risk

In our initial application, we focus on a straightforward profit-maximization model for a single-product firm with predetermined input prices and output price. We then extend the analysis to incorporate stochastic elements, introducing uncertainty into the input prices and output prices, thereby considering the stochastic case. This extension enables us to explore the implications of uncertainty on the firm's profit-maximization strategy within a geometric programming framework. The analysis begins with a simplified scenario where the production of a single output involves the combination of three inputs: labor x 3 , raw materials x 2 , and capital x 3 . The geometric programming formulation for the profit-maximization model in the deterministic case can be represented as follows min 0.14x 1 + 0.04x 2 + 0.06x 3

(2.17)

s.t. 10x -1.6 1 x -0.5 2 x -0.7 3 ≤ 1, (2.18)
0.1x 1 x -1 3 ≤ 1, (2.19)
5x -1 5 x 3 ≤ 1.
(2.20)

The dual program can be expressed as follows through max 0.14

δ 1 δ 1 0.04 δ 2 δ 2 0.06 δ 3 δ 3 10 δ 4 δ 4 0.1 δ 5 δ 5 5 δ 6 δ 6 s.t. δ 1 + δ 2 + δ 3 = 1, δ 1 -1.6δ 4 + δ 5 -δ 6 = 0, δ 2 -0.5δ 4 = 0, δ 3 -0.7δ 4 -δ 5 + δ 6 = 0, δ 1 , δ 2 , δ 3 ≥ 0.
The optimal solution is 0.6516 and we get As a straightforward approach to incorporate stochastic elements, we can consider a scenario where production needs to align with sales. In this context, sales are treated as a random variable with a normal distribution, with a mean of 1000 units and a variance of

0.1x 1 x -1 3 ≤ 1, (2.24)
5x -1 5 x 3 ≤ 1.

(2.25)

where y * denotes sales, a random variable with a normal distribution with a mean of ȳ = 1000 and variance of σ 2 = 100. The constraint (2.22) indicates that the firm aims to implement a production policy that ensures a high probability of meeting sales targets.

The tolerance measure α represents the desired probability threshold.

. Extensions of Geometric Programming

A notable limitation in the effective application of geometric programming to optimize engineering design has been the requirement that all functions involved in the problem must be posynomials. As a result, extensive research efforts have been undertaken to develop extensions and enhancements for geometric programming. These extensions aim to relax the posynomial constraint and allow for a broader range of functions to be incorporated into the optimization framework. These advancements have expanded the applicability and versatility of geometric programming in engineering design optimization.

Sigmonial geometric programming

SGP, signomial geometric programming, represents a significant expansion of GP. Several approaches have been proposed and numerous specialized algorithms have been developed to tackle SGP. However, despite these valuable contributions, the challenge of solving SGP problems remains an open issue. This difficulty primarily arises from the inherent non-convexity of SGP. Unlike GP problems, SGP problems retain their non-convex nature in both the primal and dual forms, and there is no known transformation capable of convexifying them. As a result, efficiently computing only a locally optimal solution for an SGP problem is currently achievable.

SGP is similar to Posynomial GP, with the key distinction that the term coefficients c i are no longer required to be positive. In Posynomial GP, all coefficients must be positive, whereas in SGP, the coefficients can have both positive and negative values. A general form of a signomial program can be stated as follows

min x∈R M ++ I 0 i=1 σ 0 i c 0 i M j=1 x a 0 ij j , (2.26)
s.t

I k i=1 σ k i c k i M j=1 x a k ij j ≤ σ k0 , k = 1, ..., K, c k i > 0, σ k i , σ k0 = +1 or -1.

29

The dual program is given by

max σ 00   K k=0 i∈I k c i δ i σ k i δ i K k=1 λ λ k σ k0 k   , (2.27)
s.t σ k0 i∈I k σ k i δ i = λ k , k = 1, ..., K, K k=0 i∈I k a ij δ i σ k i = 0, j = 1, ..., M, σ 00 i∈I 0 σ 0 i δ i = 1, δ i , λ k ≥ 0,
where σ 00 represents the sign of the primal objective at the optimum.

There exist a lot of approaches to solve SGP, we can state the following three major techniques

• Branch and bound: The first approach involves employing standard branch and bound techniques for general non-convex optimization. Gochet and Smeers [START_REF] Gochet | A branch-and-bound method for reversed geometric programming[END_REF] propose a branch-and-bound method for solving SGP, ensuring convergence to the global optimum. The subproblems involved in this method are convex and were combined with a cutting plane technique to generate linear subproblems. Shen and Jiao [START_REF] Shen | A new rectangle branch-and-pruning approach for generalized geometric programming[END_REF] propose an algorithm that successfully integrates the branch-andbound method with the pruning technique, leading to convergence to the global minimum. This achievement is made possible by iteratively refining the linear relaxation of the feasible region of the objective function and solving a series of linear relaxation problems. Through this process of successive refinement and solving relaxed problems, the algorithm progressively approaches the global minimum.

• Relaxations: Shen and Zhang [START_REF] Shen | Global optimization of signomial geometric programming using linear relaxation[END_REF] use an exponential variable transformation, tangential hypersurfaces and convex envelope approximations to obtain a linear relaxation of SGP. Consequently, the initial nonconvex nonlinear problem of SGP is transformed into a sequence of linear programming problems via the successive improvement of a linear relaxation of the objective function's feasible region. By iteratively solving this series of linear programming problems, the proposed algorithm demonstrates convergence to the global minimum of SGP. Chandrasekaran and Shah [START_REF] Chandrasekaran | Relative entropy relaxations for signomial optimization[END_REF] use a hierarchy of convex relaxations to obtain progressively tighter lower bounds on the optimal value of SGPs. These lower bounds are obtained by solving increasingly larger relative entropy optimization problems, which are convex programs formulated using linear and relative entropy functions. The key insight behind this approach is the observation that the relative entropy function, due to its joint convexity with respect to both arguments, offers a convex parameterization for specific sets of globally nonnegative signomials. This parameterization allows for efficiently computable nonnegativity certificates using the arithmeticgeometric-mean inequality.

• Reversed GP: The third approach involves converting an SGP problem into a Reversed GP problem and then applying iterative monomial approximation. Shen et al. [START_REF] Shen | A signomial programming approach for binary image restoration by penalized least squares[END_REF] minimize a SGP problem by converting it into a reversed GP problem, which is then solved using standard GP solvers. Aliabadi et al. [START_REF] Aliabadi | Conference: 13th International Conference on Industrial Engineering (IIEC 2017)At[END_REF] consider a nonlinear model that incorporates joint partial delayed payments, pricing, and marketing strategies in a supply chain consisting of a retailer and multiple customers. The problems are formulated using a constrained SGP with two levels of difficulty. To solve the proposed models, the researchers transform them into a reversed constraint programming framework and obtain optimal solutions in closed forms for each case.

As an example, we can take the problem of maximizing the total capacity of the wireless network from [START_REF] Haleem | On the capacity and transmission techniques of massive mimo systems[END_REF], which can be stated as follows

max p∈ K ++ min i∈U log 2   1 + p i |g H i g i | 2 j∈U ,j̸ =i p j |g H i g j | 2 +|σ i | 2    , (2.28)
s.t

P min ≤ p i ≤ P max , ∀i ∈ U, (2.29)
By taking

a ij = |g H i g j | 2 and b i = |σ i | 2 ,
a ij p j + b i , ∀i ∈ U, (2.31)
P min ≤ p i ≤ P max , ∀i ∈ U.
(2.32)

The search for the global optimal solution of problem (2.30)-(2.32) is challenging due to its non-convexity. However, in an effort to find the global optimal solution, one can still employ any of the aforementioned methods mentioned earlier to solve the problem.

Quasi geometric programming

Quasi geometric programming (QGP) is a special class of non-convex optimization that can be transformed into geometric programming by holding certain variables constant. The exploration of this specific type of optimization problem, which can be nonlinear and possibly non-smooth, is motivated by the observation that numerous engineering problems can be effectively formulated or approximated as QGPs. Nevertheless, solving a QGP remains challenging due to its inherent non-convex nature. This extended version of GP was introduced in 2010 by Toscano and Lyonnet [START_REF] Toscano | Evolutionary algorithms for solving quasi geometric programming problems[END_REF].

A nonlinear optimization problem is referred to as a quasi-geometric programming problem if it can be formulated in the following form

min x∈ p ++ ,ζ∈ q ++ ϕ 0 (x, ζ) -Q 0 (ζ), (2.33)
s.t ϕ k (x, ζ) ≤ Q k (ζ), k = 1, ..., K,
where the function

ϕ k (x, ζ), i = 0, 1, ..., K are posynomials. Problem (2.33) is called quasi geometric programming in posynomial form if Q k (ζ), k = 0, 1, ..., K are ratios of posyno- mial functions. Otherwise, if Q k (ζ), k = 0, 1, .
.., K are just assumed to be positive problem (2.33) is called quasi geometric programming in general form. It is crucial to note that in both of these cases, problem (2.33) cannot be converted into a standard form GP, and as a result, the problem is non-convex. Consequently, there is no approach available to swiftly find even a suboptimal solution using existing convex solvers. However, specific algorithms can be designed to discover suboptimal solutions for problem (2.33) in both posynomial or general form. We consider as example the following optimization problem from [START_REF] Qu | A new global optimization algorithm for signomial geometric programming via lagrangian relaxation[END_REF] min 0.5x

1 x 1 2 -x 1 -5x 1 2 , (2.34)
s.t 0.01x 2 x 1 3 + 0.01x 2 + 0.0005x 1 x 3 ≤ 1, 70 ≤ x 1 ≤ 150, 1 ≤ x 2 ≤ 30, 0.5 ≤ x 3 ≤ 21.
Problem 2.34 can be rewritten as a quasi geometric program (2.35) and can be solved using a global optimization algorithm via Lagrangian relaxation.

min λ -1 (2.35) s.t λ + 0.5x 1 x 1 2 x 1 + 5x 1 2 ≤ 1, 0.01x 2 x 1 3 + 0.01x 2 + 0.0005x 1 x 3 ≤ 1, 70 ≤ x 1 ≤ 150, 1 ≤ x 2 ≤ 30, 0.5 ≤ x 3 ≤ 21.

. Chance Constrained Optimization

In today's complex and dynamic world, decision-makers face a multitude of uncertainties and risks that can significantly impact the outcomes of their choices. Traditional optimization techniques often assume deterministic inputs and fail to capture the inherent variability and randomness present in many real-world scenarios. Chance constrained optimization, also known as probabilistic optimization or stochastic programming, offers a principled approach to handle uncertainties and make decisions that are robust and reliable. In this section, we provide an extensive survey on chance-constrained programming, including its concepts, results, and applications. We start by delving into the fundamental concepts underlying chance-constrained programming, offering a comprehensive discussion on the topic. Following that, we summarize the established results for various chance-constrained models. Subsequently, we shift our focus to conducting a survey of real-world applications where chance-constrained programming has been employed.

Chance constrained optimization was introduced for the first time by Charnes, Cooper and Symonds [START_REF] Charnes | Cost horizons and certainty equivalents: An approach to stochastic programming of heating oil[END_REF]. Since then, it has been applied to several real-world problems as aircraft systems [START_REF] Clark | Gradient-Based Optimization of Time-Dependent Aircraft Subsystems under Uncertainty[END_REF], machine learning [START_REF] Ang | A constrained teaching-learning-based optimization with modified learning phases for constrained optimization[END_REF][START_REF]Machine Learning and Wireless Communications[END_REF], water distribution networks [START_REF] Köker | Chance constrained optimization of booster chlorination in water distribution networks[END_REF][START_REF] Lee | Chance-constrained optimization for active distribution networks with virtual power lines[END_REF], radiation therapy [START_REF] Lodwick | Optimization under uncertainty: Methods and applications in radiation therapy[END_REF], power systems [START_REF] Roald | Power systems optimization under uncertainty: A review of methods and applications[END_REF], critical infrastructure protection [START_REF] Singh | Chance-constrained optimization for critical infrastructure protection[END_REF] and supply chain finance [START_REF] Tseng | Decision-making model for sustainable supply chain finance under uncertainties[END_REF].

A general form of a chance-constrained optimization problem is given below.

min x∈X f (x), subject to P (h(x, ζ) ≤ 0) ≥ p, (2.36)
where

f is the objective function, h = (h 1 , h 2 , ..., h m ) is an m-dimensional function, X ∈ R n is a deterministic set, x is the decision variable, ζ is a random variable and p ∈ [0, 1] is a
given probability level.

. Individual and joint chance constraints

The constraint in (2.36) is called a joint chance constraint since the constraints

h 1 (x, ζ) ≤ 0, h 2 (x, ζ) ≤ 0,..., h m (x, ζ) ≤ 0,
have to be satisfed simultaneously with probability p. Alternatively, each of the following m constraints can be referred to as an individual chance constraint

P (h 1 (x, ζ) ≤ 0) ≥ p 1 , P (h 2 (x, ζ) ≤ 0) ≥ p 2 , ..., P (h m (x, ζ) ≤ 0) ≥ p m (2.

37)

. Chance constraints were initially formulated as individual chance constraints [START_REF] Charnes | Cost horizons and certainty equivalents: An approach to stochastic programming of heating oil[END_REF], and further advancements were made in [START_REF] Charnes | Chance-constrained programming[END_REF][START_REF] Charnes | Chance constraints and normal deviates[END_REF][START_REF] Charnes | Deterministic equivalents for optimizing and satisficing under chance constraints[END_REF]. The concept of Joint Chance Constraints was first introduced by Miller and Wagner [START_REF] Miller | Chance constrained programming with joint constraints[END_REF] where they explore the mathematical properties of chance constrained programming problems, specifically focusing on situations where the constraint is on the joint probability of a multivariate random event. In the early 1970s, Prékopa made significant contributions [START_REF] Prekopa | ON PROBABILISTIC CONSTRAINED PROGRAMMING[END_REF][START_REF] Prékopa | A class of stochastic programming decision problems[END_REF], including the development of comprehensive convexity results for the feasible set of (2.36).

Generally, joint chance constraints are more adequate to use especially when multiple inequalities are used to define a system they lead to more robust solutions. Nevertheless, they are generally more difficult to solve. There exist several methods enabling to convert joint constraints to individual ones.

A joint chance constrained can be written as a set of individual constraints. Using Boole's inequality we have

P m i=1 h i (x, ζ) ≤ 0 ≤ m i=1 P (h i (x, ζ) ≤ 0) . (2.38) If m i=1 p i ≤ p,
then a feasible solution to (2.37) is also a feasible solution to the joint constraint in (2.36). However, the choice of parameters p i is generally a difficult exercise.

Insightful observations regarding the limitations of this approach can be found in [START_REF] Chen | From cvar to uncertainty set: Implications in joint chance-constrained optimization[END_REF].

Tractable approaches for addressing joint chance constrained programming problems can be categorized into two groups : sampling-based approximations and analytical safe approximations. Luedtke and Ahmed [START_REF] Luedtke | A sample approximation approach for optimization with probabilistic constraints[END_REF] study sample approximations of optimization problems with probabilistic constraints. Jeff Hong et al. [START_REF] Hong | Sequential convex approximations to joint chance constrained programs: A monte carlo approach[END_REF] propose sequential convex approximations to joint chance constraints based on a Monte Carlo approach. Blackmore et al. [START_REF] Blackmore | A probabilistic particle-control approximation of chance-constrained stochastic predictive control[END_REF] present a novel method for chance-constrained predictive stochastic control of dynamic systems. The method approximates the distribution of the system state using a finite number of particles. Li et al. [START_REF] Li | An analytical safe approximation to joint chance-constrained programming with additive gaussian noises[END_REF] present a safe approximation to joint chance-constrained programming, where the constraint functions are additively dependent on a normally-distributed random vector.

. An overview of solutions to chance-constrained optimization

Chance-constrained optimization models often present computational challenges, making them difficult to solve efficiently. Various existing approaches have been developed to address this issue, including decomposition methods, sample average approximations, and the formulation of deterministic equivalents.

Decomposition methods

Luedtke [START_REF] Luedtke | A branch-and-cut decomposition algorithm for solving chanceconstrained mathematical programs with finite support[END_REF] introduces a novel approach for the exact solution of chance-constrained mathematical programs that involve discrete distributions with finite support and random polyhedral constraints. The proposed approach combines decomposition techniques and integer programming methods to handle the complexity of the problem. Bai et al. Adam et al. [START_REF] Adam | Solving joint chance constrained problems using regularization and benders' decomposition[END_REF] addresse the problem of stochastic programs with joint chance constraints that involve discrete random distributions and propose a reformulation technique that introduces auxiliary variables. By adding these variables, they transform the original problem into a new formulation. The transformed problem, however, has a non-regular feasible set, which can pose challenges in finding optimal solutions. To address this issue, the authors propose a regularization approach. This involves expanding the feasible set by introducing additional constraints, which increases the solution space. To solve the regularized problem, the authors employ an iterative approach. They solve a master problem iteratively while incorporating Benders' cuts from a slave problem. This iterative process helps refine the solution and improve the convergence towards an optimal solution. [START_REF] Bai | An augmented lagrangian decomposition method for chance-constrained optimization problems[END_REF] present an augmented Lagrangian decomposition method that aims to find high-quality feasible solutions for complex optimization problems, including nonconvex chance-constrained problems. The proposed method differs from existing augmented Lagrangian approaches by allowing randomness to appear in both the left-hand-side matrix and the right-hand-side vector of the chance constraint. The method involves solving a convex subproblem and a 0-1 knapsack subproblem at each iteration. The special structure of the chance constraint enables the efficient computation of the 0-1 knapsack problem in quasi-linear time, which helps keep the computation for discrete optimization subproblems at a relatively low level.

Sample average approximations

The ample average approximation (SAA) method involves generating a set of scenario samples from the underlying random distribution and computing the empirical average of the objective function and constraints based on these samples. Chunlin and Liu [START_REF] Chunlin | Sample average approximation method for chance constrained stochastic programming in transportation model of emergency management[END_REF] propose the use of Conditional Value at Risk (CVaR) to approximate the chance constraint in stochastic programming problems. Specifically, the chance constraint is approximated by its CVaR counterpart, which provides a measure of the expected loss beyond a certain confidence level. To solve the resulting chance-constrained stochastic programming problem, the authors employ the SAA method. Cheng et al. [START_REF] Cheng | Partial sample average approximation method for chance constrained problems[END_REF] introduce a novel samplingbased method to solve chance constrained programs, which offers several advantages over traditional approaches. The key innovation is that the proposed method formulates an approximation problem that involves only continuous variables, as opposed to the standard SAA formulation that includes binary variables. Kleywegt et al. [START_REF] Kleywegt | The sample average approximation method for stochastic discrete optimization[END_REF] investigate a Monte Carlo simulation-based approach for solving stochastic discrete optimization problems. The fundamental concept behind such methods is to generate a random sample and approximate the expected value function using the corresponding sample average function. By solving the resulting sample average optimization problem, the approach aims to approximate the optimal solution of the original stochastic problem. This process is iterated multiple times until a predefined stopping criterion is met.

Deterministic equivalent formulations

Nemirovski and Shapiro [START_REF] Nemirovski | Convex approximations of chance constrained programs[END_REF] address the problem of chance constrained optimization, where the objective is to minimize a convex function subject to a set of convex constraints that are randomly perturbed. The goal is to develop a computationally tractable approximation of this problem, specifically a deterministic optimization program that can be efficiently solved and whose feasible set is contained within the original chance constrained problem. A class of convex conservative approximations was proposed for the chance constrained problem. These approximations provide a computationally tractable representation of the problem while ensuring that the feasible solutions remain within the original chance constrained set. Xie and Ahmed [START_REF] Xie | On deterministic reformulations of distributionally robust joint chance constrained optimization problems[END_REF] investigate deterministic reformulations for joint chance constraints and explores the conditions under which these reformulations are convex. Gopalakrishnan et al. [START_REF] Gopalakrishnan | Solving chanceconstrained optimization under nonparametric uncertainty through hilbert space embedding[END_REF] present a novel algorithm that efficiently solves a class of chance-constrained optimization problems under nonparametric uncertainty. The algorithm leverages the concept of Reproducing Kernel Hilbert Space (RKHS) to represent arbitrary distributions as functions allowing to formulate the chanceconstrained optimization problem as a minimization of the distance between a desired distribution and the distribution of the constraint functions within the RKHS framework.

. Application : Chance Constrained Portfolio Problem Description

Portfolio optimization is a widely researched problem in the field of finance. It involves selecting the optimal allocation of investments among a set of K possible companies. The objective is to maximize the revenue generated from these investments.

In the literature, various formulations of the portfolio optimization problem have been studied. One common approach is the risk-averse formulation, where the goal is to minimize the variance of the total return rate. Some recent applications and methods can be found in [START_REF] Faria | Scenario generation and riskaverse stochastic portfolio optimization applied to offshore renewable energy technologies[END_REF][START_REF] Edirisinghe | Meanndash;variance portfolio efficiency under leverage aversion and trading impact[END_REF][START_REF] Chen | A sparse chance constrained portfolio selection model with multiple constraints[END_REF][START_REF] Sun | Chance-constrained optimization for pension fund portfolios in the presence of default risk[END_REF].

In our model, we make certain assumptions about the problem. These assumptions will serve as the foundation for our analysis and optimization approach. We suppose that we have exact knowledge of the return rate distributions. Additionally, there are no extra costs associated with making investments, allowing for flexibility in allocating funds among companies without incurring differential costs. We also include a risk-free investment option with a return rate of 0. To address risk aversion, we formulate the problem using chance constraints, ensuring that the potential loss does not fall below a certain threshold. Our objective is to maximize the expected revenue. We consider a system represented by a sequence of random variables R i that represent the return rates of K possible investments, and our aim is to find the optimal solution to this system.

max K i=1 E[R i ]x i , (2.39) s.t P K i=1 R i x i ≤ -η ≤ ϵ (2.40) K i=1 x i ≤ 1, x i ≥ 0, i = 1, ..., K.
(2.41) Under the assumption that η > 0 and ϵ < 0.5, finding a solution to the aforementioned problem allows us to establish an investment strategy that maximizes the expected return rate while ensuring the probability of experiencing losses exceeding η remains below ϵ.

By optimizing the portfolio allocation based on these criteria, we aim to strike a balance between maximizing potential returns and managing risk within the specified threshold.

We consider that the return rates of investments follow a normal distribution and we have complete knowledge of the mean µ i and standard deviation σ i for all distributions, i.e., R i ∼ N (µ i , σ i ) for i = 1, ..., K. Additionally, we assume that the return rates of different investments are independent of each other. We have then

K i=1 R i x i ∼ N   K i=1 µ i x i , K i=1 σ 2 i x 2 i   .
(2.42)

We have then

P K i=1 R i x i ≤ -η = P    K i=1 µ i x i - K i=1 R i x i K i=1 σ 2 i x 2 i ≥ η+ K i=1 R i x i K i=1 σ 2 i x 2 i    = 1-ϕ    η+ K i=1 µ i x i K i=1 σ 2 i x 2 i   ,
where ϕ is the cumulative probability function of the standard normal distribution. Con-straint (2.40) can be written equivalently as

1 -ϕ       η + K i=1 µ i x i K i=1 σ 2 i x 2 i       ≤ ϵ, (2.43)
which is finally equivalent to

ϕ -1 (1 -ϵ) K i=1 σ 2 i x 2 i - K i=1 µ i x i ≤ η.
(2.44)

Problem (2.39)-(2.41) becomes then max K i=1 µ i x i , (2.45)
s.t ϕ -1 (1 -ϵ) K i=1 σ 2 i x 2 i - K i=1 µ i x i ≤ η, (2.46)
K i=1 x i ≤ 1, x i ≥ 0, i = 1, ..., K. (2.47)
Problem (2.45)-(2.47) being a convex second-order cone program implies that it can be solved using existing methods developed for solving convex optimization problems.

In the upcoming Chapter, we will provide a comprehensive theoretical overview of neurodynamic neural networks and biconvex optimization. This detailed exploration aims to enhance the understanding of the thesis by establishing a solid foundation in these key areas of study.

-Neurodynamic Neural Networks and Biconvex Optimization

. Dynamical Neural Networks

. Introduction

Throughout one's lifetime, the brain's networks are constantly abuzz with waves of activity. Countless signals harmonize and oscillate, forming the underlying functional architecture that shapes every facet of existence. Emotions, stress levels, musical preferences, aspirations, and dreams all stem from the dynamic behavior of these networks, complemented by adaptive memory systems that strive to better align with one's ever-changing surroundings.

The vast scale and intricate integration of networks in the human brain present formidable challenges in understanding whether traditional forms of computation occur within them and how. Evidence points to information processing taking place at various levels, encompassing modulatory proteins within individual cells, cortical microcircuits, and expansive functional networks spanning the entire brain. However, the experimental grasp of the brain's workings advances slowly. Thankfully, resourceful engineers have devised algorithms that partially simulate aspects of these networks. While no single model can fully encapsulate the sheer complexity and behavior of the human brain, these tools offer researchers a valuable window into how information might be computed and ultimately represented within the activity of distributed networks.

A significant category of logical problems that stem from real-world scenarios can be framed as optimization problems, which can be qualitatively understood as a quest for the optimal solution. These problems are prevalent in fields such as engineering and commerce, as well as in perceptual challenges that require swift resolution by nervous systems.

The fundamental concept behind the neural network approach for optimization is to create a nonnegative energy function and establish a dynamic system that mimics an artificial neural network. This dynamic system typically takes the form of first-order ordinary differential equations. The objective is for the dynamic system to converge towards a static state or equilibrium point, which corresponds to the solution of the underlying optimization problem, starting from an initial point. One notable advantage of using neural networks for solving optimization problems is their hardware implementability. In other words, these neural networks can be implemented using integrated circuits. Moreover, neural networks designed for circuit implementation offer real-time processing capabilities, further enhancing their practical applicability.

In 1985, Hopfield and Tank introduced the first neural network model for solving the Traveling Salesman Problem [START_REF] Hopfield | Neural computation of decisions in optimization problems[END_REF]. Since then, numerous neural network models have been proposed and utilized in both linear programming and nonlinear programming do-mains. Kennedy and Chua [START_REF] Kennedy | Neural networks for nonlinear programming[END_REF] extend the Tank-Hopfield network by incorporating a finite penalty parameter into a neural network framework. This extension demonstrated the versatility and adaptability of neural networks in addressing a broader range of optimization problems beyond the initial applications of the Tank-Hopfield network. In the pursuit of finding exact solutions, researchers have continued to develop additional neural network models for optimization. These new neural network approaches build upon the foundation established by earlier models and aim to improve the accuracy and efficiency of optimization solutions.

Hosseini et al. [START_REF] Hosseini | A recurrent neural network for solving a class of generalized convex optimization problems[END_REF] propose a penalty-based recurrent neural network designed to address a specific class of constrained optimization problems featuring generalized convex objective functions. The proposed model is characterized by a straightforward structure represented through a differential inclusion. It is applicable to a wide range of nonsmooth optimization problems that incorporate affine equality and convex inequality constraints. Notably, the objective function should be regular and pseudoconvex within the problem's feasible region for successful application. Qin and al. [START_REF] Qin | A neurodynamic approach to convex optimization problems with general constraint[END_REF] present a neurodynamic approach that utilizes a recurrent neural network to solve convex optimization problems with general constraints. The key contribution of the proposed approach is that for any initial point, the state of the neural network reaches the constraint set within a finite time and ultimately converges to an optimal solution of the convex optimization problem. Qin and Xue [START_REF] Qin | A two-layer recurrent neural network for nonsmooth convex optimization problems[END_REF] introduce a novel approach to tackle nonsmooth convex optimization problems with convex inequality and linear equality constraints using a two-layer recurrent neural network. In contrast to existing neural network models, the proposed approach offers several advantages. Firstly, it exhibits a low model complexity, making it computationally efficient and easier to implement. Additionally, the proposed neural network does not rely on penalty parameters, eliminating the need for parameter tuning and improving the overall simplicity of the optimization process. Liu and Zhou [START_REF] Liu | A one-layer recurrent neural network for non-smooth convex optimization subject to linear inequality constraints[END_REF] present a one-layer recurrent network designed specifically for solving non-smooth convex optimization problems subject to linear inequality constraints. Liu and al. [START_REF] Liu | A one-layer recurrent neural network for constrained pseudoconvex optimization and its application for dynamic portfolio optimization[END_REF] introduce a one-layer recurrent neural network as a solution for pseudoconvex optimization problems that involve linear equality and bound constraints. In comparison to existing neural networks used for optimization, such as projection neural networks, the proposed model exhibits enhanced capabilities in solving a broader range of pseudoconvex optimization problems with equality and bound constraints. Furthermore, it can handle constrained fractional programming problems as a special case.

Recurrent neural networks have proven to be highly effective in solving a wide range of real-world problems. Liu et al. [START_REF] Liu | A one-layer recurrent neural network for constrained pseudoconvex optimization and its application for dynamic portfolio optimization[END_REF] propose a decision-making model based on a recurrent neural network for dynamic portfolio optimization. Yan et al. [START_REF] Yan | Recurrent neural network-based nonlinear optimization for braking control of electric vehicles[END_REF] use recurrent neural networks to address the problem of electro-hydraulic braking (EHB) force allocation for electric vehicles. Kerboua and Kelaiaia [START_REF] Kerboua | Recurrent neural network optimization for wind turbine condition prognosis[END_REF] employ recurrent neural networks for the health prognosis of wind turbine operations using vibration time series data. The study explores various memory cell variations, including Long Short Time Memory (LSTM), Bilateral LSTM (BiLSTM), and Gated Recurrent Unit (GRU), integrated into different net-work architectures. Urias et al. [START_REF] Urias | Electrical microgrid optimization via a new recurrent neural network[END_REF] design a recurrent neural network for optimization in the context of optimal operation of an electrical microgrid. The microgrid is connected to the utility grid and includes various components such as batteries for energy storage and supply, as well as an electric car. The main objective of the proposed neural network is to determine the optimal power allocation over a one-week time horizon for wind, solar, and battery systems, taking into account the electric car, with the goal of minimizing the power obtained from the utility grid and maximizing the power generated from renewable energy sources.

The dynamical recurrent neural networks are based on ODE systems. We present in the following subsection some necessary theoretical results to study the stability of ODE systems.

. Stability in differential equations

Let consider the following differential equation

ẋ(t) = f (x(t)), x(t 0 ) = x0 ∈ R n . (3.1)
We present some classical results that establish the existence and uniqueness of a solution for (3.1) [START_REF] Reiszig | lefschetz, stability by liapunov's direct method with applications. vii + 134 s. new york/london 1961. academic press. preis geb. $ 5[END_REF].

Theorem 1. Assuming that f is a continuous mapping from R n to R n , the following conclusions can be made regarding the existence and uniqueness of solutions to equation (3.1).

• Existence: For any given t 0 ≥ 0 and initial condition x 0 ∈ R n , there exists a local solution x(t) defined on an interval [t 0 , s) for some s > t 0 . This means there is a time interval where the solution exists and is well-defined.

• Uniqueness: If f is locally Lipschitz continuous at x 0 , then the solution to equation (3.1) is unique within the local interval [t 0 , s). This means there is only one possible solution for a given initial condition in this interval.

• Maximal Interval of Existence: The maximal interval of existence, denoted by [t 0 , s(x 0 )),

is the largest time interval on which a local solution exists for a given initial condition x 0 . If a local solution cannot be extended to a larger interval [t 0 , s), then the local interval [t 0 , s) is considered the maximal interval of existence.

Theorem 2. Assuming that f is a continuous mapping from R n to R n , if x(t), t ∈ [t 0 , s(x 0 )), is a maximal solution such s(x 0 ) > ∞, then it implies the following lim t-→s(x 0 ) ||x(t)||= +∞. (3.2) Definition 1. A point x * ∈ R n is called an equilibrium point of (3.1) if f (x * ) = 0.
Definition 2. Let x(t) be a solution of equation (3.1). An isolated equilibrium point x * is Lyapunov stable if, for any initial condition x 0 = x(t 0 ) and any scalar ϵ > 0, there exists a δ > 0 such that if ||x(t0) -x * ||≤ δ, then ||x(t) -x * ||≤ ϵ for all t ≥ t 0 .

In other words, if we can find a small enough neighborhood around the equilibrium point such that any initial condition within that neighborhood will result in the solution staying within another specified neighborhood around the equilibrium point for all future time, then the equilibrium point is said to be Lyapunov stable. This stability condition ensures that the solution does not diverge too far from the equilibrium point for sufficiently small perturbations in the initial condition. Definition 3. An isolated equilibrium point x * is said to be asymptotically stable if in addition to being Lyapunov stable, it also has the property that if ||x(t0) -

x * ||≤ δ then lim t-→∞ x(t) = x * .
In other words, if we can find a neighborhood around the equilibrium point such that any initial condition within that neighborhood will result in the solution approaching the equilibrium point as time goes to infinity, then the equilibrium point is said to be asymptotically stable. This stability property ensures that the system converges to a stable state over time.

Definition 4.

Let Ω ⊆ R n be an open neighborhood of x. A continuously differentiable function η : R n -→ R is said to be a Lyapunov function at the state x (over the set Ω) for the differential equation (3.1) if the following conditions hold

• η(x) = 0, η(x) > 0 for x ∈ Ω, x ̸ = x. • dη(x(t)) dt
≤ 0, for all x ∈ Ω Remark 1. The above conditions define the Lyapunov function and characterize it as an energy function for the system described by the differential equation (3.1). The Lyapunov function provides a measure of the system's energy or potential, and its derivative describes the change in energy over time. By analyzing the behavior of the Lyapunov function and its derivative, we can determine the stability properties of the system and assess whether the equilibrium point x* is stable or not.

Theorem 3. We have the following results

• An isolated equilibrium point x * is Lyapunov stable if there exists a Lyapunov function over some neighborhood Ω of x * .

• An isolated equilibrium point x * is asymptotically stable if there exists a Lyapunov function over some neighborhood Ω of x * satisfying dη(x(t)) dt < 0, for all x(t) ∈ Ω, x(t) ̸ = x * .

The neurodynamic system approach is a powerful technique for tackling optimization problems. It leverages artificial recurrent neural networks to convert these problems into dynamic systems represented by first-order differential equations. By initializing the system at a specific point, it is expected to converge towards a static state or equilibrium point, which corresponds to the solution of the original optimization problem.

One notable advantage of neural networks for optimization problems is their parallel information processing capability. Neural networks can simultaneously evaluate multiple inputs and compute corresponding outputs, thanks to their inherent parallel structure. This parallelism enables efficient and concurrent information processing, resulting in faster optimization performance than sequential algorithms.

Furthermore, neural networks designed for optimization problems can be implemented in hardware using specialized technologies like integrated circuits or dedicated processing units. This hardware implementation takes advantage of the parallel nature of neural networks, further enhancing their computational speed and efficiency. By utilizing dedicated hardware resources, neural networks can be deployed in real-time applications or embedded systems, enabling efficient and rapid optimization across diverse domains.

To better comprehend the application of recurrent neural networks in optimization, we present two recurrent neural networks from the literature in the coming subsections to solve convex nonlinear optimization problems.

. A recurrent neural network for nonlinear convex optimization subject to nonlinear inequality constraints

Consider the following nonlinear convex program.

min f (x), (NCP)
s.t c(x) ≤ 0, x ≥ 0.
where x ∈ R n , f is convex and twice differentiable and c(x

) = [c 1 (x), c 2 (x), ..., c m (x)] T
is an m-dimensional vector-valued continuous function, where each component is convex and twice differentiable. A vector x is considered a feasible solution if it satisfies the constraints of (NCP). In particular, a feasible solution is deemed a regular point if the gradients of the objective function components, and the constraint functions, are linearly independent.

if a feasible regular point x exists, the Karush-Kuhn-Tucker (KKT) conditions for (NCP) can be expressed as follows

y ≥ 0, c(x) ≤ 0, x ≥ 0, (3.3)
∇f (x) + ∇c(x)y = 0, y T c(x) = 0. (3.4)
By utilizing the optimization projection technique, we can establish a relationship between the solution to (NCP) and the solution to a projection formulation given by

(x -α(∇f (x) + ∇c(x)y)) + -x = 0, (3.5)
(y + αc(x)) + -y = 0.

(3.6)

where (x) + = [(x 1 ) + , ..., (x n ) + ], (x i )+ = max(x i , 0) and α is a positive constant.

Theorem 4. [START_REF] Xia | A recurrent neural network for nonlinear convex optimization subject to nonlinear inequality constraints[END_REF] Let x be a feasible and regular point of (NCP). Then, x is an optimal solution to (NCP) if and only if there exists a nonzero m-dimensional vector y such that for any α ≥ 0 (3.5)-(3.6) holds for (x, y).

Based on the projection formulation presented in equation (3.5)-(3.6), a recurrent neural network model for solving (NCP) is proposed. The following dynamical equation describes the model [START_REF] Xia | A recurrent neural network for nonlinear convex optimization subject to nonlinear inequality constraints[END_REF].

dx dt = -x + (x -α(∇f (x) + ∇c(x)y)) + , (3.7

. A dynamic system model for solving convex nonlinear optimization problems

Let us consider the following general CNLP (Convex Nonlinear Programming) problem.

min f (x), (CNLP) s.t g(x) ≤ 0, h(x) = 0. where x ∈ R n , f : R n -→ R, g(x) = [g 1 (x), g 2 (x), ..., g m (x)]
T is an m-dimensional vectorvalued continuous function of n variables, and the functions f, g 1 , ..., g m are assumed to be convex and twice differentiable. The constraint function h(x) = Ax -b, where A ∈ R (l×n) and b ∈ R l , represents a linear equality constraint.

The Lagrange function of (CNLP) is given by [START_REF] Nazemi | A dynamic system model for solving convex nonlinear optimization problems[END_REF]

L(x, u, v) = f (x) + 1 2 m k=1 u 2 k g k (x) + l p=1 v p h p (x).
(3.9)

where u ∈ R m and v ∈ R l are the Lagrange multipliers. x ∈ R n is an optimal solution of (CNLP) if and only if there exist u ∈ R m and v ∈ R l such that the following Karush-Kuhn-Tucker (KKT) system is verified.

u ≥ 0, g(x) ≤ 0, u T g(x) = 0, (3.10)
∇f (x) + ∇g(x) T u + ∇h(x) T v = 0 (3.11) h(x) = 0.
(3.12) Lemma 7. If f and g k , k = 1, ..., m are all convex, then x * is an optimal solution of (CNLP) if and only if x * is KKT point of (CNLP). Now, consider the time-dependent variables x(t), u(t), and v(t). Nazemi [START_REF] Nazemi | A dynamic system model for solving convex nonlinear optimization problems[END_REF] designs the following neural network that settles at the saddle point of the Lagrangian function

L(x, u, v). dx dt = -∇ x L(x, u, v), (3.13)
du dt = ∇ u L(x, u, v), (3.14)
dv dt = ∇ v L(x, u, v). (3.15)
For the sake of simplicity the dynamical neural network (3.13)-(3.15) can be written as follows

dy dt = ϕ(y), (3.16
)

y(t 0 ) = y 0 , u(t 0 ) ̸ = 0.
(3.17)

where y = (x, u, v) T and y 0 is a given initial points. Nazemi [START_REF] Nazemi | A dynamic system model for solving convex nonlinear optimization problems[END_REF] proves the following results. globally converges to a KKT point y = (x, u, v) T , where x is the optimal solution of problem (CNLP), and u, v are the optimal solutions of its dual.

. Biconvex Optimization

. Introduction

Biconvex optimization refers to a class of optimization problems where the objective function is convex in each of its variables while keeping the other variables fixed. This type of optimization problem arises in various fields, including machine learning, signal processing, and information processing. Biconvex optimization models have been extensively studied due to their practical applications and mathematical properties. One popular approach for solving biconvex optimization problems is the alternating direction method of multipliers (ADMM), which is an iterative optimization algorithm that handles the separable structure of biconvex problems. ADMM has been widely used in biconvex optimization due to its efficiency and convergence guarantees.

In recent years, there have been several advancements and applications of biconvex optimization. In machine learning, biconvex optimization has been used for tasks such as matrix completion, collaborative filtering, and multi-view learning. In signal and information processing, biconvex optimization has been applied to problems like compressed sensing, image reconstruction, and channel estimation. Research in biconvex optimization has focused on developing efficient algorithms, understanding theoretical properties, and exploring applications in various domains. Different techniques, such as convex relaxation, non-convex relaxation, and convex-concave procedures, have been proposed to solve biconvex optimization problems. Additionally, research has investigated the impact of incorporating additional constraints, such as sparsity or low-rank structure, into biconvex optimization models. These extensions have led to improved performance and more accurate solutions in many practical scenarios.

. Definition and Properties

The study of convexity of functions is motivated by the property that if a function is convex, then every local minimum is guaranteed to be a global minimum. However, it is important to note that convexity is not a necessary condition for this local-global property.

Unfortunately, biconvexity, as defined, does not imply the local-global property. While biconvex functions possess certain desirable properties, they do not guarantee that all local minima will be global minima. Therefore, the study and understanding of biconvexity require additional considerations beyond the local-global property. Definition 5. Let X ⊆ R n and Y ⊆ R m be two non-empty, convex sets and let B ⊆ X × Y.

We define

B x = {y ∈ Y|(x, y) ∈ B} (3.18) and B y = {x ∈ X |(x, y) ∈ B}. (3.19)
The set B is said to be biconvex on X × Y, if B x is biconvex for every x ∈ X and B y is biconvex for every y ∈ Y.

Definition 6. A function f : R n × R m -→ R is called a biconvex function on R n × R m , if f x (.) := f (x, .) : R m -→ R (3.20) is convex on R m for every fixed x ∈ R n and f y (.) := f (., y) : R n -→ R (3.21)
is convex on R n for every fixed y ∈ R m . Theorem 12. [START_REF] Goh | Biaffine matrix inequality properties and computational methods[END_REF] Let X ⊆ R n and Y ⊆ R m , a function f (x, y) is biconvex on X × Y if and only if for all quadruples (x 1 , y 1 ), (x 1 , y 2 ), (x 2 , y 1 ), (x 2 , y 2 ) ∈ X × Y and for every

λ, β ∈ [0, 1] f ((1 -λ)x 1 + λx 2 , (1 -β)y 1 + βy 2 ) ≤ (1-λ)(1-β)f (x 1 , y 1 )+β(1-λ)f (x 1 , y 2 )+λ(1-β)f (x 2 , y 1 )+λβf (x 2 , y 2 Theorem 13. [54] If f is biconvex, then its level sets L c = {(x, y) ∈ X × Y|f (x, y) ≤ c} (3.22)
are all convex for every c ∈ R. Lemma 14. [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF] Let f and g be two biconvex functions on X × Y and let µ ∈ R + , then the functions f + g and µf are biconvex. Lemma 15. [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF] Let f a biconvex function on X × Y and let ϕ : R -→ R a e a convex, non-decreasing function. Then the function ϕ(f (., .)) is biconvex on X × Y.

. Biconvex optimization problems

The existing research highlights the challenging nature of solving biconvex programming problems, making it a highly valuable area of study that necessitates the exploration of new theories and solution methods. In this section, we present two significant theoretical results related to the partial optimum of biconvex programming, utilizing the objective penalty function.

The first result establishes the equivalence between the partial Karush-Kuhn-Tucker (KKT) condition and the partial exactness property of the objective penalty function in biconvex programming. This finding provides an important insight into the relationship between the KKT condition and the behavior of the objective penalty function, enabling a deeper understanding of the optimization process.

The second result demonstrates the equivalence between the partial stability condition and the partial exactness property of the objective penalty function in biconvex programming. This result sheds light on the significance of the stability condition in relation to the objective penalty function, offering valuable insights into the convergence behavior of algorithms aiming to solve the partial optimum of biconvex programming.

These theoretical results offer crucial guarantees for the convergence of algorithms designed to tackle the partial optimum of biconvex programming problems. By establishing the equivalence between key conditions and properties of the objective penalty function, this research provides a solid foundation for developing efficient and effective algorithms that can reliably converge towards optimal solutions for biconvex programming.

We consider the following biconvex programming min f (x, y), subject to

g i (x, y) ≤ 0, i = 1, ..., K, (BP)
where f g i , i = 1, ..., K are biconvex functions on R n × R m . the feasible set for (BP) is denoted by

B = {(x, y) ∈ R n × R m | g i (x, y) ≤ 0, i = 1, ..., m} .
We define

B x = {y ∈ R m | g i (x, y) ≤ 0, i = 1, ..., m} and B y = {x ∈ R n | g i (x, y) ≤ 0, i = 1, ..., m} .
We additionally define two suboptimal problems

min f (x, y), subject to x ∈ B y , (BP1)
and

min f (x, y), subject to y ∈ B x . (BP2) Definition 7. Let (x * , y * ) ∈ B, (x * , y * ) is called a partial optimum of (BP) if f (x * , y * ) ≤ f (x, y * ), ∀x ∈ B y and f (x * , y * ) ≤ f (x * , y), ∀y ∈ B x .
. Remark 2. It is straightforward that an optimal solution for (BP) is a partial optimum of (BP).

Remark 3. In biconvex optimization problems, there always exists a partial optimal solution [START_REF] Jiang | Partial exactness for the penalty function of biconvex programming[END_REF].

Definition 8. Let (x * , y * ) ∈ B, if there exist α i , i = 1, ..., K such that ∇f (x * , y * ) + K i=1 α i ∇g i (x * , y * ) = 0, α i g i (x * , y * ) = 0, α i ≥ 0, i = 1, ..., K, then (x * , y * ) is called a KKT point of (BP). Definition 9. Let (x * , y * ) ∈ B, if there exist α 1 i , α 2 i , i = 1, ..., K such that ∇ x f (x * , y * ) + K i=1 α 1 i ∇ x g i (x * , y * ) = 0, ∇ y f (x * , y * ) + K i=1 α 2 i ∇ y g i (x * , y * ) = 0, α 1 i g i (x * , y * ) = 0, α 2 i g i (x * , y * ) = 0, α 1 i ≥ 0, α 2 i ≥ 0, i = 1, ..., K, then (x * , y * ) is called a partial KKT point of (BP).
Definition 10. Let (x * , y * ) ∈ B, if there exist (x, ỹ) such that g i (x * , ỹ) < 0, g i (x, y * ) < 0, i = 1, ..., K, then (BP) is said to satisfy partial Slater constraint qualification at (x * , y * ). Theorem 16. [START_REF] Jiang | Partial exactness for the penalty function of biconvex programming[END_REF] Let (x * , y * ) ∈ B. If (BP) satisfies partial Slater constraint qualification at (x * , y * ), (x * , y * ) is a partial optimum of (BP) if and only if (x * , y * ) is a partial KKT point of (BP). Corollary 17. [START_REF] Jiang | Partial exactness for the penalty function of biconvex programming[END_REF] Let (x * , y * ) ∈ B be a partial optimum of (BP). If (BP) satisfies partial Slater constraint qualification at (x * , y * ), then (x * , y * ) ∈ B is a KKT point of (BP) if and only if

α 1 i = α 2 i , i = 1, ..., K.

. Algorithms

In this subsection, we focus on various methods and algorithms specifically designed to solve biconvex minimization problems represented by (BP). These methods take advantage of the inherent biconvex structure of the problem to achieve efficient solutions. For each solution approach, we provide concise algorithmic descriptions and discuss their convergence properties and limitations.

Alternate convex search

Firstly, we explore the Alternate Convex Search (ACS) method, which is a special case of Block-Relaxation Methods. This approach involves iteratively optimizing each block of variables while keeping the other blocks fixed. By alternating between the blocks and iteratively refining the solution, this method aims to converge to a biconvex minimizer. We delve into the details of this algorithm, discuss its convergence properties, and highlight any limitations or challenges associated with its implementation. The pseudocode of the alternate convex search algorithm is given as follows.

Algorithm 1

The alternate convex search 1-Select an arbitrary starting point, denoted as z 0 = (x 0 , y 0 ) ∈ B. 2-Initialize the iteration counter i = 0. 3-For a fixed y i , solve the convex problem (BP1). if There exists an optimal solution x * ∈ B y i for (BP1) then

Set

x i+1 = x * , else STOP. end if

4-For a fixed x i , solve the convex problem (BP2).

if There exists an optimal solution y * ∈ B x i for (BP2) then

Set

y i+1 = y * , else STOP. end if 5-Set z i+1 = (x i+1 , y i+1 ). if A stopping criterion is satisfied, i.e, ||z i+1 -z i ||≤ ϵ then STOP. else i = i + 1.

end if

Theorem 18. [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF] Assuming that B ⊆ R n × R m , where f : B -→ R is bounded from below, and the optimization problems (BP1) and (BP2) are solvable, we can conclude that the sequence {f (z i )} generated Algorithm 1 converges monotonically.

The statement of Theorem 27 may appear weak because the boundedness of the objective function f guarantees the convergence of the sequence {f (z i )} generated by the ACS algorithm, but it does not automatically guarantee the convergence of the sequence {z i } itself. It is possible to encounter scenarios where the sequence {f (z i )} converges while the sequence {z i } diverges, as shown in Example 4.3 in [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF].

In other words, while the objective function values approach a limit as the ACS algorithm iterates, the sequence of points {z i } may not converge to a specific point or region in the optimization space. This behavior can occur even for biconvex functions, where the convergence of {f (z i )} does not imply the convergence of {z i }.

The existence of such cases highlights the importance of considering additional conditions or assumptions to ensure the convergence of the sequence {z i }. Simply relying on the boundedness of f may not be sufficient to guarantee the convergence of the optimization variables.

Therefore, it is crucial to carefully examine the specific properties of the optimization problem, such as additional constraints or assumptions, to ascertain the convergence of the sequence {z i }.

Definition 11. [57] Let B ⊆ R n × R m , let z 1 = (x 1 , y 1 ) and z 2 = (x 2 , y 2 ) in B.
The map A : B -→ P(B), where P(B) is the power set of (B) defined by z 2 ∈ A(z 1 ) if and only if

f (x 2 , y 1 ) ≤ f (x, y 1 ) ∀x ∈ B y 1 and f (x 2 , y 2 ) ≤ f (x 2 , y) ∀y ∈ B x 2
is called the algorithmic map of the ACS algorithm.

In other words,A(z 1 ) = z 2 in Definition 21 implies that the ACS algorithm, based on the algorithmic map A, identifies the point z 2 as a candidate solution that satisfies the given conditions. This map plays a crucial role in guiding the ACS algorithm by determining the next point to explore in the optimization process.

By utilizing the algorithmic map A, the ACS algorithm dynamically selects points in B, allowing for an iterative exploration of the optimization space and providing information about the neighborhood regions where local minima are identified, helping to refine the search and improve the convergence properties of the ACS algorithm. Theorem 19. [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF] Let X ⊆ R n and Y ⊆ R m be closed sets and let f : X × Y -→ R be a continuous function. If the sequence {z i } generated by the ACS algorithm converges to z * ∈ X × Y , then z * is a partial optimum of (BP).

To sum up, to seek the global optimum of a biconvex minimization problem using the ACS algorithm, a multistart version of ACS can be employed, as suggested by Goh et al. [START_REF] Goh | Biaffine matrix inequality properties and computational methods[END_REF]. However, it is important to note that even with the multistart approach, there is no guarantee of finding the global optimum within a reasonable amount of time or being certain that the identified best minimum is indeed the global minimum. The multistart version of ACS involves running the ACS algorithm from multiple different starting points in order to explore different regions of the optimization space. By initiating the algorithm from various initial points, it increases the chances of discovering better solutions and potentially finding the global minimum. However, due to the inherent complexity of biconvex optimization problems, it is challenging to ensure global optimality. The presence of multiple local minima and the intricate interplay between the convex components make it difficult to guarantee that the identified solution is indeed the global minimum. Additionally, the computational time required to exhaustively search the entire optimization space can be prohibitive. Therefore, while the multistart version of ACS improves the likelihood of finding better solutions, it does not provide a foolproof guarantee of locating the global optimum within a reasonable timeframe. It is crucial to carefully analyze the problem characteristics, employ suitable search strategies, and evaluate the obtained results to make informed decisions regarding the optimality of the solution and the trade-off between computational resources and solution quality. A detailed discussion with examples can be found in [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF].

The global optimization algorithm

In this subsection, we present a comprehensive review of the Global Optimization Algorithm (GOP), developed by Floudas and Visweswaran [START_REF] Visweswaran | A global optimization algorithm (gop) for certain classes of nonconvex nlps-ii. application of theory and test problems[END_REF], for solving constrained biconvex minimization problems. This algorithm leverages the convex substructure of the problem through a primal-relaxed dual approach, building upon decomposition concepts introduced by Benders [START_REF] Benders | Partitioning procedures for solving mixed-variables programming problems[END_REF] and Geoffrion [START_REF] Geoffrion | Generalized benders decomposition[END_REF].

The GOP algorithm follows a two-step process similar to the ACS method. In the first step, the constrained problem is solved for a fixed value of the variable y, leading to an upper bound on the solution of the biconvex problem. This initial problem is referred to as the primal problem.

To obtain a lower bound for the solution, the algorithm utilizes duality theory and applies linear relaxation. The resulting relaxed dual problem is solved for every possible combination of bounds in a subset of the variables x.

By iteratively alternating between the primal and relaxed dual problems, the GOP algorithm achieves finite ϵ-convergence to the global optimum. This means that as the iterations progress, the algorithm converges towards the global optimum within a specified tolerance ϵ.

The primal-relaxed dual approach of the GOP algorithm combines decomposition techniques, duality theory, and linear relaxation to effectively explore the solution space of constrained biconvex minimization problems. Through this iterative process, the algorithm provides a robust framework for finding globally optimal solutions. The finite ϵ-convergence guarantees that the algorithm will eventually converge to the global optimum, although the convergence rate and the number of iterations required may vary depending on the specific problem instance and input parameters. Theorem 20. [START_REF] Floudas | Deterministic global optimization: Theory, methods, and applications[END_REF] If the conditions of Theorem 4.11 [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF] hold, then the GOP algorithm will terminate at the global optimum of the biconvex minimization problem.

Advantages of the GOP algorithm

• Convexity of the primal problem: The fact that the primal problem in the first step of each iteration is a convex problem is a major advantage. This means that every local optimum is also the global minimum of the subproblem. This property simplifies the optimization process and ensures that the obtained solution is globally optimal within the subproblem.

• Simplification of constraints: The set of constraints for the convex subproblem often simplifies to linear or quadratic constraints in the x-variables. This simplification enables the use of conventional non-linear local optimization solvers to efficiently solve the primal problem.

• Reduction of variables: The relaxed dual problem only needs to be solved for the connected x-variables, reducing the number of variables involved. This reduction in variables can lead to improved computational efficiency and reduced computational resources.

Drawbacks of the GOP algorithm

• Number of subproblems: In each iteration of the algorithm, a potentially large number of non-linear subproblems have to be solved to obtain a new lower bound for the problem. This can be computationally expensive, especially when dealing with a large number of connected x-variables.

• Total enumeration: The algorithm requires a total enumeration of all possible assignments of the connected x-variables to their lower and upper bounds in each iteration. This exhaustive search can be time-consuming and computationally demanding, particularly for problems with a large search space.

• Potential complexity: Depending on the structure of the given biconvex problem, the number of relaxed dual problems may still be significant even with the improvements suggested in the literature. The complexity of solving these problems can impact the overall efficiency and scalability of the algorithm.

For a detailed survey on biconvex sets and optimization, we refer the reader to the work of Gorski, Pfeuffer and Klamroth [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF].

-Joint Chance Constrained Geometric Optimization With Independent Row Vectors

In this chapter, we investigate a dynamical neural network approach to address nonconvex geometric programs with joint probabilistic constraints. These programs involve coefficients that follow a normal distribution and independent matrix row vectors. We establish the stability and convergence of our neural network based on Lyapunov analysis. Additionally, we demonstrate the equivalence between the optimal solution of the deterministic equivalent problem of geometric programs with joint probabilistic constraints and the solution of the dynamical system. Finally, we present numerical experiments to showcase the performance of our approach compared to state-of-the-art methods.

. Geometric programs

This section is published in Mathematics and Computers in Simulation Journal [START_REF] Tassouli | A neural network approach to solve geometric programs with joint probabilistic constraints[END_REF].

We consider the following general form of geometric programs with joint probabilistic constraints

min t∈R M ++ i∈I 0 c i M j=1 t a ij j , (GPJC) s.t P   i∈I k c i M j=1 t a ij j ≤ 1, k = 1, ...., K   ≥ 1 -ϵ,
where {c i } i∈I k , k ∈ {1, .., K} are pairwise independent normally distributed random variables i.e c i ∼ N (µ i , σ 2 i ), i ∈ I k where µ i ≥ 0 is the mean value and σ i is the standard deviation of c i . The coefficients a ij , i ∈ I k , j = 1, ..., M are deterministic parameters, and 1 -ϵ is a given probability level with ϵ ∈ (0, 0.5].

. Deterministic equivalent problem

Using the pairwise independence between {c i } i∈I k , the joint constraint of (GPJC) can be simplified to be equivalent to

K k=1 P   i∈I k c i M j=1 t a ij j ≤ 1   ≥ 1 -ϵ (4.1)
Consequently, problem (GPJC) can be reformulated as follows [START_REF] Liu | Stochastic geometric optimization with joint probabilistic constraints[END_REF]

min t∈R M ++ E   i∈I 0 c i M j=1 t a ij j   , s.t P   i∈I k c i M j=1 t a ij j ≤ 1   ≥ y k , k = 1, ..., K, (4.2) K k=1 y k ≥ 1 -ϵ, 0 ≤ y k ≤ 1 , k = 1, ..., K.
where y k , k = 1, ..., K, are nonnegative auxiliary variables.

A deterministic equivalent problem of (4.2) is given by [START_REF] Cheng | A second-order cone programming approach for linear programs with joint probabilistic constraints[END_REF] min

t∈R M ++ i∈I 0 µ i M j=1 t a ij j , s.t i∈I k µ i M j=1 t a ij j + ϕ -1 (y k ) i∈I k σ 2 i M j=1 t a ij j ≤ 1, k = 1, ..., K, (4.3)
K k=1 y k ≥ 1 -ϵ, 0 ≤ y k ≤ 1 , k = 1, ..., K,
where ϕ -1 (.) is the quantile of the standard normal distribution and can be represented using the inverse error function ϕ -1 (y k ) = √ 2erf -1 (2y k -1). We apply the log-transformation r j = log(t j ), j = 1, ..., M and x k = log(y k ), k = 1, ..., K, to (4.3) which leads to the following biconvex optimization problem.

min i∈I 0 µ i exp    M j=1 a ij r j    , s.t i∈I k µ i exp    M j=1 a ij r j    (4.4) + i∈I k σ 2 i exp    M j=1 (2a ij r j + log(ϕ -1 (e x k ) 2 ))    ≤ 1, k = 1, ..., K, K k=1 x k ≥ log(1 -ϵ), x k ≤ 0, k = 1, ..., K.

. Convex approximations

We introduce in this section two convex approximations proposed by Liu et al. [START_REF] Liu | Stochastic geometric optimization with joint probabilistic constraints[END_REF] to solve (4.4). Later, to evaluate the performance of our proposed neurodynamic method we compare it with those approximations.

The expression log(ϕ -1 (e x k ) 2 ) is approximated by the following piecewise linear function

F (x k ) = max s=1,...,S F s (x k ), (4.5)
where

F s (x k ) ≤ log(ϕ -1 (e x k ) 2 ), ∀x k ∈ [log(1 -ϵ), 0) and F s (x k ) = d s x k + b s , s = 1, ..., S.
One can choose the tangent lines of the expression log(ϕ -1 (e x k ) 2 ) at various points within the interval [log(1 -ϵ), 0), denoted as ζ 1 , ζ 2 , ..., ζ S and define

d s = 2e ζs (ϕ -1 ) ′ (e ζs ) ϕ -1 (e ζs ) , (4.6
)

b s = -d s ζ s + log(ϕ -1 (e ζs ) 2 ). (4.7)
We then obtain the following convex approximation of (4.4)

min i∈I 0 µ i exp    M j=1 a ij r j    , s.t i∈I k µ i exp    M j=1 a ij r j    (4.8)
+ i∈I k σ 2 i exp    M j=1 (2a ij r j + d s x k + b s )    ≤ 1, s = 1, ..., S, k = 1, ..., K, K k=1 x k ≥ log(1 -ϵ), x k ≤ 0, k = 1, ..., K.
The optimal value of problem (4.8) represents a lower bound for problem (4.4). Additionally, as the number of tangent lines S approaches infinity, (4.8) converges to an equivalent reformulation of problem (4.4) [START_REF] Liu | Stochastic geometric optimization with joint probabilistic constraints[END_REF].

To obtain an upper bound for the joint probabilistic constrained problem (GPJC), the widely used technique of sequential convex approximation was applied. The underlying principle of this approximation involves decomposing the original problem into subproblems, where a subset of variables is fixed alternately. In our specific problem, we begin by fixing y = y n and proceed to update t by solving the following subproblem

min t∈R M ++ i∈I 0 µ i M j=1 t a ij j , s.t i∈I k µ i M j=1 t a ij j + ϕ -1 (y n k ) i∈I k σ 2 i M j=1 t a ij j ≤ 1, k = 1, ..., K, (4.9)
and then fix t = t n , where t n is the optimal solution of (4.9) and update y by solving the subsequent subproblem

min y∈R K + K k=1 ψ k y k , s.t y k ≤ ϕ       1 - i∈I k µ i M j=1 (t n j ) a ij i∈I k σ 2 i M j=1 (t n j ) a ij       , k = 1, ..., K, (4.10)
K k=1 y k ≥ 1 -ϵ, 0 ≤ y k ≤ 1 , k = 1, ..., K,
where ψ k represents a selected searching direction, i.e,

ψ k = θ n k (ϕ -1 ) ′ (y n k ) i∈I k σ 2 i M j=1 (t n j ) a ij
with t n and θ n represent an optimal solution and an optimal solution of the Lagrangian dual variable θ of problem (4.9), respectively.

. A dynamical neural network approach

Let f (r) =

i∈I 0 µ i exp M j=1 a ij r j , h(r, x) = log(1 -ϵ) - K k=1 x k , g k (r, x) = x k and l k (r, x) = i∈I k µ i exp M j=1 a ij r j + i∈I k σ 2 i exp M j=1 (2a ij r j + log(ϕ -1 (e x k ) 2 )) -1.
We can write problem (4.4) as follows:

min f (r), s.t l k (r, x) ≤ 0, k = 1, .., K, h(r, x) ≤ 0, (4.11)
g k (r, x) ≤ 0, k = 1, .., K.
The feasible set of (4.11) is given by

U = (r, x) ∈ R m × R k | l k (r, x) ≤ 0, h(r, x) ≤ 0 and g k (r, x) ≤ 0, k = 1, .., K .
We define (4.11), then there exists λ (r) , λ (x) , β

U(r) = x ∈ R k | l k (r, x) ≤ 0, h(r, x) ≤ 0 and g k (r, x) ≤ 0, k = 1, .., K, and
U(x) = {r ∈ R m | l k (r, x) ≤ 0, h(r, x) ≤ 0 and g k (r, x) ≤ 0, k = 1, .., K} . Let (r * , x * ) ∈ R m × R K a partial KKT point of
i and γ (x)

i , i = 1, .., K such that (4.12)

∇f (r * ) + K i =1 β (r) i ∇ r l i (r * , x * ) + λ (r) ∇ r h(r * , x * ) + K i =1 γ (r) i ∇ r g i (r * , x * ) = 0, (4.13)
K i =1 β (x) i ∇ x l i (r * , x * ) + λ (x) ∇ x h(r * , x * ) + K i =1 γ (x) i ∇ x g i (r * , x * ) = 0, β i (r) ≥ 0, β (r)
i l i (r * , x * ) = 0, λ (r) ≥ 0 λ (r) h(r * , x * ) = 0, i = 1, .., K, (4.14) γ (r) i ≥ 0, γ (r)
i g i (r * , x * ) = 0, β (x) i ≥ 0, β (x)
i l i (r * , x * ) = 0, i = 1, .., K, (4.15) λ (x) ≥ 0, λ (x) h(r * , x * ) = 0, γ (x) i ≥ 0, γ (x)
i g i (r * , x * ) = 0, i = 1, .., K, (4.16)

Corollary 21. Let (r * , x * ) ∈ R M × R K be
(r) = λ (x) , β (r) = β (x)
and γ (r) = γ (x) , then (r * , x * ) is a KKT point of (4.11).

Let l(r, x) = (l 1 (r, x), .., l K (r, x)) T and g(r, x) = (g 1 (r, x), .., g K (r, x)) T , we propose a recurrent dynamical neural network model to solve problem (4.11). r(.), x(.), β(.), λ(.) and γ(.) are time continuous variables. The dynamical neural network is driven by the following system.

dr dt = -(∇f (r) + ∇ r l(r, x) T (β + l(r, x)) + + ∇ r h(r, x) T (λ + h(r, x)) + + ∇ r g(r, x) T (γ + g(r, x)) + ), (4.17)

dx dt = -(∇ x l(r, x) T (β + l(r, x)) + + ∇ x h(r, x) T (λ + h(r, x)) + + ∇ x g(r, x) T (γ + g(r, x)) + ), (4.18
)

dβ dt = (β + l(r, x)) + -β, (4.19
)

dλ dt = (λ + h(r, x)) + -λ, (4.20)
dγ dt = (γ + g(r, x)) + -γ. (4.21)
Let z = (r, x, β, λ, γ), we can write the dynamical system (4.17)-(4.21) as follows

dz dt = κΦ(z), z(t 0 ) = z 0 , (4.22)
where z 0 is the given initial point and κ is a scale parameter indicating the convergence rate of the neural network (4.17)-(4.21). For the sake of simplicity, we set κ = 1. Figure 4.12 illustrates the interconnections between the various inputs of the neural network (4.17)-(4.21).

In the following sections, we investigate the convergence and stability properties of the proposed dynamical neural network. Proof. Let (r * , x * , β * , λ * , γ * ) an equilibrium point of the neural network defined by (4.17)-(4.21). Then, dr * dt = 0, dx * dt = 0, dβ * dt = 0, dλ * dt = 0 and dγ * dt = 0. We have

-(∇f (r * ) + ∇ r l(r * , x * ) T (β * + l(r * , x * )) + + ∇ r h(r * , x * ) T (λ * + h(r * , x * )) + + ∇ r g(r * , x * ) T (γ * + g(r * , x * )) + ) = 0, -(∇ x l(r * , x * ) T (β * + l(r * , x * )) + + ∇ x h(r * , x * ) T (λ * + h(r * , x * )) + + ∇ x g(r * , x * ) T (γ * + g(r * , x * )) + ) = 0, (β * + l(r * , x * )) + -β * = 0, (λ * + h(r * , x * )) + -λ * = 0, (γ * + g(r * , x * )) + -γ * = 0.
Observe that (β * +l(r * , x * )) + -β * = 0 if and only if β * ≥ 0, l(r * , x * ) ≤ 0 and β * T l(r * , x * ) = 0. Similarly, we have λ * ≥ 0, h(r * , x * ) ≤ 0, λ * T h(r * , x * ) = 0, γ * ≥ 0, g(r * , x * ) ≤ 0 and γ * T g(r * , x * ) = 0. By substitution in (4.17)-(4.21), we obtain the partial KKT system (4.12)-(4.16) with λ

(r) = λ (x) = λ, β (r) = β (x) = β and γ (r) = γ (x) = γ .
The reverse implication of the theorem is straightforward.

In order to establish the stability and convergence of the neural network, our first step is to demonstrate that the Jacobian matrix ∇Φ(z) is a negative semidefinite matrix.

Lemma 23. The Jacobian matrix ∇Φ(z) is negative semidefinite matrix.

Proof. We assume the existence of 0 < p, q < K such that (β + l) + = (β 1 + l 1 (r, x), β 2 + l 2 (r, x), ....., β p + l p (r, x), 0, ...., 0

K-p
), (γ + g) + = (γ 1 + g 1 (r, x), γ 2 + g 2 (r, x), ....., γ q + g q (r, x), 0, ...., 0

K-q

).

Assuming without loss of generality that λ + h(r, x) ̸ = 0, we represent the Jacobian matrix ∇Φ as follows

∇Φ(z) =       A 1 A 2 A 3 A 4 A 5 B 1 B 2 B 3 B 4 B 5 C 1 C 2 C 3 C 4 C 5 D 1 D 2 D 3 D 4 D 5 E 1 E 2 E 3 E 4 E 5      
, where

A 1 = -(∇ 2 f (r) + p i=1 ((β i + l i )∇ 2 r l p i (r, x)) + ∇ r l p (r, x) T ∇ r l p (r, x)) + (λ + h)∇ 2 r h(r, x) + ∇ r h(r, x) T ∇ r h(r, x)) + q i=1 ((γ i + g i )∇ 2 r g q i (r, x)) + ∇ r g q (r, x) T ∇ r g q (r, x)), A 2 = -( p i=1 ((β i + l i )∇ x ∇ r l p i (r, x)) + ∇ r l p (r, x) T ∇ x l p (r, x)) + (λ + h)∇ x ∇ r h(r, x) + ∇ x h(r, x) T ∇ r h(r, x)) + q i=1 ((γ i + g i )∇ x ∇ r g q i (r, x)) + ∇ r g q (r, x) T ∇ x g q (r, x)), B 1 = -( p i=1 ((β i + l i )∇ r ∇ x l p i (r, x)) + ∇ r l p (r, x) T ∇ x l p (r, x)) + (λ + h)∇ r ∇ x h(r, x) + ∇ x h(r, x) T ∇ r h(r, x)) + q i=1 ((γ i + g i )∇ r ∇ x g q i (r, x)) + ∇ x g q (r, x) T ∇ r g q (r, x)), B 2 = -( p i=1 (β i + l i )∇ 2 x l p i (r, x) + ∇ x l p (r, x) T ∇ x l p (r, x)) + (λ + h)∇ 2 x h(r, x) + ∇ x h(r, x) T ∇ x h(r, x) + q i=1 ((γ i + g i )∇ 2 x g q i (r, x)) + ∇ x g q (r, x) T ∇ x g q (r, x)), A 3 = -C 1 = -∇ r l p (r, x) T , A 4 = -D 1 = -∇ r h(r, x) T , A 5 = -E 1 = -∇ r g q (r, x) T , B 3 = -C 2 = -∇ x l p (r, x) T , B 4 = -D 2 = -∇ x h(r, x) T , B 5 = -E 2 = -∇ x g q (r, x) T , C 3 = -S p = O p×p O p×(K-p) O (K-p)×p I (K-p)×(K-p) , E 5 = -S q = O q×q O q×(K-q) O (K-q)×q I (K-q)×(K-q)
, C 4 = 0, C 5 = 0, D 3 = 0, D 4 = 0, D 5 = 0, E 3 = 0 and E 4 = 0. Since g, h and l are twice differentiable, by Schwarz's theorem, we have ∇ r ∇ x g q i (x, y) = ∇ x ∇ r g q i (x, y), ∀i ∈ {1, .., q}, ∇ r ∇ x l p i (x, y) = ∇ x ∇ r l p i (x, y), ∀i ∈ {1, .., p} and ∇ r ∇ x h(x, y) = ∇ x ∇ r h(x, y). It follows that A 2 = B T 1 and ∇Φ(z) becomes

∇Φ(z) =       A 1 B T 1 A 3 A 4 A 5 B 1 B 2 B 3 B 4 B 5 -A 3 -B 3 -S p 0 0 -A 4 -B 4 0 0 0 -A 5 -B 5 0 0 -S q      
.

Demonstrating that -S p and -S q are negative semidefinite is straightforward. Hence, we

can conclude that S =   -S p 0 0 0 0 0 0 0 -S q   is negative semidefinite.
Given that the function f is convex and twice differentiable, we can conclude that the Hessian matrix ∇ 2 f (x) is positive semidefinite. Furthermore, since g, h, and l are biconvex and twice differentiable, it follows from [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF] that the Hessian matrices ∇ 2 x g q i (x, y), ∇ 2 y g q i (x, y), ∇ 2 x h(x, y), ∇ 2 y h(x, y), ∇ 2 x l p i (x, y), and ∇ 2 y l p i (x, y) are positive semidefinite for i = 1, . . . , q. Consequently, we can conclude that A 1 and B 2 are negative semidefinite

and then A = A 1 B T 1 B 1 B 2 is negative semidefinite. Let B = A 3 A 4 A 5 B 3 B 4 B 5 , ∇Φ(z) can be written as ∇Φ(z) = A B -B T S
. Since S and A are negative semidefinite, the conclusion follows [START_REF] Foias | Positive Definite Block Matrices[END_REF]. The proof when λ + h(r, x) = 0 follows along the same lines.

Before proving the stability of our neural network, we introduce the following definition and lemma. Definition 12. [START_REF] Rockafellar | Variational Analysis[END_REF] A mapping F : R n -→ R n is said to be monotonic if

(x -y) T (F (x) -F (y)) ≥ 0, ∀x, y ∈ R n .
Lemma 24. [START_REF] Rockafellar | Variational Analysis[END_REF]. A differentiable mapping F : R n -→ R n is monotonic, if and only if the jacobian matrix ∇F (x), ∀x ∈ R n , is positive semidefinite.

Theorem 25. The neural network (4.17 Proof. The proof of this theorem follows a similar approach to that of Theorem 3 in [START_REF] Xia | A recurrent neural network for solving linear projection equations[END_REF].

Consider y * = (r * , x * , β * , λ * , γ * ) as an equilibrium point for the neural network described by equations (4.17)-(4.21). We define the following Lyapunov function Hence, we can conclude that

V (y) = ||Φ(y)|| 2 + 1 2 ||y -y * || 2 .
dV (y(t)) dt ≤ 0.
Consequently, the neural network defined by equations (4.17)-(4.21) is stable at the equilibrium point y in the sense of Lyapunov. Since V (y) ≥ 1 2 ∥y -y * ∥ 2 , there exists a convergent subsequence (y(t k ) k≥0 ) such that lim k-→∞ y(t k ) = ŷ and dV (ŷ(t)) dt = 0. Let M = {y(t)| dV (y(t)) dt = 0}, using LaSalle's invariance principle [START_REF] Slotine | Applied nonlinear control[END_REF] we can see that any solution starting from a certain y 0 converges to the largest invariant set contained in M .

Notice that

             dr dt = 0 dx dt = 0 dβ dt = 0 dλ dt = 0 dγ dt = 0 ⇔ dV (y) dt = 0.
It follows that ŷ is an equilibrium point of the neural network (4.17)-(4.21).

Now we consider a new Lyapunov function defined by V (y) = ∥Φ(y)∥ 2 + 1 2 ∥y -ŷ∥ 2 . Since V is continuously differentiable, V (ŷ) = 0 and lim k-→∞ y(t k ) = ŷ then lim t-→∞ V (y(t)) = V (ŷ) = 0. Furthermore, since 1 2 ∥y -ŷ∥ 2 ≤ V (y(t)) then lim t-→∞ ∥y -ŷ∥ = 0 and lim t-→∞ y(t) = ŷ. We conclude that the neural network (4.17

. Numerical experiments

To assess the performance of our approach, we investigate a transportation problem that involves optimizing the shape of a transportation box while adhering to geometric constraints. Initially, we focus on a three-dimensional scenario and subsequently extend it to larger instances to showcase the effectiveness of our approach. For implementation, we utilize Python to develop our dynamical neural network. We employ numpy.random to randomly generate problem instances and use solve_ivp from the scipy.integrate package to solve the system of ordinary differential equations (ODEs). The deterministic equivalent programs are solved using the gekko package, while autograd.grad and autograd.jacobian are employed to compute gradients and partial derivatives.

To evaluate the performance, we compare our neural network approach with a state-ofthe-art convex approximation method proposed by Liu et al. [START_REF] Liu | Stochastic geometric optimization with joint probabilistic constraints[END_REF] presented in Section 4.1.2. Our assessment primarily focuses on the quality of the solution, and we do not record the CPU time as the current ODE solvers tend to be time-consuming. All experiments are conducted on an Intel(R) Core(TM) i7-10610U CPU @ 1.80 GHz.

. A three-dimension shape optimization problem

We first consider a transportation problem involving the transfer of grain from a warehouse to a factory. The grain is transported using an open rectangular box with dimensions of length x 1 meters, width x 2 meters, and height x 3 meters, as illustrated in Figure 4.2. The primary objective is to maximize the volume of the rectangular box, which is given by x 1 x 2 x 3 . There are two additional constraints imposed on the floor area A floor and the wall area A wall of the box to ensure it conforms to the specified shape of a given truck. To account for the stochastic nature of the problem, we assume that the wall area A wall and the floor area A floor are independent random variables. The minimization problem Figure 4.2: 3D-box shape [1] can be formulated as follows

min x∈R 3 ++ x -1 1 x -1 2 x -1 3 , s.t P 1 A wall (2x 3 x 2 + 2x 1 x 3 ) ≤ 1, 1 A f loor x 1 x 2 ≤ 1 ≥ 1 -ϵ. αx -1 1 x 2 ≤ 1, (4.24
)

βx -1 2 x 2 ≤ 1, γx -1 2 x 3 ≤ 1, δx -1 3 x 2 ≤ 1,
The deterministic equivalent problem of (4.24) can be written as

minimize x∈R 3 ++ x -1 1 x -1 2 x -1 3 , s.t P 1 A wall (2x 3 x 2 + 2x 1 x 3 ) ≤ 1 ≥ y 1 , P 1 A f loor x 1 x 2 ≤ 1 ≥ y 2 , αx -1 1 x 2 ≤ 1, βx -1 2 x 2 ≤ 1, γx -1 2 x 3 ≤ 1, (4.25
) In addition to solving the stochastic problem (4.24), we solve also the corresponding deterministic problem where the random variables 1

δx -1 3 x 2 ≤ 1, y 1 y 2 ≥ 1 -ϵ, 0 ≤ y 1 , y 2 ≤ 1.
A wall and 1

A floor are replaced by their means. We compare the results of the individual and joint probabilistic constraints formulations.

To evaluate the robustness of the deterministic formulation compared to the individual and joint probabilistic formulations, we have randomly generated a set of 100 instances of the stochastic variables 1

A wall and 1

A floor according to the mentioned means and variances.

We have checked whether the solutions of the three formulations are feasible for problem (6.33) for each of the 100 scenarios. If the solutions are not feasible for a scenario, we label it as a violated scenario (VS).

The numerical results are summarized in Table 4.1. The first column represents the objective value obtained by the deterministic approach. The second column shows the number of violated scenarios out of the 100 scenarios. Columns three to six present the numerical results for the individual and joint chance-constrained formulations, respectively. The table indicates that our neural network successfully solves both the individual and joint chance-constrained 3D-transportation problems. It is worth noting that the joint formulation exhibits the highest level of robustness, as none of the 100 scenarios are violated, whereas the deterministic approach violates all 100 scenarios, as depicted in Figure 4.3.

In order to assess the performance of our approach, we compare our solution with the sequential convex approximation algorithm presented in Section 4.1.2. The numerical results are provided in Table 4.2. The first and second columns represent the objective value and the number of violated scenarios obtained using our neural network, respectively. The third and fourth columns display the objective value and the number of violated scenarios obtained using the sequential algorithm, respectively. The fifth column presents the gap between the two objective values, denoted as GAP, which is calculated Table 4.2 demonstrates that the upper bounds obtained by our neural network and the sequential algorithm are identical. However, the sequential algorithm violates 4 scenarios, while our neural network solution remains feasible for all scenarios, as depicted in

. Multidimensional transportation problem

In order to evaluate the performance of our approach on large-size instances, we extend the transportation problem studied in the previous subsection to higher-dimensional problems. The problem formulation for the m-dimensional case can be written as follows.

min x∈R M ++ m i=1 x -1 i , s.t P( m-1 j=1 ( m -1 A wallj x 1 m i=1,i̸ =j x i ) ≤ 1, 1 A f loor m j=2 x j ≤ 1, (4.26) 1 γ i,j x i x -1 j ≤ 1, 1 ≤ i ̸ = j ≤ m) ≥ 1 -ϵ.
For the numerical experiments, we set the following parameters: ϵ = 0.05, 1

A f loor follows a normal distribution with a mean of 1.0/20.0 and a standard deviation of 0.1, 1

A wall j
is drawn from a normal distribution with a mean randomly selected from the interval [1.0/60.0, 1.0/40.0] and a standard deviation randomly selected from the interval [0.3, 0.5],

and 1 γ i,j follows a normal distribution with a mean of 0.5 and a standard deviation of 0.1. The results for different values of m for problem (4.26) are given in Table 4. [START_REF] Adasme | A stochastic geometric programming approach for power allocation in wireless networks[END_REF] where column one gives the data of the problem, e.g., [START_REF]Machine Learning and Wireless Communications[END_REF][START_REF] Adasme | A stochastic geometric programming approach for power allocation in wireless networks[END_REF] means the problem is composed of two variables and four constraints. Columns two and three show the optimal value obtained by the deterministic approach and the number of VS over 100 generated scenarios, respectively. Columns four to six give the objective value of the individual chance constrained problem, the number of VS and the gap with the deterministic problem, respectively. We note that the gap is computed compared to the deterministic approach. The last three columns give the optimal value obtained by the joint chance-constrained problem, the number of VS and the gap with the deterministic approach, respectively. Our joint chance-constrained neural network demonstrates greater robustness compared to the individual chance constraint and the deterministic models. On average, our joint chance constraint neural network has 10 violated scenarios (VS), while the average numbers of VS for the individual chance constraint and the deterministic approach are 40 and 92, respectively. Figure 4.5 illustrates the distributions of VS for the three formulations. Please note that we have not recorded the CPU time as our neural network problems are solved using standard ODE solvers, which are currently more time-consuming compared to standard optimization methods.

In a similar way as for the previous subsection, we compare the results obtained by the neural network for solving (4.26) with those obtained using the sequential algorithm. The results are presented in Table 4.4, where column one indicates the number of variables, column two indicates the number of constraints, column three shows the objective value obtained using the neural network, column four shows the corresponding number of violated scenarios (VS), column five shows the objective value obtained by the sequential

. Rectangular geometric programs

In this Section, we study the following stochastic rectangular programming problem

min t∈R M ++ E   i∈I 0 c i M j=1 t a ij j   , (4.27)
s.t. P   α k ≤ i∈I k c i M j=1 t a ij j ≤ β k , k = 1, ..., K   ≥ 1 -ϵ, (4.28)
where t = (t j ) 1≤j≤M is the decision variable, c i , i ∈ I k , k = 0, 1, ..., K are uncorrelated normally distributed random variables , i.e., c i ∼ N (c i , σ 2 i ), ci ≥ 0. The coefficients a ij , i ∈ I k , j = 1, ..., M are deterministic, 0 < α k < β k , and 1 -ϵ is a given probability level. Liu et al. [START_REF] Lisser | Rectangular chance constrained geometric optimization[END_REF] propose convex approximations based on the variable transformation to solve problem (4.27)-(4.28) with an elliptical distribution. They give upper and lower bounds for the optimal solution.

. Deterministic biconvex equivalent problem

Problem (4.27)-(4.28) is a joint constrained program. To transform the joint constraints into deterministic ones, we assume that the row vector constraints are mutually independent. Then, we introduce auxiliary variables y k , k = 1, ..., K and we rewrite the constraint (4.28) equivalently as

P   α k ≤ i∈I k c i M j=1 t a ij j ≤ β k   ≥ y k , k = 1, ..., K, (4.29) K k=1 y k ≥ 1 -ϵ, 0 ≤ y k ≤ 1, k = 1, ..., K. (4.30)
The rectangular constraints (6.76) are equivalent to

P   i∈I k c i M j=1 t a ij j ≥ α k   + P   i∈I k c i M j=1 t a ij j ≤ β k   -1 ≥ y k , k = 1, ..., K. (4.31)
Then, we introduce two additional K-dimensional auxiliary variables z + , z -∈ R K + [START_REF] Lisser | Rectangular chance constrained geometric optimization[END_REF] such that (4.31) is equivalent to

P   i∈I k c i M j=1 t a ij j ≥ α k   ≥ z + k , k = 1, ..., K, (4.32)
P   i∈I k c i M j=1 t a ij j ≤ β k   ≥ z - k , k = 1, ..., K, (4.33)
z + k + z - k -1 ≥ y k , 0 ≤ z - k ; z + k ≤ 1, k = 1, ..., K. (4.34)
Deterministic reformulations of constraints (4.32) and (4.33) are given as follows [START_REF] Cheng | A second-order cone programming approach for linear programs with joint probabilistic constraints[END_REF] -

i∈I k ci M j=1 t a ij j + ϕ -1 (z + k ) i∈I k σ 2 i M j=1 t 2a ij j ≤ -α k , k = 1, ..., K, (4.35)
i∈I k ci M j=1 t a ij j + ϕ -1 (z - k ) i∈I k σ 2 i M j=1 t 2a ij j ≤ β k , k = 1, ..., K. (4.36)
A biconvex equivalent formulation of constraint (4.79) can be obtained using a simple logarithmic transformation. Constraint (4.35) can be reformulated as follows

ϕ -1 (z + k ) 2 i∈I k σ 2 i M j=1 t 2a ij j ≤   i∈I k ci M j=1 t a ij j -α k   2 , k = 1, ..., K. (4.37)
We write (4.37) equivalently for a given k as

2α k i∈I k ci M j=1 t a ij j - i∈I k p∈I k ci cp M j=1 t (a ij +a pj ) j + ϕ -1 (z + k ) 2 i∈I k σ 2 i M j=1 t 2a ij j ≤ α 2 k ,
which can be reformulated as

2α k i∈I k ci M j=1 t a ij j + ϕ -1 (z + k ) 2 i∈I k σ 2 i M j=1 t 2a ij j i∈I k p∈I k ci cp M j=1 t (a ij +a pj ) j + α 2 k ≤ 1. (4.38)
We approximate the denominator of constraint (4.38) with a monomial function by applying the arithmetic-geometric mean inequality [START_REF] Adasme | A stochastic geometric programming approach for power allocation in wireless networks[END_REF].

i∈I k p∈I k ci cp M j=1 t (a ij +a pj ) j + α 2 k ≥ i∈I k p∈I k      ci cp M j=1 t (a ij +a pj ) j δ ip      δ ip α 2 k δ 0 δ 0 , (4.39)
where δ 0 and δ ip are nonnegative parameter ∀i ∈ I k and δ 0 + i,p∈I k

δ ip = 1.
Then, we approximate the denominator of constraint (4.38) with

i∈I k p∈I k    ci cp M j=1 t (a ij +a pj ) j δ ip    δ ip α 2 k δ 0 δ 0
. We write consequently problem (2.1)-(2.36) as

min t∈R M ++ ,y,z + ,z - i∈I 0 ci M j=1 t a ij j , (4.40)
s.t.   2α k i∈I k ci M j=1 t a ij j + ϕ -1 (z + k ) 2 i∈I k σ 2 i M j=1 t 2a ij j   × i∈I k p∈I k      ci cp M j=1 t (a ij +a pj ) j δ ip      -δ ip α 2 k δ 0 -δ 0 ≤ 1, k = 1, ..., K, (4.41)
i∈I k ci M j=1 t a ij j + ϕ -1 (z - k ) i∈I k σ 2 i M j=1 t 2a ij j -β k ≤ 0, k = 1, ..., K, z + k + z - k -1 ≥ y k , 0 ≤ z - k ; z + k ≤ 1, k = 1, ..., K, (4.42)
K k=1 y k ≥ 1 -ϵ, 0 ≤ y k ≤ 1, k = 1, ..., K. (4.43)
We apply a logarithmic transformation of the problem by introducing r j = log(t j ), j = 1, ..., M , x k = log(y k ), k = 1, ..., K. We obtain the following biconvex equivalent problem

min r,x,z + ,z - i∈I 0 ci exp    M j=1 a ij r j    , (4.44)
s.t. i∈I k p∈I k ci cp δ ip -δ ip α 2 k δ 0 -δ 0 ×   2α k i∈I k ci exp    M j=1 a ij r j    + ϕ -1 (z + k ) 2 i∈I k σ 2 i exp    M j=1 2a ij r j      exp    i∈I k p∈I k -δ ip M j=1 (a ij + a pj )r j    ≤ 1, k = 1, ..., K, (4.45
)

i∈I k ci exp    M j=1 a ij r j    + ϕ -1 (z - k ) i∈I k σ 2 i exp    M j=1 2a ij r j    -β k ≤ 0, k = 1, ..., K, (4.46)
log(1 -ϵ) -

K k=1 x k ≤ 0, x k ≤ 0, k = 1, ..., K, (4.47) exp(x k ) -z + k -z - k + 1 ≤ 0, k = 1, ..., K, (4.48)
z + k -1 ≤ 0, k = 1, ..., K, (4.49)
z - k -1 ≤ 0, k = 1, ..., K, (4.50)
-z + k ≤ 0, k = 1, ..., K, (4.51)
-z - k ≤ 0, k = 1, ..., K.

(4.52) Let z = (z + , z -) T , f (r) = i∈I 0 ci exp M j=1 a ij r j , h(x) = log(1 -ϵ) - K k=1 x k , x 1 , . . . , x K T , l(x, z) = exp(x 1 ) -z + 1 -z - 1 + 1, ..., exp(x K ) -z + K -z - K + 1 T , w(z) = (z 1 + -1, . . . , z K + -1, z 1 --1, . . . , z K --1, -z 1 + , . . . , -z K + , -z 1 -, . . . , -z K -) T and g(r, z) =                                   i∈I 1 p∈I 1 ci cp δ ip -δ ip α 2 1 δ 0 -δ 0 (2α 1 i∈I 1 ci exp M j=1 a ij r j + ϕ -1 (z + 1 ) 2 × i∈I 1 σ 2 i exp M j=1 2a ij r j )exp i∈I 1 p∈I 1 -δ ip M j=1 (a ij + a pj )r j -1 . . . i∈I K p∈I K ci cp δ ip -δ ip α 2 K δ 0 -δ 0 (2α K i∈I K ci exp M j=1 a ij r j + ϕ -1 (z + K ) 2 × i∈I K σ 2 i exp M j=1 2a ij r j )exp i∈I K p∈I K -δ ip M j=1 (a ij + a pj )r j -1 i∈I 1 ci exp M j=1 a ij r j + ϕ -1 (z - 1 ) i∈I 1 σ 2 i exp M j=1 2a ij r j -β 1 . . . i∈I K ci exp M j=1 a ij r j + ϕ -1 (z - K ) i∈I K σ 2 i exp M j=1 2a ij r j -β K                                   ,

. Optimality conditions

Now, we study the optimality conditions for problem (4.53). Since the problem is biconvex, we do not talk about a KKT system but rather a partial KKT system [START_REF] Jiang | Partial exactness for the penalty function of biconvex programming[END_REF].

Definition 13. Let (r * , z * , x * ) ∈ R M × R 2K × R K , if
there exist µ (1) , µ (2) , λ (1) , λ (2) , γ and ω such that

∇f (r * ) + µ (1) T ∇ r g(r * , z * ) = 0, (4.54) µ (1) ≥ 0, µ (1) T g(r * , z * ) = 0, (4.55) µ (2) T ∇ z g(r * , z * ) + λ (1) T ∇ z l(x * , z * ) + γ T ∇ z w(z * ) = 0, (4.56) µ (2) ≥ 0, µ (2) T g(r * , z * ) = 0, λ (1) ≥ 0, λ (1) T l(x * , z * ), γ ≥ 0, γ T w(z * ) = 0, (4.57) λ (2) T ∇ x l(x * , z * ) + ω T ∇ x h(x * ) = 0, (4.58) λ (2) ≥ 0, λ (2) T l(x * , z * ) = 0, ω ≥ 0, ω T h(x * ) = 0 (4.59)
Then (r * , z * , x * ) is a partial KKT point of (4.53).

The following theorem is driven by the equivalence between a partial optimum and a partial KKT point of a biconvex program. Theorem 26. Let (r * , z * , x * ) ∈ R M × R 2K × R K be a partial optimum of (4.53), with respect to partial Slater constraints qualification [START_REF] Jiang | Partial exactness for the penalty function of biconvex programming[END_REF] at (r * , x * ). Then (r * , z * , x * ) is a KKT point of (4.53) if and only if the partial KKT system (4.54)-(4.59) holds with µ (1) = µ (2) and λ (1) = λ (2) . Furthermore, if µ (1) = µ (2) and λ (1) = λ (2) then (r * , z * , x * ) is a KKT point of (4.53).

. A dynamical neural network approach

Based on the partial KKT system (4.54)-(4.59) obtained in the previous section, we construct a dynamical neural network system that converges to a partial KKT point of (4.53). The dynamical neural network is driven by the following system, where r, z, x, µ, λ, γ, and ω are time-dependent variables

dr dt = -(∇f (r) + ∇ r g(r, z) T (µ + g(r, z)) + ), (4.60
)

dz dt = -(∇ z g(r, z) T (µ + g(r, z)) + ∇ z l(x, z) T (λ + l(x, z)) + + ∇ z w(z) T (γ + w(z)) + ), (4.61)
dx dt = -(∇ x l(x, z) T (λ + l(x, z)) + ∇ x h(x) T (ω + h(x)) + ), (4.62)
dµ dt = (µ + g(r, z)) + -µ, (4.63
)

dλ dt = (λ + l(x, z)) + -λ, (4.64
)

dγ dt = (γ + w(z)) + -γ, (4.65
)

dω dt = (ω + h(x)) + -ω. (4.66)
where (x) + = max(x, 0). For the sake of simplicity, let y = (r, z, x, µ, λ, γ, ω) we rewrite the dynamical system (4.60)-(4.66) equivalently as follows Proof. (x + y) + = x means that x ≥ 0. Furthermore, (x + y) + = x implies that (x T x + x T y) + = x T x which means that x T y = 0 since x T x ̸ = 0. Remember that x and y are n-dimensional vectors, then (x + y) + = x means that ∀ 1 ≤ i ≤ n (x i + y i ) + = x i . If x i = 0 then (y i ) + = 0 which means y i ≤ 0, else if x i > 0 then y i = 0. In conclusion y ≤ 0. The conclusion follows.

Theorem 28. Let y = (r, z, x, µ, λ, γ, ω) ∈ R M × R 2K × R K × R 2K × R K+1 × R K × R 4K , y
is an equilibrium point of (4.60)-(4.66) if and only if (r, z, x) is a KKT point of (4.53).

Proof. Let (r, z, x, µ, λ, γ, ω) an equilibrium point of (4.60)-(4.66), then dr dt = 0, dz dt = 0, dx dt = 0, dµ dt = 0, dλ dt = 0, dγ dt = 0 and dω dt = 0.

By Lemma 27, dµ dt = 0 ⇔ (µ + g(r, z)) + = µ ⇔ µ ≥ 0 and g(r, z) ≤ 0 and µ T g(r, z) = 0 ⇔ (4.56) . We use the same approach to obtain (4.57) and (4.59).

dr dt = 0 ⇔ ∇f (r) + ∇ r g(r, z) T (µ + g(r, z)) + = 0 ⇔ ∇f (r) + ∇ r g(r, z) T µ = 0 ⇔ (4.54) . We obtain (4.56) and (4.58) following the same steps. We conclude that (r, z, x) is a partial KKT system of (4.53). It is easy to check the converse part of the theorem. Now, to show the stability and the convergence of our neural network, we need first to prove the negative semidefiniteness of the Jacobian matrix ∇Φ(y).

Theorem 29. The Jacobian matrix ∇Φ(y) is negative semidefinite.

Proof. Let p, q, s, t ∈ N such that (µ + g) + = (µ 1 + g 1 (r, z), µ 2 + g 2 (r, z), ....., µ p + g p (r, z), 0, ...., 0

2K-p

), (λ + l) + = (λ 1 + l 1 (x, z), λ 2 + l 2 (x, z), ....., λ q + l q (x, z), 0, ...., 0

K-q
), (γ + w) + = (γ 1 + w 1 (z), γ 2 + w 2 (z), ....., γ s + w s (z), 0, ...., 0

4K-s

),

(ω + h) + = (ω 1 + h 1 (x), ω 2 + h 2 (x), ....., ω t + h t (x), 0, ...., 0

K+1-t

).

We write

∇Φ(z) =           A 1 A 2 A 3 A 4 A 5 A 6 A 7 B 1 B 2 B 3 B 4 B 5 B 6 B 7 C 1 C 2 C 3 C 4 C 5 C 6 C 7 D 1 D 2 D 3 D 4 D 5 D 6 D 7 E 1 E 2 E 3 E 4 E 5 E 6 E 7 F 1 F 2 F 3 F 4 F 5 F 6 F 7 G 1 G 2 G 3 G 4 G 5 G 6 G 7          
, where,

A 1 = -(∇ 2 f (r) + p i=1 ((µ i + g i )∇ 2 r g p i (r, z)) + ∇ r g p (r, z) T ∇ r g p (r, z)), A 2 = -( p i=1 ((µ i + g i )∇ z ∇ r g p i (r, z)) + ∇ z g p (r, z) T ∇ r g p (r, z)), A 4 = -∇ r g p (r, z) T , A 3 = 0, A 5 = 0, A 6 = 0, A 7 = 0, B 1 = -( p i=1 ((µ i + g i )∇ r ∇ z l p i (r, z)) + ∇ z g p (r, z) T ∇ r g p (r, z)), B 2 = -( p i=1 ((µ i + g i )∇ 2 z g p i (r, z)) + ∇ r g p (r, z) T ∇ z g p (r, z) + q i=1 ((λ i + l i )∇ 2 z l q i (x, z)) + ∇ z l q (x, z) T ∇ z l q (x, z) + s i=1 ((γ i + w i )∇ 2 z w s i (z)) + ∇ z w s (z) T ∇ z w s (z)), B 3 = -( q i=1 ((λ i + l i )∇ x ∇ z l q i (x, z)) + ∇ z l q (x, z) T ∇ x l q (x, z)), B 4 = -∇ z g p (r, z) T , B 5 = -∇ z l q (x, z) T , B 6 = -∇ z w s (z) T , B 7 = 0, C 1 = 0, C 2 = -( q i=1 ((λ i + l i )∇ z ∇ x l q i (x, z)) + ∇ x l q (x, z) T ∇ z l q (x, z)), C 3 = -( q i=1 ((λ i + l i )∇ 2 x l q i (x, z)) + ∇ x l q (x, z) T ∇ x l q (x, z) + t i=1 ((ω i + h i )∇ 2 x h t i (x)) + ∇ x h t (x) T ∇ x h t (x)), C 4 = 0, C 6 = 0, C 5 = -∇ x l q (x, z) T , C 7 = -∇ x h t (x) T , D 1 = ∇ r g p (r, z) T , D 2 = ∇ z g p (r, z) T , D 3 = 0, D 4 = S p = - O p×p O p×(2K-p) O (2K-p)×p I (2K-p)×(2K-p)
, D 5 = 0, D 6 = 0,D 7 = 0,

E 1 = 0, E 2 = ∇ z l q (x, z) T , E 3 = ∇ x l q (x, z) T , E 4 = 0, E 5 = S q = - O q×q O q×(K-q) O (K-q)×q I (K-q)×(K-q) , , E 6 = 0, E 7 = 0, F 1 = 0, F 2 = ∇ z w s (z) T , F 3 = 0, F 4 = 0, F 5 = 0, F 6 = S s = - O s×s O s×(4K-s) O (4K-s)×s I (4K-s)×(4K-s) , , F 7 = 0, G 1 = 0, G 2 = 0, G 3 = ∇ x h t (x) T , G 4 = 0, G 5 = 0, G 6 = 0, G 7 = S t = - O t×t O t×(K+1-t) O (K+1-t)×t I (K+1-t)×(K+1-t) , .
We rewrite the Jacobian matrix ∇Φ as follows,

∇Φ(z) =           A 1 B T 1 0 A 4 0 0 0 B 1 B 2 B 3 B 4 B 5 B 6 0 0 B T 3 C 3 0 C 5 0 C 7 -A 4 -B 4 0 S p 0 0 0 0 -B 5 -C 5 0 S q 0 0 0 -B 6 0 0 0 S s 0 0 0 -C 7 0 0 0 S t           , We can represent ∇Φ as ∇Φ(z) = A B -B T C , where A =   A 1 B T 1 0 B 1 B 2 B 3 0 B T 3 C 3   , B =   A 4 0 0 0 B 4 B 5 B 6 0 0 C 5 0 C 7   , and C =     S p 0 0 0 0 S q 0 0 0 0 S s 0 0 0 0 S t     ,
Since w and h are convex and twice differentiable, there follows that ∇ 2 z w s i (z) and ∇ 2 x h t i (x) are positive semidefinite. Furthermore, g and l are biconvex and twice differentiable, then we have ∇ 2 z g p i (r, z), ∇ 2 z l q i (x, z) ∇ 2 x l q i (x, z) are positive semidefinite [START_REF] Gorski | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF]. It is easy to see that ∇ r g p (r, z) T ∇ z g p (r, z), ∇ z l q (x, z) T ∇ z l q (x, z) and ∇ x l q (x, z) T ∇ x l q (x, z) are positive semidefinite. We conclude that B 2 and C 3 are negative semidefinite and hence

B 2 B 3 B T 3 C 3
is negative semidefinite [START_REF] Foias | Positive Definite Block Matrices[END_REF]. Following the same steps, we show that A 1 is negative semidefinite, and we conclude that A is negative semidefinite. We can easily see that C is negative semidefinite. We conclude that ∇Φ is negative semidefinite.

Theorem 30. The neural network (4.60)-(4.66) is stable and converges to y * = (r * , z * , x * , µ * , λ * , γ * , ω * ) where (r * , z * , x * ) is a KKT point of (4.53).

Before giving the proof of Theorem 30, we need to introduce the relationship between the monotocity of mapping and the semidefiniteness of its Jacobian matrix. Definition 14. [START_REF] Rockafellar | Variational Analysis[END_REF] A mapping F : R n -→ R n is said to be monotonic if

(x -y) T (F (x) -F (y)) ≥ 0, ∀x, y ∈ R n
Lemma 31. [START_REF] Rockafellar | Variational Analysis[END_REF] A differentiable mapping F : R n -→ R n is monotonic, if and only if the Jacobian matrix ∇F (x), ∀x ∈ R n , is positive semidefinite.

Proof. of Theorem 30. There follows that the neural network (6.47)-(4.66) is stable in the sense of Lyapunov.

Notice that V 1 (y) ≥ 1 2 ∥y -y * ∥ 2 , consequently there exists a convergent subsequence (y(t k ) k≥0 ) such that lim k-→∞ y(t k ) = ỹ and dV 1 (ỹ) dt = 0. Starting from a certain y 0 , we have by LaSalle's invariance principle that the neural network converges to the largest invariant set contained in M which is defined by M = {y(t)∥ dV 1 (y) dt = 0}.

Observe that dy dt = 0 ⇔ dV 1 (y) dt = 0, we have then ỹ is an equilibrium point of (4.60)-(4.66). Let show now that the neural network converges to ỹ. For this, we consider the following Lyapunov function V 2 (y) = ∥Φ(y)∥ 2 + 1 2 ∥y -ỹ∥ 2 . We have that V 2 is continuously differentiable, V 2 (ỹ) = 0 and lim k-→∞ y(t k ) = ỹ, then

lim t-→∞ V 2 (y(t)) = V 2 (ỹ) = 0.
Additionally, since 1 2 ∥y -ỹ∥ 2 ≤ V 2 (y) then lim t-→∞ ∥y -ỹ∥ = 0 and lim t-→∞ y(t) = ỹ. There follows that the neural network converges to an equilibrium point ỹ = (r, z, x, μ, λ, γ, ω) where (r, z, x) is a KKT point of (4.53).

. Numerical experiments

In order to test the performances of our proposed neural network, we first consider a transportation problem. Then, in Subsection 4.2.6, we study a generalized shape optimization problem to analyze the behavior of the neural network for different size instances and we compare it to a state-of-the-art alternate convex search. All the numerical experiments are done using Python. To compute the partial derivatives and the Jacobians, we use the package autograd. To generate the random instances, we use the package numpy.random. The ODEs of the recurrent dynamical neural networks are solved using the function solve_ivp of scipy.integrate library. We run our algorithms on Intel(R) Core(TM) i7-10610U CPU @ 1.80GHz.

. Minimizing transport cost problem

In order to shift V m 3 grains from a warehouse to a factory, we can use an open rectangular box of length x 1 meters, of width x 2 meters, and of height x 3 meters see Figure 4.7. We use the parameters c1, c2, c3 and c4 for the bottom costs, for each side costs and each end costs, and for each round trip of the box costs, respectively. We aim to find the minimum cost of transporting V m 3 of grain. We use a transporter to carry the box into the truck. The floor area of the box x 1 x 2 must be less than β f loor A f loor and larger than α f loor A f loor in order to avoid wasting in capacity, where A f loor is the floor area and α f loor and β f loor are the minimum and the maximum occupancy rates, respectively. The same thing applies for the wall area 2x 1 x 3 + 2x 2 x 3 which must be less than β wall A wall and larger then α wall A wall . We assume that the floor and the wall areas of the transporter are random.

We reformulate then our minimization problem as

min x∈R 3 ++ c 1 x 1 x 2 + 2c 2 x 1 x 3 + 2c 3 x 2 x 3 + c 4 V x 1 x 2 x 3 , s.t. P(α wall A wall ≤ 2x 1 x 3 + 2x 2 x 3 ≤ β wall A wall , (4.67)
α f loor A f loor ≤ x 1 x 2 ≤ β f loor A f loor ) ≥ 1 -ϵ.
To solve problem (4.67) using our proposed neural network, we set α wall = α f loor = 50%,

β wall = β f loor = 95%, c 1 = 80, c 2 = 20, c 3 = 30, c 4 = 1, V = 80m 3 , 1
A wall ∼ N (1.0/6.0, 0.01) and 1

A f loor ∼ N (3.0, 0.01). The neural network converges to a minimum of 260.81 at x 1 = 0.68m, x 2 = 0.46m and x 3 = 2.01m. We follow the convergence of x 1 , x 2 and x 3 in Figure 4.8.

. Stochastic shape optimization problem

In order to evaluate the performances of the proposed dynamical network on different instances, we introduce the following shape optimization problem taken from [START_REF] Lisser | Rectangular chance constrained geometric optimization[END_REF].

We remind that A wallj and A f loor are random and defined as in the previous Subsection.

The generalized problem is defined as follows.

min x∈R M ++ m i=1 x -1 i , s.t. P(α wall ≤ m-1 j=1 ( m -1 A wallj x 1 m i=2,i̸ =j x i ) ≤ β wall , (4.68
)

α f loor ≤ 1 A f loor m j=2 x j ≤ β f loor ) ≥ 1 -ϵ.
For the sake of comparison, we additionally solve the problem with individual constraints. In order to check the robustness of our approach, we generate a set of 100 scenarios of the stochastic constraints, and we visualize the number of violated scenarios (VS) for each problem.

min x∈R M ++ m i=1 x -1 i , s.t. P   α wall ≤ m-1 j=1 ( m -1 A wallj x 1 m i=2,i̸ =j x i ) ≤ β wall   ≥ 1 -ϵ, (4.69)
P   α f loor ≤ 1 A f loor m j=2 x j ≤ β f loor   ≥ 1 -ϵ.
The numerical results are represented in Table 4 Finally, we test the impact of the confidence level 1 -ϵ. We solve then problems (4.68) and (4.69) for different value of ϵ and we fix m = 5. The obtained results are recapitulated in Table 4.6. We observe that as the value of ϵ increases, the value of the objective function decreases though the number of VS increases. In fact, a higher value of ϵ means a larger risk area and a less restrictive minimization problem.

To the best of our knowledge, the only work dealing with rectangular programs with joint probabilistic constraints is .20 [START_REF] Lisser | Rectangular chance constrained geometric optimization[END_REF]. They propose new convex approximations based on the variable transformation and piecewise linear approximation methods to come up with lower and upper bounds for the optimal solutions. Since in the worst case, our neural network converges to a partial optimum of the minimization problem, we converge then at worst to an upper bound of the optimal solution. The advantage of our approach is that we don't use any convex approximation to approximate the optimal solution though the neurodynamic approach takes more time to solve the stochastic rectangular programs as the size of the problem increases. Therefore, we must note that our approach doesn't replace the existing convex approximations but gives promising results and opens the way for a new vision of the joint probabilistic problems.

In order to evaluate the quality of the upper bound obtained using the dynamical neural network, we compare the solutions with the values given by the alternate convex search method in [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF]. We set ϵ = 0.05 and solve problem (4.68) for different values of m.

The obtained results are presented in Figure 4.11. The orange and bleu curves represent the optimal solution obtained by the dynamical neural network and the alternate convex search method, respectively. We notice that the results of our approach are comparable or even relatively better than the alternate convex search.

. A case study : Maximizing Signal to Interference Noise Ratio for Massive MIMO

The content of this Chapter was accepted in The 19 th International Conference on Mobile Web and Intelligent Information Systems.

Massive Multiple Input Multiple Output (MaMIMO) is a promising technology that is gaining attention in the field of communication systems and the Internet of Things (IoT). It involves the utilization of a large number of antennas that can simultaneously transmit and receive signals, resulting in increased capacity and improved performance. MaMIMO is considered as a potential solution for 5G networks and is expected to succeed the existing 4G LTE and LTE-A technologies. The key advantage of MaMIMO lies in its ability to handle interference among a large number of antennas. By leveraging spatial multiplexing and beamforming techniques, MaMIMO enables higher connectivity and improved network efficiency. It can adapt to high-density environments with a large number of devices, making it suitable for scenarios such as crowded stadiums, urban areas, and IoT deployments. Furthermore, MaMIMO offers reduced transmission latency, which is crucial for real-time applications like augmented reality and virtual reality. It also addresses energy efficiency concerns by optimizing the use of resources, thereby aligning with green communications guidelines. Additionally, MaMIMO provides better signal quality and enhanced security by exploiting multiple signal paths and employing advanced signal processing algorithms. We propose in this section a stochastic geometric formulation to solve a problem of maximizing signal to interference noise ratio for massive MIMO.

. Stochastic geometric formulation

We consider a single cell area, as depicted in Figure 4.12, where there exists a group of users denoted by U = 1, ..., K. Each user in this set employs a single antenna to receive data from the base station. On the other hand, the base station is equipped with T antennas. Our objective is to maximize the worst user's Signal to Interference Noise Ratio (SINR) while respecting some power constraints imposed on each user. The SINR value for user i can be mathematically expressed as follows [START_REF] Adasme | Maximizing signal to interference noise ratio for massive mimo: A mathematical programming approach[END_REF]

SINR i = p i |g H i g i | 2 j∈U ,j̸ =i p j |g H i g j | 2 +|σ i | 2 (4.70)
We formulate our optimization problem as follows

max p∈ K ++ min i∈U p i |g H i g i | 2 j∈U ,j̸ =i p j |g H i g j | 2 +|σ i | 2 , (4.71)
s.t

P min ≤ p i ≤ P max , ∀i ∈ U, (4.72)
where p i represents the power assigned to user i from the base station. The terms

g i ∈ C T ×1 , g H i ∈ C 1×T
, and σ 2 i correspond to the beam domain channel vector associated with user i, its Hermitian transpose, and the Additive White Gaussian Noise (AWGN), respectively. We assume that the AWGN follows an independent complex Gaussian distribution with zero mean and unit variance (σ i ∼ CN (0, 1)). The channel vectors g i and g H i represent the quasi-static independent and identically distributed Rayleigh fading channels, where each entry is a complex number. The variables P min and P max define the lower and upper bounds, respectively, for the power assigned to each user, ensuring that the power remains within a certain range.

Let

a ij = |g H i g j | 2 |g H i g i | -2 and b i = |σ i | 2 |g H i g i | -
a ij p j p -1 i w + b i p -1 i w ≤ 1, ∀i ∈ U, (4.74)
P min ≤ p i ≤ P max , ∀i ∈ U. (4.75)
An equivalent geometric minimization problem is given by We consider the case where the coefficients a ij and b i are not completely known and follow independent normal distributions, i.e., a ij ∼ N (µ ij , σ 2 ij ) and b i ∼ N (µ i , σ 2 i ). In this case, we replace the deterministic constraint (4.77) with the following chance constraint, as proposed in [START_REF] Adasme | A stochastic geometric programming approach for power allocation in wireless networks[END_REF].

min p∈ K ++ ,w∈ ++ w -1 , (4.76) s.t j∈U ,j̸ =i a ij p j p -1 i w + b i p -1 i w ≤ 1, ∀i ∈ U,
P    j∈U ,j̸ =i a ij p j p -1 i w + b i p -1 i w ≤ 1, ∀i ∈ U    ≥ 1 -ϵ. (4.79)
with 1 -ϵ, 0 < ϵ ≤ 0.5, is a given confidence level. The joint chance constraint, as opposed to individual constraints, is used because it ensures that the constraint as a whole is satisfied to a certain confidence level. While individual chance constraints are easier to solve, they only guarantee that each constraint is satisfied to a certain confidence level independently, without considering the overall satisfaction of the constraints.

Using the independence between the coefficients and by introducing auxiliary variables y i ∈ + , i ∈ U, we give the following deterministic equivalent for the joint constraint (4.79) [START_REF] Tassouli | A neural network approach to solve geometric programs with joint probabilistic constraints[END_REF] j∈U ,j̸ =i

µ ij p j p -1 i w + µ i p -1 i w + ϕ -1 (y i )    j∈U ,j̸ =i σ 2 ij p 2 j p -2 i w 2 + σ 2 i p -2 i w 2    ≤ 1, ∀i ∈ U, (4.80)
i∈U

y i ≥ 1 -ϵ, (4.81) 0 ≤ y i ≤ 1, ∀i ∈ U, (4.82)
We write then (6.76)-(6.78) equivalently as

min p∈ K ++ ,w∈ ++ w -1 , s.t j∈U ,j̸ =i µ ij p j p -1 i w + µ i p -1 i w + ϕ -1 (y i )    j∈U ,j̸ =i σ 2 ij p 2 j p -2 i w 2 + σ 2 i p -2 i w 2    ≤ 1, ∀i ∈ U, 1 -ϵ - i∈U y i ≤ 0, (GP) -y i ≤ 0, y i -1 ≤ 0, ∀i ∈ U, P min -p i ≤ 0, p i -P max ≤ 0, ∀i ∈ U.
The obtained equivalent deterministic problem (GP) is nonconvex, we apply then the log-

arithmic transformation r i = log(p i ), x i = log(y i ),
µ ij exp(r j -r i + t) + µ i exp(t -r i ) (4.84) +ϕ -1 (e x i )    j∈U ,j̸ =i σ 2 ij exp(2r j -2r i + 2t) + σ 2 i exp(2t -2r i )    ≤ 1, ∀i ∈ U, log(1 -ϵ) - i∈U x i ≤ 0, x i ≤ 0, i ∈ U, (4.85)
log(P min ) -r i ≤ 0, r i -log(P max ) ≤ 0, ∀i ∈ U.

(4.86)

Let z = (r, t) T , for the sake of simplicity we write the optimization problem as

min f (z), (4.87) s.t g i (z, x) ≤ 0, ∀i ∈ U, (4.88) l(x) ≤ 0, h i (x) ≤ 0, i ∈ U, (4.89) v i (z) ≤ 0, w i (z) ≤ 0, ∀i ∈ U. (4.90)
where

f (z) = exp(-t),l(x) = log(1 -ϵ) - i∈U x i , h i (x) = x i , v i (z) = log(P min ) -r i , w i (z) = r i -log(P max ) and g i (z, x) = j∈U ,j̸ =i µ ij exp(r j -r i + t) + µ i exp(t -r i ) + ϕ -1 (e x i )×    j∈U ,j̸ =i σ 2 ij exp(2r j -2r i + 2t) + σ 2 i exp(2t -2r i )    -1
The dynamical neural network is given in this case by

dz dt = -(∇ z f (z) + ∇ z g(z, x) T (α + g(z, x)) + + ∇ z v(z) T (γ + v(z)) + + ∇ z w(z) T (λ + w(z)) + ), (4.91
)

dx dt = -(∇ x l(x) T (β + l(x)) + + ∇ x g(z, x) T (α + g(z, x)) + + ∇ x h(x) T (ζ + h(x)) + ), (4.92
)

dα dt = (α + g(z, x)) + -α, (4.93
)

dγ dt = (γ + v(z)) + -γ, (4.94
)

dλ dt = (λ + w(z)) + -λ, (4.95
) where α (1) = (α

dβ dt = (β + l(x)) + -β, ( 4
(1) 1 , ..., α (1)
K ) T , α (2) = (α (2) 1 , ..., α (2)
K ) T , γ = (γ 1 , ..., γ K ) T , λ = (λ 1 , ..., λ K ) T , ζ = (ζ 1 , ..., ζ K ) T , g = (g 1 , ..., g K ) T , v = (v 1 , ..., v K ) T , w = (w 1 , ..., w K ) T and h = (h 1 , ..., h K ) T
are time-continuous vectors.

. Numerical analysis

For the numerical experiments, we set P min = 0.1 and P max = 0.5. The confidence level is ϵ = 0.1. We generate complex vectors g i ∈ C T ×1 and g H i ∈ C 1×T for each i ∈ U from an independent complex Gaussian distribution with zero mean and unit variance. Each of these vectors is then multiplied by a factor randomly chosen from the set 3.0, 4.0, 5.0, 7.0. We also generate the parameter σ i for each i ∈ U from an independent complex Gaussian distribution with zero mean and unit variance. The variables a ij and b i are computed as explained in Section 4.3.1. We assume that µ ij = a ij and µ i = b i , and we vary the values of σ ij and σ i within the set 0.1, 0.2, 0.3.

We first solve (GP) for K = 5 with different feasible initial points y 0 and observe the convergence process of the neural network for each case. As shown in Figure 4.13, we observe that the neural network converges to the same final value regardless of the starting point.

In order to demonstrate the advantage of using joint constraints instead of individual constraints to handle the uncertainty in constraints (4.77), we solve (GP) for different numbers of users, ranging from K = 2 to K = 20, using both joint and individual chance constraints. We generate 100 instances of the stochastic variables a ij and b i and track the number of times the constraints (4.77) are violated. The results are summarized in Table 4.7.

Column one represents the number of users K, while columns two and three present the optimal solution and the number of VS obtained using individual constraints, respectively. Similarly, columns four and five display the optimal solution and the number of VS obtained using joint constraints.

We observe that the number of VS when using individual constraints is consistently bigger than the number of VS when using joint constraints. This difference becomes more significant as the value of K increases. By employing joint chance constraints, we ensure a better coverage of the risk area and achieve a more robust solution. For comparison purposes, we solve problem (GP) using both the neurodynamic approach and the sequential algorithm. The results obtained are summarized in Table 6.1.

Column one represents the number of users K, while columns two and three provide the optimal solution and the number of VS obtained using the sequential algorithm, respectively. Similarly, columns four and five present the optimal solution and the number of VS obtained using the neurodynamic approach. Finally, column six displays the gap between the two solutions. We observe that the neurodynamic approach yields better solutions compared to the sequential algorithm. Additionally, the number of violated scenarios for the solutions obtained using the neurodynamic approach is slightly lower than those obtained using the sequential algorithm.

We consider the case where K = 5 and vary the value of ϵ in the range of [0.05, 0.4].

The results obtained are summarized in Table 4.9.

We observe that as ϵ increases, the problem becomes less conservative. This means that the chance constraint becomes less strict, allowing for a higher probability of violating the constraint. Additionally, we observe that the gap between the two approaches (neurodynamic and sequential) increases as ϵ increases, as shown in Figure 4.14. This indicates that the neurodynamic approach provides better solutions compared to the sequential algorithm when dealing with larger values of ϵ. Furthermore, we notice that the number of violated scenarios increases as ϵ increases, as shown in Figure 4.15. The difference in the number of violated scenarios becomes more significant as ϵ increases, indicating that the neurodynamic approach ensures better robustness in handling uncertainties.

-Dependent Joint Chance Constrained Geometric Optimization

In this chapter, we delve into the application of copula theory for modeling dependencies in joint probabilistic constrained geometric programs with dependent rows. We also study in the second Section linear programs as a particular case. We focus on scenarios where the random variables are assumed to follow an elliptical distribution. To lay the foundation, we introduce the fundamental concepts of copula theory and refresh our understanding of the key properties associated with elliptical distributions. This groundwork will provide the necessary background for our subsequent analysis.

. Copula Theory

The use of copula theory in statistical modeling has gained significant interest in recent years. Copulas provide a means to capture the dependencies between random variables and derive joint distributions. The concept of copulas was introduced by Fisher in 1959 [START_REF] Sklar | Fonctions de répartition à n dimensions et leurs marges[END_REF], and since then, it has found applications in various fields. One of the key advantages of copulas is that they provide scale-free measures of statistical dependence. This means that the measure of dependence is independent of the marginal distributions of the variables involved. This property makes copulas useful in capturing and analyzing the dependence structure between variables. The application of copulas extends to fields such as finance and economics. Jouini and Clemen [START_REF] Jouini | Copula models for aggregating expert opinions[END_REF] discuss the use of copulas in aggregating information from different sources by considering multivariate distributions that are functions of their marginals. Patton [START_REF] Patton | Copula-based models for financial time series[END_REF] provides a comprehensive review of the applications of copulas in finance and economics, highlighting their usefulness in modeling dependencies and risks. In the context of joint chance-constrained optimization, Houda and Lisser [START_REF] Houda | On the use of copulas in joint chance-constrained programming[END_REF] use copulas to model the dependence between random variables. Copulas offer a flexible framework to capture the joint distribution of uncertain variables and enable the formulation of joint chance constraints. Furthermore, Liu et al. [START_REF] Liu | A collective neurodynamic approach to distributed constrained optimization[END_REF] propose a collective neurodynamic approach using multiple interconnected recurrent neural networks for distributed constrained optimization. This approach leverages the power of neural networks to tackle optimization problems with complex constraints.

In this section, we introduce the notion of a copula and the correspondent properties we use in this paper.

∀x ∈ d F(x) = C(F 1 (x 1 ), F 2 (x 2 ), .., F d (x d )).
If F 1 , F 2 , .., F d are all continuous then C is unique and given by

C(u) = F(F 1 -1 (u 1 ), F 2 -1 (u 2 ), .., F d -1 (u d )),
otherwise C is uniquely determined on rang(F 1 ) × rang(F 2 ), .., rang(F d ).

Definition 16. If a copula C has a density, then it is expressed as follows

C(u 1 , u 2 , .., u d ) = ∂ d C(u 1 , u 2 , .., u d ) ∂u 1 u 2 ..u d
To capture the dependencies between stochastic parameters, we often rely on copulas. While many copulas lack explicit analytical expressions, such as the Gaussian copula, there exists a class of copulas known as Archimedean copulas that offer explicit mathematical formulations. Archimedean copulas derive their name from the Greek mathematician Archimedes and are defined based on specific generator functions. These generator functions characterize the underlying dependence structure of the copula.

Definition 17. A copula C is Archimedean if it is represented as follows C(u 1 , u 2 , .., u d ; θ) = ψ [-1] (ψ(u 1 ; θ), ψ(u 2 ; θ), .., ψ(u d ; θ); θ),
where ψ : [0, 1] × Θ -→ [0, ∞) a strictly continuously declining function called generator of C, such that ψ(1; θ) = 0. ψ [-1] is a pseudo-inverse of ψ defined by ψ

[-1] (t; θ) = ψ -1 (t; θ) if 0 ≤ t ≤ ψ(0; θ) 0 if t ≥ ψ(0; θ)
and θ is a dependency parameter.

(-1) k d k dt k ψ -1 (t) ≥ 0, k = 0, 1, ..K. (5.1)
Then ψ is a strict copula generator.

Theorem 34. All copula generators are convex.

We give in Table 5.1 some examples of Archimedean copulas often used in the literature.

Name of copula

Parameter θ

Generator ψ θ (t)

Frank θ > 0 -ln( e -θt-1 e -θ -1 )
Gumbel-Hougaard

θ ≥ 1 e -t 1 θ Joe θ > 1 -ln(1 -(1 -t) θ )
Table 5.1: Examples of some commonly used Archimedean copulas

. Elliptically symmetric random vectors

The concept of elliptically symmetric random vectors emerged in the field of probability theory in the early 1970s as an extension of the class of multivariate normal distributions. Extensive research on elliptical distributions and their fundamental properties has been conducted, and a comprehensive survey of these results can be found in the book by Fang et al. [START_REF] Fang | Spherically and elliptically symmetric distributions[END_REF]. This survey provides a comprehensive overview of the basic concepts, key findings, and properties of elliptical distributions, making it an invaluable resource for further exploration and study in this area. We consider in this chapter that the random vectors are elliptically distributed. We then recall some definitions and properties of the elliptical family. Definition 19. An n-dimensional vector X follows an elliptical symmetric distribution if there exists a function Ψ such that its characteristic function is given by

ϕ(z) = E[e iz T X ] = e iz T µ Ψ(z T Σz). (5.2)
where µ is the location parameter, and Σ is a positive-definite matrix. The function Ψ is called a characteristic generator of the elliptical distribution. We note X ∼ Ellip(µ, Σ, ϕ). Definition 20. When the density function of an elliptical distribution exists, it must have the structure

f (x) = C √ detΣ g (x -µ) T Σ -1 (x -µ) . (5.3)
where g : R + → R ++ is a so-called radial density and c > 0 is a normalization factor ensuring that f integrates to one.

Remark 4. Table 5.2 gives a selection of some multivariate elliptical distributions, together with their characteristic generators, radial densities, quantile functions of the standard distribution, and the convexity conditions of the quantile function. The quantile functions of the selected distributions are given in Figure 5.2. In this section, we study the following joint chance-constrained geometric program

. Geometric programs

. Deterministic equivalent formulation

Law

Characteristic generator Radial density Quantile Function Convexity Condition Normal

e -1 2 t e -1 2 t 2 √ 2erf -1 (2α -1) ≤ α ≤ 1 Laplace (1 + 1 2 t) -1 e - √ 2|t| ln(2α), if 0 ≤ α ≤ 1 2 ≤ α ≤ 1 -ln(2(1 -α)), if 1 2 ≤ α ≤ 1 Cauchy e - √ t (1 + t 2 ) -n+1 2 tan(π(α -1 2 )) ≤ α ≤ 1 Logistic 2π √ t e π √ t -e -π √ t e -t 2 (1+e -t 2 ) 2 ln( α 1-α ) ≤ α ≤ 1 Table 5.2: Table of selected elliptical distributions
min t∈R M ++ E   I 0 i=1 c i M j=1 t a ij j   , s.t P   I k i=1 c i M j=1 t a ij j ≤ 1, k = 1, ...., K   ≥ α, (5.4)
with α ∈ [0.5, 1).

Assumption 35.

We assume that the row vectors are dependent and elliptically distributed and that the dependence is captured by a Gumbel-Hougaard copula, i.e., C 0 = (c 1 , ..., c I 0 ) ∼ Ellip(µ 0 , Σ 0 , ϕ) and

C k = (c 1 , ..., c I k ) ∼ Ellip(µ k , Σ k , ϕ).
Where µ k is a mean vector and

Σ k =    σ 11 σ 12 . . . σ 1I k . . . . . . . . . . . . σ I k 1 σ I k 2 . . . σ I k I k    is a positive definite matrix, k = 1, ..., K.
Assumption 36. In this research, we consider only the elliptical distributions present in Table 5.2.

Using the theory of copula and Sklar's properties, a deterministic equivalent of (5.4) is given by [START_REF] Cheng | Chance constrained 0-1 quadratic programs using copulas[END_REF]

min t∈R M ++ I 0 i=1 µ i M j=1 t a ij j , s.t I k i=1 µ i M j=1 t a ij j + ϕ -1 ((α y k ) 1 θ ) I k i=1 I k l=1 σ il M j=1 t a ij +a lj j ≤ 1, k = 1, ..., K, (5.5) K k=1 y k ≥ α, 0 < y k ≤ 1 , k = 1, ..., K,
where ϕ is one-dimensional standard elliptical distribution. Problem (5.5) is not convex, so we apply a logarithmic transformation by setting t j = e r j . This transformation allows us to reformulate problem (5.5) as follows min

I 0 i=1 µ i exp    M j=1 a ij r j    , s.t I k i=1 µ i exp    M j=1 a ij r j    + ϕ -1 ((α y k ) 1 θ ) i∈I k l∈I k σ il exp    M j=1 (a ij + a lj )r j    ≤ 1, k = 1, ..., K, K k=1 -ln(y k ) ≤ -ln(α), (5.6) -y k < 0 , k = 1, ..., K, y k ≤ 1 , k = 1, ..., K.
Theorem 37. The resulting problem (5.6) is biconvex.

Proof. The convexity on r is straightforward. We have for Gumbel-Hougaard copula that θ ≥ 1, there follows that 0

< 1 θ ≤ 1. Since 0 < y k ≤ 1, we obtain 0 < y k 1 θ ≤ 1. We have then, (α y k ) 1 θ ≥ α ≥ 0.5. From Table 5.2, y k → ϕ -1 ((α y k ) 1 θ ) is convex. The conclusion follows. Theorem 38. Let (r * , y * ) ∈ R M × R K and g k (r, y k ) = ϕ -1 ((α y k ) 1 θ ) i∈I k l∈I k σ kil exp M j=1 (a kji + a kjt )r j . If there exist γ r , γ y , β + , β -and λ such that      i∈I 0 µ i a i1 exp M j=1 a ij r j . . . i∈I 0 µ i a iM exp M j=1 a ij r j      + K k=1 γ rk      i∈I k a i1 exp M j=1 µ i a ij r j + ∇ r 1 g k (r, y k ) . . . i∈I k µ i a iM exp M j=1 a ij r j + ∇ r M g k (r, y)      = 0, i∈I k γ y k ∇ y k g k (r, y) + β + k -β -k - λ y k = 0, k = 1, ..., K (5.7)
γ r ≥ 0, i∈I k γ rk    i∈I k µ i exp    M j=1 a ij r j    + g k (r, y)    = 0, γ y ≥ 0, i∈I k γ y k    i∈I k µ i exp    M j=1 a ij r j    + g k (r, y) β + ≥ 0, β + k (y k -1) = 0, β -≥ 0, β -k (y k ) = 0, λ ≥ 0, λ ≥ 0, λ(ln(α)- K k=1 ln(y k )) = 0, k = 1, ..., K
then (r * , y * ) is a partial KKT point of (5.5).

Let (r * , y * ) ∈ R M × R K . We recall that if (5.5) is satisfied with partial Slater constraint qualification at (r * , y * ), then (r * , y * ) is a partial optimum of (5.5) if and only if (r * , y * ) is a partial KKT point of (5.5). Furthermore if γ r = γ y , then (r * , y * ) is a KKT point of (5.5).

. The dynamical neural network

We propose the following dynamical neural network to solve program (5.5)

dr j dt = - i∈I 0 a ij exp    M j=1 a ij r j    + K k=1   γ k +    i∈I k µ i exp    M j=1 a ij r j    + g k (r, y)      + ×   i∈I k a ij µ i exp    M j=1 a ij r j    + ∇ r j g k (r, y k )   , j = 1, ..., M, dy k dt = - i∈I k ∇ y k g k (r, y) ×   γ k +    i∈I k µ i exp    M j=1 a ij r j    + g k (r, y)      + + β + k (β + k + y k -1) + + β -k (β -k -y k ) + - λ y k λ + ln(α) - K k=1 ln(y k ) + , k = 1, ..., K, dγ k dt =   γ k +    i∈I k µ i exp    M j=1 a ij r j    + g k (r, y)      + -γ k , k = 1, ..., K, (DNN)
dβ + k dt = β + k + y k -1 + -β + k , k = 1, ..., K, dβ -k dt = β -k -y k + -β -k , k = 1, ..., K, dλ dt = λ + ln(α) - K k=1 ln(y k ) + -λ.
Let r = (r 1 , ..., r M ) T , y = (y 1 , ..., r K ) T , γ = (γ 1 , ..., γ K ) T , β + = (β + 1 , ..., β + K ) T and β -= (β -1 , ..., β -K ) T . For the sake of simplicity, we can write (DNN) as

dr dt = A(r, y, γ), (5.8)
dy dt = B(r, y, γ, β + , β -, λ), (5.9)
dγ dt = C(r, y, γ), (5.10)
dβ + dt = D(y, β + ), (5.11
)

dβ - dt = E(y, β -), (5.12)
dλ dt = F (y, λ).

(5.13)

A simplified circuit implementation of (5.8)-(5.13) is given in Figure 5.2.

Theorem 39. Suppose that (r * , y * ) is a partial optimum of (5.6) and γ * , β * + , β * -and λ * the corresponding Lagrange multipliers, then (r * , y * , γ * , β * + , β * -, λ * ) is a an equilibrium point of (DNN). Additionally, every equilibrium point of (DNN) is a KKT point of (5.6). Proof. Let µ = (r, y, γ, β + , β -, λ), we consider the following Lyapunov function

V (µ) = ||Φ(µ)|| 2 2 + 1 2 ||µ -µ * || 2 2 ,
where

Φ(y) =         A(r, y, γ) B(r, y, γ, β + , β -, λ) C(r, y, γ) D(y, β + ) E(y, β -) F (y, λ)        
. Similar to the analysis of Theorem 46, the Jacobian matrix ∇Φ is negative semidefinite. We have dV (µ)

dt = ( dΦ dt ) T Φ + Φ T dΦ dt + (µ -µ * ) T dµ dt . Observe that dΦ dt = dΦ dµ dµ dt = ∇Φ(µ)Φ(µ), we write then dV (µ) dt = Φ T (∇Φ(µ) T + ∇Φ(µ))Φ ≤0 since ∇Φ is negative semidfinite + (µ -µ * ) T (Φ(µ) -Φ(µ * ))
≤0 by Lemma 66

.

There follows that dV (µ) dt ≤ 0, thus V (µ) is a global Lyapunov function for the dynamical system (DNN) and (DNN) is stable in the sense of Lyapunov.

According to LaSalle's invariance principle, the trajectories µ(t) of system (DNN) converge to the largest invariant set M = µ | dV (µ) dt = 0 . It is easy to see that dµ dt = 0 if and only if dV dt = 0, therefore, the proposed neural network converges globally to the solution set of the problem (5.5).

. Numerical experiments

To assess the performance of our dynamical neural network, we employed the multidimensional shape optimization problem with joint chance constraints from [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF].

min x∈R M ++ m i=1 x -1 i , s.t P   m-1 j=1 ( m -1 A wallj x 1 m i=1,i̸ =j x i ), 1 A f loor m j=2 x j ≤ 1   ≥ 1 -ϵ, (5.14) 1 γ i,j x i x -1 j ≤ 1, 1 ≤ i ̸ = j ≤ m.
In our numerical experiments, we fixed the following parameters: m = 10, 1 γ i,j = 0.5, ϵ wall = 0.15, ϵ f loor = 0.15 and ϵ = 0.15. The inverse of floor's area ( 1A f loor ) and the inverse of wall area ( 1

A wall j
) for each j = 1, ..., m were considered as dependent variables. The mean value of 1 A f loor was set to 1.0/20.0 and the mean values of 1 A wall j , j = 1, ..., m were randomly selected from the range [1.0/60.0, 1.0/40.0]. The value of the coefficients of the covariance matrix is set to 0.01. We first compare the results of our approach with those of the sequential convex approximation algorithm from [START_REF] Shiraz | Copula theory approach to stochastic geometric programming[END_REF] for θ = 10.0 (i.e, high dependency). We test the robustness of the different approaches by creating 100 random samples of the variables 1 A wall j and 1

A f loor using the same mentioned moments. We then examine if the solutions from the three methods meet the constraints of (5.14) for all 100 cases. If the solutions are not viable for a particular case, it is referred to as a violated scenario (VS). The numerical results are displayed in Table 5. of the normal distribution, Laplace distribution and the logistic distribution, respectively.

For the three tables, column one gives the number of variables m. Columns two and three give the objective value and the number of VS when θ = 1, respectively. Columns four and five show the objective value and the number of VS when θ = 2, respectively.

Finally, columns six and seven present the objective value and the number of VS when θ = 10, respectively. We observe that for the different distributions, the objective value decreases as the dependency between the random variables becomes higher. The last observation is coherent since the joint chance constraint becomes less restrictive as θ increases. Nevertheless, the number of violated scenarios increases for the problems with high dependency.

. Linear programs as a particular case

This section was published in Results in Control and Optimization Journal [START_REF] Tassouli | Solving linear programs with joint probabilistic constraints with dependent rows using a dynamical neural network[END_REF]. We study now the following stochastic linear programming problem min c T x,

s.t. P(T k x ≤ D k , k = 1, .., K) ≥ 1 -α, (5.15)
x ≥ 0, where T k , k = 1, .., K are random vector rows, D k , k = 1, .., K are constants and 1 -α is a given confidence level. In our study, the dependence among the random rows Tk, k = 1, .., K, is modeled using the Gumbel-Hougaard copula with a parameter θ. Each row vector Tk is characterized by a known mean vector µ k = (µ k 1 , ...., µ kn ) T and a covariance m θ = 1 θ = 2 θ = 10

(

. Deterministic formulation equivalent

Problem (5.15) can be equivalently formulated as [START_REF] Cheng | Chance constrained 0-1 quadratic programs using copulas[END_REF] min c T x,

s.t. µ T k x + ϕ -1 ((1 -α) y k 1 θ ) x T Σ k x ≤ D k , k = 1, .., K, K k=1 y k = 1, (5.16
)

y k ≥ 0, k = 1, .., K, x ≥ 0,
where ϕ is one-dimensional standard normal distribution.

The convexity of the feasible set of problem (5.16) is guaranteed under certain conditions, as stated in the following Theorem. The detailed proof of this Theorem can be found in the work of Cheng et al. [START_REF] Cheng | Chance constrained 0-1 quadratic programs using copulas[END_REF].

Theorem 41. [27]. If 1 -α ≥ max(ϕ(u * , v * )) where u * = max k r k + 1 r k -1 (λ k min ) -1 2 ||µ k || and v * = max k r k + 1 (r k + 1) + ( 1-θ 2pln(p) f ϕ ( √ r k + 1)) 2 + 1-θ 2pln(p) f ϕ ( √ r k + 1)
, for some r k ≥ 1 with λ min is the smallest eigenvalue of Σ k , ϕ is the standard normal distribution and f ϕ is the associated density, then the feasible set of ( 8)

X(α) = x ≥ 0 | ∀k = 1, .., K, ∃y k ≥ 0 : µ T k x + ϕ -1 ((1 -α) y k 1 θ ) √ x T Σx ≤ D k , K k=1 y k = 1 is convex.
Corollary 42. By Theorem 41, problem (5.16) is biconvex.

Let f (x) = c T x, g θ k (x, y) = µ T k x + ϕ -1 ((1 -α) y k 1 θ ) x T Σ k x -D k , l(y) = ( K k=1 y k - 1, 1 -K k=1 y k ) T , h k (y) = -y k and w(x) = -x.
We write problem (5.15) as

min f (x), s.t. g θk (x, y) ≤ 0, k = 1, .., K, l(y) ≤ 0, (5.17)
h k (y) ≤ 0, k = 1, .., K, w(x)
≤ 0, Theorem 43. Let U the feasible set of (5.16), let (x * , y * ) ∈ U. If there are µ (1) , µ (2) , γ, λ, ω such that

∇ x f (x * ) + µ (1) T ∇ x g θ (x * , y * ) + ω T ∇ x w(x * ) = 0, (5.18)
µ (1) ≥ 0, µ (1) T g θ (x * , y * ) = 0, ω ≥ 0, ω T w(x * ) = 0, (5.19)
µ (2) T ∇ y g θ (x * , y * ) + γ T ∇ y l(y * ) + λ T ∇ y h(y * ) = 0, (5.20)
µ (2) ≥ 0, µ (2) T g θ (x * , y * ) = 0, γ ≥ 0, γ T l(y * ) = 0, (5.21)
λ ≥ 0, λ T h(y * ) = 0, (5.22)
with g θ (x, y) = (g θ 1 (x, y), . . . , g θ K (x, y)) T and h(y) = (h 1 (y), . . . , h K (y)) T , then (x * , y * ) is a partial KKT point of (5.17).

As stated in the previous Section, if (5.17) is satisfied with partial Slater constraints qualificationt (x * , y * ), then (x * , y * ) is a partial optimum of (5.17) if and only if (x * , y * ) is a partial KKT point of (5.17). Moreover, if µ (1) = µ (2) then (x * , y * ) is a KKT point of (5.17).

. The dynamical neural network

The neural network to solve (5.17) is driven by the following differential system

dx dt = -(∇ x f (x) + ∇ x g θ (x, y) T (µ + g θ (x, y)) + + ∇ x w(x) T (ω + w(x)) + ), (5.23)
dy dt = -(∇ y g θ (x, y) T (µ + g θ (x, y)) + + ∇ y l(y) T (γ + l(y)) + + ∇ y h(y) T (λ + h(y)) + ), (5.24)
dµ dt = (µ + g θ (x, y)) + -µ, (5.25)
dγ dt = (γ + l(y)) + -γ, (5.26
)

dλ dt = (λ + h(y)) + -λ, (5.27
)

dω dt = (ω + w(x)) + -ω.
(5.28)

We rewrite the dynamical system as

dz dt = Φ(z) z(t 0 ) = z 0
, where z = (x, y, µ, γ, λ, ω) and

Φ(z) =         -(∇ x f (x) + ∇ x g θ (x, y) T (µ + g θ (x, y)) + + ∇ x w(x) T (ω + w(x)) + ) -(∇ y g θ (x, y) T (µ + g θ (x, y)) + + ∇ y l(y) T (γ + l(y)) + + ∇ y h(y) T (λ + h(y)) + ) (µ + g θ (x, y)) + -µ (γ + l(y)) + -γ (λ + h(y)) + -λ (ω + w(x)) + -ω        
. The neural network system given by equations (5.23)-(5.28) can be implemented as an autonomous differential system using a generalized circuit. Figure 5.3 illustrates a generalized circuit implementation of the neural network. To begin the computation, we initialize the values of the variables x, y, µ, γ, λ, and ω with the initial values x 0 , y 0 , µ 0 , γ 0 , λ 0 , and ω 0 , respectively. We then solve equations (5.23)-(5.28) iteratively until convergence is achieved.

The following Theorem establishes the equivalence between an equilibrium point of the neural network system (5.23)-(5.28) and a Karush-Kuhn-Tucker (KKT) point of the optimization problem (5.16).

Theorem 44. Let (x * , y * , µ * , γ * , λ * , ω * ) an equilibrium point of (5.23)-(5.28), then (x * , y * ) is a partial KKT point of (5.16) with µ (1) = µ (2) . If (x * , y * ) is a partial optimum of (5.16), then there exist (µ * , γ * , λ * , ω * ) such that (x * , y * , µ * , γ * , λ * , ω * ) is an equilibrium point of

∇ x f (x * ) + ∇ x g θ (x * , y * ) T (µ * + g θ (x * , y * )) + + ∇ x w(x * ) T (ω * + w(x * )) + = 0, ∇ y g θ (x * , y * ) T (µ + g θ (x * , y * )) + + ∇ y l(y * ) T (γ * + l(y * )) + + ∇ y h(y * ) T (λ * + h(y * )) + = 0, (µ * + g θ (x * , y * )) + -µ * = 0, (γ * + l(y * )) + -γ * = 0, (λ * + h(y * )) + -λ * = 0, (ω * + w(x * )) + -ω * = 0.
We observe that (µ * + g θ (x * , y * )) + -µ * = 0 if and only if µ * ≥ 0, g θ (x * , y * ) ≤ 0 and µ * T g θ (x * , y * ) = 0. The same observation applies for (γ * +l(y * )) + -γ * = 0, (λ * +h(y * )) +λ * = 0 and (ω * + w(x * )) + -ω * = 0. By substitution, we find the KKT system (5.18)- (5.22).

The proof of the converse is straightforward and the conclusion follows.

Theorem 45. For any initial point (x(t 0 ), y(t 0 ), µ(t 0 ), γ(t 0 ), λ(t 0 ), ω(t 0 )), there exists an unique continuous solution (x(t), y(t), µ(t), γ(t), λ(t), ω(t)) for (5.23)-(5.28).

Proof. Since ∇ x f (x), ∇ x g θ (x, y), ∇ y g θ (x, y), ∇ x w(x), ∇ y h(y) and ∇ y l(y) are continuously differentiable on open sets, then all the second terms of the differential equations (5.23)-(5.28) are locally Lipschitz continuous. According to the local existence of ordinary differential equations also known as Picard-Lindelöf Theorem [START_REF] Miller | Ordinary Differential Equations[END_REF], the neural network (5.23)-

(5.28) has a unique continuous solution (x(t), y(t), µ(t), γ(t), λ(t), ω(t)).

To show the stability of the neural network (5.23)-(5.28), we need to show the semidefiniteness of the Jacobian matrix ∇Φ(z).

Theorem 46. The Jacobian matrix ∇Φ(z) is negative semidefinite.

Proof. We consider the dynamical neural network (5.23)-(5.28). We assume that there

exist 0 ≤ p ≤ K, 0 ≤ q ≤ 2, 0 ≤ s ≤ K and 0 ≤ t ≤ n such that (µ + g θ ) + = (µ 1 + g θ1 (x, y), µ 2 + g θ2 (x, y), ....., µ p + g θp (x, y), 0, ...., 0 K-p ), (γ + l
) + = (γ 1 + l 1 (y), γ 2 + l 2 (y), ....., γ q + l q (y), 0, ...., 0

), (λ + h) + = (λ 1 + h 1 (y), λ 2 + h 2 (y), ....., λ s + h s (y), 0, ...., 0

K-s
),

(ω + w) + = (ω 1 + w 1 (x), ω 2 + w 2 (x), ....., ω t + w t (x), 0, ...., 0 n-t

).

We write the Jacobian matrix of ϕ as follows

∇Φ(z) =         A 1 A 2 A 3 0 0 A 6 B 1 B 2 B 3 B 4 B 5 0 C 1 C 2 C 3 0 0 0 0 D 2 0 D 4 0 0 0 E 2 0 0 E 5 0 F 1 0 0 0 0 F 6        
, where

A 1 = -(∇ 2 x f (x) + p i=1 ((µ i + g θi )∇ 2 x g θ p i (x, y)) + ∇ x g θ p (x, y) T ∇ x g θ p (x, y) + t i=1 ((ω i + w i )∇ 2 x w t i (x)) + ∇ x w t (x) T ∇ x w t (x)), A 2 = -( p i=1 ((µ i + g θi )∇ y ∇ x g θ p i (x, y)) + ∇ y g θ p (x, y) T ∇ x g θ p (x, y)), B 1 = -( p i=1 ((µ i + g θi )∇ x ∇ y g θ p i (x, y)) + ∇ x g θ p (x, y) T ∇ y g θ p (x, y)), B 2 = -( p i=1
((µ i + g θi )∇ 2 y g θ p i (x, y)) + ∇ y g θ p (x, y) T ∇ y g θ p (x, y) + q i=1 ((γ i + l)∇ 2 y l q i (y))

+ ∇ y l q (y) T ∇ y l q (y) + s i=1

((λ i + h)∇ 2 y h s i (y)) + ∇ y h s (y) T ∇ y h s (y)), A 3 = -∇ x g p θ (x, y) T , A 6 = -∇ x w t (x) T , B 3 = -∇ y g p θ (x, y) T B 4 = -∇ y l q (y) T , B 5 = -∇ y h s (y) T , C 1 = ∇ x g p θ (x, y), C 2 = ∇ y g p θ (x, y), D 2 = ∇ y l q (y), E 2 = ∇ y h s (y), F 1 = ∇ x w t (x), C 3 = S p = - O p×p O p×(K-p) O (K-p)×p I (K-p)×(K-p) , D 4 = S q = - O q×q O q×(2-q) O (2-q)×q I (2-q)×(2-q) , E 5 = S s = - O s×s O s×(K-s) O (K-s)×q I (K-s)×(K-s)
and

F 6 = S t = - O t×t O t×(n-t) O (n-t)×t I (n-t)×(n-t) . The Jacobian matrix becomes ∇Φ(z) =         A 1 B T 1 A 3 0 0 A 6 B 1 B 2 B 3 B 4 B 5 0 -A T 3 -B T 3 S p 0 0 0 0 -B T 4 0 S q 0 0 0 -B T 5 0 0 S s 0 -A T 6 0 0 0 0 S t         . Let A = A 1 B T 1 B 1 B 2 , B = A 3 0 0 A 6 B 3 B 4 B 5 0 and S =    
S p 0 0 0 0 S q 0 0 0 0 S s 0 0 0 0 S t     , we have then

∇Φ(z) = A B -B T S .
It is clear that S is negative semidefinite. Furthermore, since g θ is biconvex, l, h and w are convex then ∇ 2

x g θ p i (x, y), ∇ 2 x w p i (x), ∇ 2 y g θ p i (x, y), ∇ 2 y l q i (y) and ∇ 2 y h s i (y) are positive semidefinite. It is easily proved that ∇ x g θ p (x, y) T ∇ x g θ p (x, y), ∇ x w t (x) T ∇ x w t (x), ∇ y g θ p (x, y) T ∇ y g θ p (x, y), ∇ y l q (y) T ∇ y l q (y) and ∇ y h s (y) T ∇ y h s (y) are positive semidefinite. There follows that A 1 and B 2 are negative semidefinite; hence, A is negative semidefinite. Since A and S are negative semidefinite, we conclude that ∇Φ is negative semidefinite [START_REF] Foias | Positive Definite Block Matrices[END_REF].

Theorem 47. The neural network (5.23)-(5.28) is stable and converges to a partial optimum of (5.16).

Proof. Let z * = (x * , y * , µ * , γ * , λ * , ω * ) an equilibrium point for the neural network (5.23)- (5.28). We consider the following Lyapunov function

E 1 (z) = ∥Φ(z)∥ 2 + 1 2 ∥z -z * ∥ 2 we have that dΦ dt = ∇Φ ∇z dz dt = ∇Φ(z)Φ(z), then dE 1 (z(t)) dt = ( dΦ dt ) T Φ + Φ T dΦ dt + (z -z * ) T dz dt = Φ T (∇Φ(z) T + ∇Φ(z))Φ + (z -z * ) T Φ(z).
By Theorem 46, we have Φ T (∇Φ(z) T + ∇Φ(z))Φ ≤ 0. By Theorem 46 and Lemma 24, we have (z As E 1 (z) ≥ 1 2 ∥z -z * ∥ 2 , then there exists a convergent subsequence (z(t k ) k≥0 ) where t k -→ ∞ when k -→ ∞, such that lim k-→∞ z(t k ) = ẑ, where ẑ satisfies dE 1 (z(t)) dt = 0. Notice that ẑ is a w-limit point of {z(t)} t≥t 0 . We have by LaSalle's invariance principle [128] that the neural network converges to the largest invariant set contained in M which is defined by M = {z(t)| dE 1 (z) dt = 0}. We observe that ( dx dt = 0, dy dt = 0, dµ dt = 0, dγ dt = 0, dλ dt = 0 and dω dt = 0) ⇔ dE 1 (z) dt = 0, therefore, ẑ is an equilibrium point for the neural network (5.23)- (5.28).

-z * ) T Φ(z) = (z -z * ) T (Φ(z) -Φ(z * )) ≤ 0 then dE 1 (z(t)) dt ≤ 0. E 1 (z) is positive
We now define a new Lyapunov function

E 2 (z) = ∥Φ(z)∥ 2 + 1 2 ∥z -ẑ∥ 2 E 2 (z) is continuously differentiable, E 2 (ẑ) = 0 and lim k-→∞ z(t k ) = ẑ then lim t-→∞ E 2 (z(t)) = E 2 (ẑ). We have also dE 2 (z) dt ≤ 0, then 1 2 ∥z -ẑ∥ 2 ≤ E 2 (z(t))
. Hence, lim t-→∞ ∥z -ẑ∥ = 0 and lim t-→∞ z(t) = ẑ. Therefore, the neural network (5.23)-(5.28)) is convergent in the sense of Lyupanov to an equilibrium point ẑ = (x, ŷ, μ, γ, λ, ω) where (x, ŷ) is a partial optimum of problem (5.16).

. Numerical experiments

To assess the performance of our approach, we conducted computational tests on a standard profit maximization problem. We compared the obtained results with a tangent approximation method proposed by Cheng et al. [START_REF] Cheng | A second-order cone programming approach for linear programs with joint probabilistic constraints[END_REF].

A manufacturing firm operates three machines to produce three different products. The manufacturing time for each unit of the products follows a normal distribution. Table 5.7 provides the mean and standard deviation values for the manufacturing time of each product, as well as the daily capacity of the machines. The first column lists the machine names, while columns 2 to 7 display the mean and standard deviation values for each product. The machine capacities are listed in column 8.

Our objective is to find the optimal daily production quantities for each product, taking into account the available machining times. The profit per unit for product 1, 2, and 3 is given by p 1 = 60, p 2 = 65, and p 3 = 62, respectively. We can formulate the minimization problem as follows min -60x

t21 x 1 + t22 x 2 + t23 x 3 + ϕ -1 ((0.95) y 2 1 θ ) σ 2 21 x 2 1 + σ 2 22 x 2 2 + σ 2 23 x 2 3 ≤ 660, t31 x 1 + t32 x 2 + t33 x 3 + ϕ -1 ((0.95) y 3 1 θ ) σ 2 31 x 2 1 + σ 2 32 x 2 2 + σ 2 33 x 2 3 ≤ 710, y 1 + y 2 + y 3 = 1, x 1 , x 2 , x 3 , y 1 , y 2 , y 3 ≥ 0.
We consider that [t 11 , t 12 , t 13 ] T , [t 21 , t 22 , t 23 ] T and [t 31 , t 32 , t 33 ] T are dependent and the dependence is driven by Gumbel-Hougaard copula. We solve problem (5.29)-(5.31) using the proposed neural network for θ = 2 and θ = 10. We compare the results obtained with the dynamical neural network with those obtained using the piecewise tangent approximation with 20 tangent points [START_REF] Cheng | A second-order cone programming approach for linear programs with joint probabilistic constraints[END_REF]. In order to evaluate the robustness of our approach, we generate a set of 100 instances of the stochastic vectors [t 11 , t 12 , t 13 ] T , [t 21 , t 22 , t 23 ] T and [t 31 , t 32 , t 33 ] T with the data in Table 5.7 and the function numpy.random.multivariate_normal.

Then, we check whether the two solutions are feasible for the constraints (5.29)-(5.31) for the 100 scenarios. If the solutions are not feasible for a given scenario, we call such a scenario a violated scenario denoted (VS). 5.9.

We observe that for all the values of α the dynamical neural network converges using the three values of θ. The solutions are all robust as they don't exceed the number of permitted violations. We notice that the solutions are more conservative when α increases.

We fix α = 0.10 and we modify the number of machines k and the number of products n. We compare the results obtained using the neural network with the tangent approximation [START_REF] Cheng | Chance constrained 0-1 quadratic programs using copulas[END_REF]. We use 20 tangent points for the tangent approximation and we recapitulate the results in Table 5.10. We notice that the gaps remain tight i.e they do not exceed 0.15% in the same way as in the independent case.

-Distributionally Robust Optimization

. Introduction

Uncertainty in parameters is a common characteristic of real-world decision problems in engineering and management. This uncertainty can arise from factors such as limited data observability, noisy measurements, implementation variations, and prediction errors. To address this uncertainty, stochastic optimization, (SO) and robust optimization (RO) frameworks have been traditionally employed in decision-making processes.

While SO assumes that the decision-maker possesses complete knowledge about the underlying uncertainty through a known probability distribution and minimizes a cost function based on this distribution, RO assumes that the decision-maker lacks detailed distributional knowledge about the uncertainty, except for its support. In RO, the goal is to minimize the worst-case cost over an uncertainty set, as studied in the works of El Ghaoui and Lebret [START_REF] Ghaoui | Robust solutions to least-squares problems with uncertain data[END_REF], Ben-Tal and Nemirovski [START_REF] Ben-Tal | Robust convex optimization[END_REF], and Bertsimas and Sim [START_REF] Bertsimas | The price of robustness[END_REF].

Distributionally Robust Optimization (DRO) is a framework used to address optimization problems in the presence of uncertain probability distributions. The goal of DRO is to find robust solutions that perform well across a range of possible distributions, rather than assuming a single fixed distribution.

Recent papers have considered the use of distributionally robust approaches in transportation network optimization problems [START_REF] Dai | A distributionally robust chance-constrained approach for modeling demand uncertainty in green port-hinterland transportation network optimization[END_REF], multistage distribution system planning [START_REF] Zare | A distributionally robust chanceconstrained milp model for multistage distribution system planning with uncertain renewables and loads[END_REF], portfolio optimization problems [START_REF] Fonseca | Robust international portfolio management[END_REF][START_REF] Wang | Distributionally robust optimization with multivariate second-order stochastic dominance constraints with applications in portfolio optimization[END_REF], planning and scheduling [START_REF] Shang | Distributionally robust optimization for planning and scheduling under uncertainty[END_REF], risk measures [START_REF] Postek | Computationally tractable counterparts of distributionally robust constraints on risk measures[END_REF], multimodal demand problems [START_REF] Hanasusanto | Distributionally robust multiitem newsvendor problems with multimodal demand distributions[END_REF], appointment scheduling [START_REF] Zhang | Distributionally robust appointment scheduling with moment-based ambiguity set[END_REF], vehicle routine problems [START_REF] Ghosal | The distributionally robust chance-constrained vehicle routing problem[END_REF] and energy and reserve dispatch [START_REF] Ordoudis | Energy and reserve dispatch with distributionally robust joint chance constraints[END_REF].

. Geometric Programs

We consider the following geometric program. In this paper, we consider the case where the coefficients c i are not known. In conse-quence, we reformulate the optimization problem (6.1) as follows

min t∈R M ++ sup F 0 ∈D 0 E F 0   I 0 i=1 c 0 i M j=1 t a 0 ij j   , (JCP) s.t inf F ∈D P F   I k i=1 c k i M j=1 t a k ij j ≤ 1, k = 1, ...., K   ≥ 1 -ϵ,
where F 0 is the probabilistic distribution of vector C 0 = (c 0 1 , .., c 0 I 0 ) T , F is the joint distribution for C 1 = (c 1 1 , .., c 1 I 1 ) T , ..., C k = (c k 1 , .., c k I k ) T , D 0 is the uncertainty set for the probability distribution F 0 , D is the uncertainty set for the probability distribution F and 1 -ϵ, ϵ ∈ (0, 0.5], is the confidence parameter for the joint constraint. This paper examines the distributionally robust geometric programs (JCP) using two different groups sets of uncertainty. The first group focuses on uncertainties in distributions, considering both known and unknown first two order moments. The second group of uncertainty sets incorporates first order moments along with nonnegative support constraints.

. Uncertainty Sets with First Two Order Moments

We first consider that the mean vector of C k , k = 0, 1, ..., K lies in an ellipsoid of size γ k 1 ≥ 0 with center µ k and that the covariance matrix of C k , k = 0, 1, ..., K lies in a positive semidefinite cone of center Σ k = σ k i,j , i, j = 1, ..., I k . We define for every

k = 0, 1, ..., K, D 2 k (µ k , Σ k ) = F k (E F k [C k ] -µ k ) T Σ -1 k (E F k [C k ] -µ k ) ≤ γ k 1 COV F k (C k ) ⪯ γ k 2 Σ k
, where

F k is the probability distribution of C k , γ k 2 ≥ 0 and COV F k is a covariance operator under probability distribution F k of C k .
Based on whether the row vectors C k , k = 1, .., K are mutually independent or dependent, we have two cases.

Case (JCP) with Jointly Independent Row Vectors.

Assumption 48. We assume that D = {F|F = F 1 F 2 ...F K }, where F is the joint distribution for mutually independent random vectors C 1 , C 2 , ..., C K with marginals F 1 , F 2 , ..., F K . Theorem 49. Given Assumption 48, (JCP) is equivalent to

(JCP ind ) min t∈R M ++ ,y∈R K + I 0 i=1 µ 0 i M j=1 t a 0 ij j + γ 0 1 I 0 i=1 I 0 l=1 σ 0 i,l M j=1 t a 0 ij +a 0 lj j , (6.2)
s.t

I k i=1 µ k i M j=1 t a k ij j + γ k 1 I k i=1 I k l=1 σ k i,l M j=1 t a k ij +a k lj j + y k 1 -y k γ k 2 I k i=1 I k l=1 σ k i,l M j=1 t a k ij +a k lj j ≤ 1, k = 1, ..., K, (6.3) K k=1 y k ≥ 1 -ϵ, 0 < y k ≤ 1 , k = 1, ..., K. (6.4)
Proof. As the row vectors C k , k = 1, ..., K are mutually independent, (JCP) is written equivalently by introducing K nonegative auxiliary variables y k as [START_REF] Tassouli | A neural network approach to solve geometric programs with joint probabilistic constraints[END_REF].

min t∈R M ++ sup F 0 ∈D 0 E F 0   I 0 i=1 c 0 i M j=1 t a 0 ij j   , s.t inf F k ∈D k P F k   I k i=1 c k i M j=1 t a k ij j ≤ 1   ≥ y k , , k = 1, ...., K, K k=1 y k ≥ 1 -ϵ, 0 < y k ≤ 1 , k = 1, ..., K.
By Theorem 1 in [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF], we conclude that (JCP) is equivalent to (JCP ind ).

Problem (JCP ind ) is not convex. By applying the logarithmic transformation r j = log(t j ), j = 1, ..., M and x k = log(y k ), k = 1, ..., K, we have the following equivalent reformulation of (JCP ind )

(JCP log ind ) min r∈R M ,x∈R K I 0 i=1 µ 0 i exp    M j=1 a 0 ij r j    + γ 0 1 I 0 i=1 I 0 l=1 σ 0 i,l exp    M j=1 (a 0 ij + a 0 lj )r j    , (6.5)
s.t

I k i=1 µ k i exp    M j=1 a k ij r j    + γ k 1 I k i=1 I k l=1 σ k i,l exp    M j=1 (a k ij + a k lj )r j    + γ k 2 I k i=1 I k l=1 σ k i,l exp    M j=1 (a k ij + a k lj )r j + log e x k 1 -e x k    ≤ 1, k = 1, ..., K, K k=1 x k ≥ log(1 -ϵ), x k ≤ 0 , k = 1, ..., K.

Case (JCP) with Jointly Dependent Row Vectors.

In this case, (JCP) is equivalent to [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF] (JCP dep ) min

t∈R M ++ ,y∈R K + I 0 i=1 µ 0 i M j=1 t a 0 ij j + γ 0 1 I 0 i=1 I 0 l=1 σ 0 i,l M j=1 t a 0 ij +a 0 lj j , (6.7)
s.t

I k i=1 µ k i M j=1 t a k ij j + γ k 1 I k i=1 I k l=1 σ k i,l M j=1 t a k ij +a k lj j + y k 1 -y k γ k 2 I k i=1 I k l=1 σ k i,l M j=1 t a k ij +a k lj j ≤ 1, k = 1, ..., K, (6.8) K k=1 y k ≥ K -ϵ, 0 < y k ≤ 1 , k = 1, ..., K. (6.9)
As for the independent case, we obtain the following biconvex equivalent problem for (JCP dep )

(JCP log dep ) min r∈R M ,x∈R K I 0 i=1 µ 0 i exp    M j=1 a 0 ij r j    + γ 0 1 I 0 i=1 I 0 l=1 σ 0 i,l exp    M j=1 (a 0 ij + a 0 lj )r j    , (6.10)
s.t

I k i=1 µ k i exp    M j=1 a k ij r j    + γ k 1 I k i=1 I k l=1 σ k i,l exp    M j=1 (a k ij + a k lj )r j    + γ k 2 I k i=1 I k l=1 σ k i,l exp    M j=1 (a k ij + a k lj )r j + log y k 1 -y k    ≤ 1, k = 1, ..., K, K k=1 y k ≥ K -ϵ, 0 < y k ≤ 1 , k = 1, ..., K.
(6.11) Theorem 51. [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF] If ϵ ≤ 0.5 and σ k i,l ≥ 0 for all i, l and k, problem (JCP log dep ) is a convex programming problem.

. Uncertainty Sets with Known First Order Moment and Nonnegative Support

In this section, we consider uncertainty sets with nonnegative supports and known first-order moments. The uncertainty sets for (JCP) can be formulated as follows

D 3 k (µ k , Σ k ) = F k E[C k ] = µ k P F k [C k ≥ 0] = 1
, k = 0, 1, ..., K, where µ k > 0.

Case (JCP) with Jointly Independent Row Vectors.

We first consider the case when the marginal distributions in the uncertainty set are jointly independent. Using the strong duality [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF], (JCP) can be reformulated as follows

(JCP ind N S ) min t∈R M ++ ,λ,β,π I 0 i=1 µ 0 i M j=1 t a 0 ij j , (6.12)
s.t

K k=1 y k ≥ 1 -ϵ, 0 ≤ y k ≤ 1 , k = 1, ..., K, (6.13
)

y k λ -1 k -λ -1 k β kT µ k ≤ 1, k = 1, ...., K, (6.14)
β k ≤ 0, 0 < λ ≤ 1, k = 1, ...., K, (6.15)
λ -1 k π k ≥ 1, k = 1, ...., K, (6.16)
(-β k ) -1 π k M j=1 t a k ij j ≤ 1, i = 1, ..., I k , k = 1, ...., K, (6.17)
(JCP) can be reformulated as a convex problem using a logarithmic transformation x j = log(y j ), t j = log(r j ), λk = log(λ k ), βk = log(-β k ), π = log(π). Problem (JCP N S ) becomes, (JCP log N S-ind ) min

x,r, λ, β,π

I 0 i=1 µ 0 i exp    M j=1 a 0 ij r j    , (6.18)
s.t

K k=1 x k ≥ log(1 -ϵ), x k ≤ 0, k = 1, ..., K, (6.19)
exp(x k -λk ) + Case (JCP) with Jointly Dependent Row Vectors.

I k i=1 exp -λk + βk i + logµ k i ≤ 1, k = 1, ...., K, λk ≤ 0, k = 1, ...., K, (6.20
In the case where the constraints of (JCP) are jointly dependent, we have the following deterministic equivalent

(JCP dep N S ) min t∈R M ++ ,λ,β,π I 0 i=1 µ 0 i M j=1 t a 0 ij j , (6.23)
s.t

K k=1 y k ≥ K -ϵ, 0 ≤ y k ≤ 1 , k = 1, ..., K, (6.24
)

y k λ -1 k -λ -1 k β kT µ k ≤ 1, k = 1, ...., K, (6.25
)

β k ≤ 0, 0 < λ ≤ 1, k = 1, ...., K, (6.26)
λ -1 k π k ≥ 1, k = 1, ...., K, (6.27)
(-β k ) -1 π k M j=1 t a k ij j ≤ 1, i = 1, ..., I k , k = 1, ...., K, (6.28)
The log transformation fails to convert (JCP For (JCP where z(.) and γ(.) are two time-dependent variable, κ is a given convergence rate and (x) + = max(x, 0).

log ind ), z = (r, x) T , f (z) = I 0 i=1 µ 0 i exp M j=1 a 0 ij r j + γ 0 1 I 0 i=1 I 0 l=1 σ 0 i,l exp M j=1 (a 0 ij + a 0 lj )r j and g(z) =                                I 1 i=1 µ 1 i exp M j=1 a 1 ij r j + γ 1 1 I 1 i=1 I 1 l=1 σ 1 i,l exp M j=1 (a 1 ij + a 1 lj )r j + e x 1 1-e x 1 γ 1 2 I 1 i=1 I 1 l=1 σ 1 i,l exp M j=1 (a 1 ij + a 1 lj )r j -1 . . . I K i=1 µ K i exp M j=1 a K ij r j + γ K 1 I K i=1 I K l=1 σ K i,l exp M j=1 (a K ij + a K lj )r j + e x K 1-e x K γ K 2 I K i=1 I K l=1 σ K i,l exp M j=1 (a K ij + a K lj )r j -1 log(1 -ϵ) -K k=1 x k x 1 log dep ), z = (r, x) T , f (z) = I 0 i=1 µ 0 i exp M j=1 a 0 ij r j + γ 0 1 I 0 i=1 I 0 l=1 σ 0 i,l exp M j=1 (a 0 ij + a 0 lj )r j and g(z) =                                I 1 i=1 µ 1 i exp M j=1 a 1 ij r j + γ 1 1 I 1 i=1 I 1 l=1 σ 1 i,l exp M j=1 (a 1 ij + a 1 lj )r j + e x 1 1-e x 1 γ 1 2 I 1 i=1 I 1 l=1 σ 1 i,l exp M j=1 (a 1 ij + a 1 lj )r j -1 . . . I K i=1 µ K i exp M j=1 a K ij r j + γ K 1 I K i=1 I K l=1 σ K i,l exp M j=1 (a K ij + a K lj )r j + e x K 1-e x K γ K 2 I K i=1 I K l=1 σ K i,l exp M j=1 (a K ij + a K lj )r j -1 log(K -ϵ) -K k=1 x k x 1       log(1 -ϵ) - K k=1 x k x 1 . . . x K exp(x 1 -λ1 ) + I 1 i=1 exp -λ1 + β1 i + logµ 1 i -1 . . . exp(x K -λK ) + I K i=1 exp -λK + βK i + logµ K i -1 λ1 . . . λ1 -π1 . . . λK -πK π1 + M j=1 a 1 ij r j -β1 i ≤ 0, i = 1, ..., I 1 (
Theorem 52. If (z * , γ * ) is an equilibrium point of (6.37)-(6.38) if and only if z * is an optimal solution of (6.34) where γ * is the correspondent Lagrange multiplier.

Proof. Let (z * , γ * ) is an equilibrium point of (6.37)-(6.38), then dz * dt = 0 and dγ * dt = 0.

dz * dt = 0 ⇔ ∇f (z * ) + ∇g(z * ) T (γ * + g(z * )) + = 0, (6.39)
dγ * dt = 0 ⇔ -γ * + (γ * + g(z * )) + = 0 (6.40)

Observe that γ * = (γ * + g(z * )) + if and only if γ * ≥ 0, g(z * ) ≤ 0 and γ * T g(z * ) = 0, we obtain then (6.130) of the KKT system (6.130)- (6.36). Furthermore, we replace (γ * + g(z * )) + by γ * in the right hand side of (6.39) we obtain then ∇f (z * ) + ∇g(z * ) T γ * = 0 which is equation (6.36) of the KKT system (6.130)- (6.36). For the converse part of the theorem, it is straightforward that if z * is an optimal solution of (6.34) where γ * is the correspondent Lagrange multiplier then (z * , γ * ) is an equilibrium point of (6.37)- (6.38).

Lemma 53. For any initial point (z(t 0 ), γ(t 0 )), there exists a unique continuous solution (z(t), γ(t)) for (6.37)- (6.38).

Proof. The right-hand side of system (6.37)-( 6

Proof. Let ζ = (z, γ), we define U (ζ) = -∇f (z) + ∇g(z) T (γ + g(z)) + -γ + (γ + g(z)) + .
First, consider the following Lyapunov function

E(ζ) = ||U (ζ)|| 2 + 1 2 ||ζ -ζ * ||, (6.41)
where ζ * = (z * , γ * ) is an equilibrium point of (6.37)- (6.38).

dE(ζ(t)) dt = dU dt T U + U T dU dt + (ζ -ζ * ) T dζ dt . Observe that dU dt = dU dζ × dζ dt = ∇U (ζ)U (ζ).
Without loss of generality suppose that there exists p ∈ N such that (γ + g(z)) + = (γ 1 + g 1 (z)), ..., (γ p + g p (z)), 0, ..., 0), and we define g p = (g 1 , ..., g p ).

We have ∇U

(ζ) = -∇ 2 f (z) + p i=1 ∇ 2 g p (z) (γ p + g p (z)) + ∇g(z) T ∇g(z) -∇g p (z) T ∇g p (z) S p .
where

S p = O p×p O p×(N -p) O (N -p)×q I (N -p)×(N -p)
, where N is the length of vector γ.

Since f and g are convex, then the hessian matrices ∇ 2 f and ∇ 2 g p are positive semidefinite. Furthermore ∇g T ∇g is positive semidefinite, we conclude that ∇U is negative semidefinite.

Back to the expression of

dE(ζ(t)) dt , we have dE(ζ(t)) dt = U T (∇U + ∇U T )U ≤0 since ∇U is negative semidfinite + (ζ -ζ * ) T (U (ζ) -U (ζ * ))
≤0 by Lemma 4 in [START_REF] Tassouli | A neural network approach to solve geometric programs with joint probabilistic constraints[END_REF] ≤ 0. So the neural network (6.37)-(6.38) is globally stable in the sense of Lyapunov.

Next similarly to the proof of Theorem 5 in [START_REF] Tassouli | A neural network approach to solve geometric programs with joint probabilistic constraints[END_REF], we prove that the dynamical neural network (6.37)-(6.38) is globally convergent to (z * , γ * ) where z * is the optimal solution of (6.34).

. A two-time scale neurodynamic duplex for (JCP

, z = (r, λ, β, π) T , f (z) = I 0 i=1 µ 0 i M j=1 t a 0 ij      K -ϵ - K k=1 y k -y 1 . . . -y K y 1 -1 . . . y K -1 y 1 exp(-λ1 ) + I 1 i=1 exp -λ1 + β1 i + logµ 1 i -1 . . . y K exp(-λK ) + I K i=1 exp -λK + βK i + logµ K i -1 λ1 . . . λ1 -π1 . . . λK -πK π1 + M j=1 a 1 ij r j -β1 i ≤ 0, i = 1, ..., I 1      .
We denote U = {z, y | g(z, y) ≤ 0} the feasible set of (6.42).The Lagrangian function of problem (6.42) is defined as follows:

L(z, y, ω) = f (z) + ω T g(z, y). (6.43)
For any (z, y) ∈ U, the KKT conditions are stated as follows:

∇L(z, y, ω) = 0, (6.44) ω ≥ 0, ω T g(z, y) = 0.

f (z) ≤ f (z), ∀z ∈ U y , (6.46)
where U y = {z | g(z, y) ≤ 0}.

Theorem 55. The KKT system (6.44)-(6.45) is equivalent to the following system ∇f (z) + ∇ z g(z, y) T (ω + g(x, z)) + = 0 (6.47)

∇ y g(z, y) T (ω + g(z, y)) + = 0 (6.48)

(ω + g(x, z)) + -ω = 0 (6.49)

Proof. The proof of Theorem 55 follows the same lines of the proof of Theorem 52.

Based on the equations (6.47)-(6.49), we consider the following two-time-scale recurrent neural network model

κ 1 dz dt = -∇f (z) + ∇ z g(z, y) T (ω + g(x, z)) + , (6.50)
κ 2 dy dt = -∇ y g(z, y) T (ω + g(z, y)) + , (6.51)
κ 2 dω dt = -ω + (λ + g(z, y)) + , (6.52)
where (z, y, ω) are now time-dependent variables and κ 1 and κ 2 are two time scaling constants with κ 1 ̸ = κ 2 . We propose a duplex of two two-time-scale recurrent neural network (6.50)-(6.52) for solving (6.34) one with κ 1 > κ 2 and the second with κ 1 < κ 2 as shown in When the value of µ approaches 1, the mutated element of the particle tends to move towards the maximum value of Λ i (j + 1). On the other hand, as µ approaches -1, the mutated element tends to move towards the minimum value of x i (j + 1). The magnitude of |µ| determines the size of the search space for x i (j + 1), with larger values indicating a wider search space. Conversely, smaller values of |µ| result in a smaller search space, allowing for fine-tuning.

To achieve fine-tuning, the dilation parameter a is adjusted based on the current iteration j relative to the total number of iterations T . Specifically, a is set as a function of j/T , with a = e 10 j T . Additionally, ϕ is randomly generated from the interval [-2.5a, 2.5a]. The algorithm details are given in Algorithm 2 where Λ = (z, y, ω) Algorithm 2 The neurodynamic duplex -Let Λ 1 (0) and Λ 2 (0) be randomly generated in the feasible region.

-Let Λ(0) = Λ(0) = y(0) the initial best previous position and best position, respectively.

-Set the convergence error ϵ. while ||Λ(j + 1) -Λ(j)||≥ ϵ do Compute the equilibrium points Λ1 (j) and Λ2 (j) of RNN 1 and RNN

2 . if f (z 1 (j)) < f (z(j)) then Λ(j + 1) = Λ1 (j) else Λ(j + 1) = Λ(j) end if if f (z 2 (j)) < f (z(j)) then Λ(j + 1) = Λ2 (j) else Λ(j + 1) = Λ(j) end if if f (z(j)) < f (ẑ(j)) then Λ(j + 1) = Λ(j + 1) else Λ(j + 1) = Λ(j) end if
Compute the value of Λ(j + 1) following (6.

κ dω dt = -ω + (λ + g(z, y)) + , (6.63)
converges to an equilibrium point, then the state of the neurodynamic model with two timescales, as described by equations (6.50)-(6.52), globally converges to a partial optimum of problem (6.42).

Proof. We recall the Lagrangian function of (6.42)

L(z, y, ω) = f (z) + ω T g(z, y). (6.64)
An equilibrium point (z * , y * , ω * ) of (6.61)-(6.63) corresponds to a KKT point of (6.42). We fix y * , and take z ∈ U y * , (6.42) becomes a convex optimization problem and we have

L(z * , y * , ω * ) ≤ L(z, y * , ω * ), (6.65)
which is equivalent to

f (z * ) + ω * T g(z * , y * ) ≤ f (z) + ω * T g(z, y * ). (6.66)
As ω * T g(z, y * ) ≤ ω * T g(z * , y * ) = 0, we have f (z * ) ≤ f (z). By Definition 21, (z * , y * ) is a partial optimum of 6.42.

Theorem 59. The duplex of two two-timescale neural networks in Figure 6.6 system is globally convergent to a global optimal solution of problem (6.34) with probability one.

Proof. By Theorem 58, the two-timescale neurodynamic models RNN 1 and RNN 2 are proven to converge to a partial optimum. From Algorithm 2, the solution sequence is generated as follows

Λ(j + 1) = Λ(j + 1) if f (z(j)) < f (ẑ(j)), Λ(j + 1) = Λ(j) else.
We observe that the generated solution sequence is monotonically increasing {f ( Λ(j))} ∞ j=1 .

Let Mi, j represent the supporting set of the initial state of RNNi at iteration j. According to equation (6.59), the mutation operation ensures that the initial states of the recurrent 123 neural networks are constrained to the feasible region U. Therefore, for every iteration index J ≥ 1, the supporting sets satisfy the following condition:

U ⊆ M = J j=1 2 i=1 M i,j . (6.67)
Consequently, we have v(U) = v(M) > 0.

By Lemma 66, we have

lim j->∞ P( Λ(j) ∈ Φ) = 1 (6.68)
where Φ is the set of the global optimal solutions of (6.34). The conclusion follows.

. Numerical experiments

We consider three geometric problems to evaluate the performance of our neurodynamic approaches. All the algorithms in this Section are implemented in Python. We run our algorithms on Intel(R) Core(TM) i7-10610U CPU @ 1.80GHz. The random instances are generated with numpy.random, and we solve the ODE systems with solve_ivp of scipy.integrate. The deterministic equivalent programs are solved with the package gekko and the gradients and partial derivatives are computed with autograd.grad and autograd.jacobian. For the following numerical experiments, we set γ k 1 = 2, γ k 2 = 2 and the error tolerance for the neurodynamic duplex ζ = 10 -4 . In the second subsection, we evaluate the quality of our neurodynamic duplex by comparing the obtained solutions with the ones given by the Convex Alternate Search (CAR) from [START_REF] Gorski Jochen | Biconvex sets and optimization with biconvex functions: a survey and extensions[END_REF]. The gap between the two solutions is computed as follows GAP = Sol CAR -Sol Duplex Sol CAR , where Sol CAR and Sol Duplex are the solutions obtained using the CAR and the neurodynamic duplex, respectively. For the neurodynamical duplex, we take κ 1 κ 2 = 0.1 for the first dynmical neural network and κ 1 κ 2 = 10.0 for the second one.

. Uncertainty Sets with First Two Order Moments

A three-dimension shape optimization problem

We first consider a transportation problem involving the shifting of grain from a warehouse to a factory. The grain is transported within an open rectangular box, with dimensions of length x 1 meters, width x 2 meters, and height x 3 meters, as illustrated in Figure 6.2. The objective of the problem is to maximize the volume of the rectangular box, given by the product of its length, width, and height (x 1 x 2 x 3 ). However, two constraints must be satisfied. The first constraint relates to the floor area of the box, and the second constraint relates to the wall area. These constraints are necessary to ensure that the shape of the box aligns with the requirements of a given truck. In our analysis, we assume that the wall area A wall and the floor area A f loor are random variables. We formulate our shape optimization problem as follows

++ x -1 1 x -1 2 x -1 3 , (6.69)
s.t inf

F ∈D P 1 A wall (2x 3 x 2 + 2x 1 x 3 ) ≤ 1, 1 A f loor x 1 x 2 ≤ 1 ≥ 1 -ϵ.
where F is the joint distribution for 1

A wall and 1

A f loor and D is the uncertainty set for the probability distribution F. We solve problem (6.69) when the uncertainty set is equal to D 2 using the dynamical neural network (6.37)-(6.38). For the numerical experiments, we take the mean and the covariance describing the uncertainty sets for 1 A wall m wall = 0.05, σ wall = 0.01, respectively and for 1

A f loor m f loor = 0.5, σ f loor = 0.1, respectively.

We recapitulate the obtained results in Table 6.1. Columns one, two and three give the optimal value, the CPU time and the number of violated scenarios (VS) in the independent case, respectively. Columns four, five and six show the optimal value, the CPU time and the number VS in the dependent case, respectively. We observe that the dependent case is more conservative compared to the independent and that the dynamic neural network covers well the risk region in both cases. Figure 6.3 show the convergence of the state variables.

A generalized shape optimization problem

To further assess the performance of our dynamical neural network, we employed the multidimensional shape optimization problem with joint chance constraints from [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF].

125 x -1 i ,

s.t inf F ∈D P F   m-1 j=1 ( m -1 A wallj x 1 m i=1,i̸ =j x i ), 1 A f loor m j=2 x j ≤ 1   ≥ 1 -ϵ, (6.70) 1 γ i,j x i x -1 j ≤ 1, 1 ≤ i ̸ = j ≤ m.
In our numerical experiments, we fixed the following parameters 1 γ i,j = 0.5 and ϵ = 0.15.

The inverse of floor's area ( 1 A f loor ) and the inverse of wall area ( 1

A wall j
) for each j = 1, ..., m were considered as random variables. We test the robustness of the different approaches by creating 100 random samples of the variables 1

A wall j and 1 A f loor . We then examine if the solutions meet the constraint of (6.70) for all 100 cases for the gaussian distribution for example. If the solutions are not feasible for a particular case, it is referred to as a violated scenario (VS).

We first solve (6.70) for m = 5 and when the uncertainty set is D 2 in the independent case for different initial points, we observe that the dynamical neural network (6.37)-(6.38) converges to the same final value independently from the starting value as shown in Figure 6.4. Now we solve (6.70) when we know the first-order moments of 1 A f loor and 1

A wall j for both the dependent and the independent case. We present the obtained results in Table 6.6. We observe again that the dependent case is more conservative compared to the independent one.

. Uncertainty Sets with Known First Order Moment and Nonnegative Support

A generalized shape optimization problem

We solve (6.70) when the uncertainty set is D 3 for both the independent and the dependent case. For the numerical experiments, we take ϵ = 0.2. We solve problem (6.70) using the neurodynamic duplex in the dependent case. We recapitulate the obtained results in Table 6.3. Column one gives the number of variables m. Columns two, three and four give the objective value, the CPU time and the number of VS in the independent case, respectively. Columns five, six and seven give the objective value, the CPU time and the number of VS for the dependent case, respectively. We observe that the problem with dependent variables is more conservative. Nevertheless, the solution, in this case, covers well the risk area as the number of VS is equal to 0 for all the values of m.

Now we additionally solve problem (6.70) using the assumption that the random variables follow a normal distribution [START_REF] Tassouli | A neural network approach to solve geometric programs with joint probabilistic constraints[END_REF] for m = 5. In order to compare the solutions obtained with the stochastic and the robust approaches, we evaluate the robustness of the solutions for different hypotheses on the true distribution of the random parameters, i.e, the uniform distribution, the normal distribution, the log-normal distribution, the logistic distribution and Gamma distribution. The obtained results are presented in Table 6.4 which gives the number of violated scenarios for both the normal solutions and the robust ones and the objective value obtained by each solution. We can infer that the distributionally robust approaches are a conservative approximation of the stochastic programs. We observe that the solutions obtained by the nonnegative support are more conservative compared to the stochastic ones. Notice that the distributionally robust solutions are

Maximizing the worst user signal-to-interference noise ratio

We consider the problem of maximizing the worst user signal-to-interference noise ratio (SINR) for Massive Multiple Input Multiple Output (MaMIMO) systems subject to antenna assignment and mul-tiuser interference constraints taken from [START_REF] Adasme | A stochastic geometric programming approach for power allocation in wireless networks[END_REF] and given by

max p∈ K ++ min i∈U p i |g H i g i | 2 j∈U ,j̸ =i p j |g H i g j | 2 +|σ i | 2 , (6.71)
s.t

P min ≤ p i ≤ P max , ∀i ∈ U, (6.72)
where p i is the power to be assigned for each user i ∈ U. g i ∈ C T ×1 , g H i ∈ C 1×T and σ 2 i are the beam domain channel vector associated to user i ∈ U, its Hermitian transpose and Additive White Gaussian Noise (AWGN), respectively.

Let

a ij = |g H i g j | 2 |g H i g i | -2 and b i = |σ i | 2 |g H i g i | -2 ,
a ij p j p -1 i w + b i p -1 i w ≤ 1, ∀i ∈ U, (6.74)
P min ≤ p i ≤ P max , ∀i ∈ U.

(6.75)

We assume that the coefficients a ij and b i are independent random variables and we propose the following optimization problem with individual and joint chance constraints

min p∈ K ++ ,w∈ ++ w -1 , s.t inf F i ∈D-i P F i    j∈U ,j̸ =i a ij p j p -1 i w + b i p -1 i w ≤ 1    ≥ 1 -ϵ i , ∀i ∈ U, (POI) P min ≤ p i ≤ P max , ∀i ∈ U.
and

min p∈ K ++ ,w∈ ++ w -1 , s.t inf F ∈D P F    j∈U ,j̸ =i a ij p j p -1 i w + b i p -1 i w ≤ 1, ∀i ∈ U    ≥ 1 -ϵ, (POJ)
P min ≤ p i ≤ P max , ∀i ∈ U.

We assume that the uncertainty set for the distributionally robust problems (POI) and

(POJ) is D 3 . We fix ϵ = 0.2. We first solve problem (POJ) for K = 10. Figure 6.5 shows the convergence of the power variables. Next, we solve (POI) and (POJ) for different values of the number of users K. Table 6.5 presents the obtained results. Column one gives the number of users K. Columns two and three give the optimal value and the number of VS for (POI), respectively. Columns four and five show the optimal value and the number of VS for (POJ), respectively. As observed in the previous section, the use of joint constraints leads to a more conservative minimization problem but covers well the risk area compared to the problem with individual constraints since the number of VS is lower. We consider in this section the case of linear programs. The distributionally robust problem can be formulated as follows.

. Linear programs as a particular case

min x∈ n + sup F 0 ∈D 0 E F 0 ζ0 T x , (6.76)
s.t inf

F ∈D P F ζk x ≤ b k , k = 1, ..., K ≥ α. (6.77)
where ζ0 ∈ n is an uncertain parameter, [ ζ1 , ζ2 , ... ζK ] T is a K × n set of pairwise independent random vectors in n and b ∈ K is a deterministic vector. We consider the case where the probability distribution F 0 of ζ0 belongs to a certain uncertainty set D 0 and the probability distributions F k of ζk , k = 1, ..., K are not completely known and belong to D k . Thus, we take the worst-case where constraints (6.77) are jointly satisfied for all possible distributions in a given distributional uncertainty set D with a given probability level α. Based on the pairwise independence between the vectors ( ζk ) k∈{1,..,K} , we introduce nonnegative auxiliary variables z k , k = 1, ..., K and rewrite the constraint (6.77) as

               inf F k ∈D k P F k ζk x ≤ b k ≥ α z k , k = 1, ..., K K k=1 z k = 1, z k ≥ 0, k = 1, ..., K. (6.78)
In the following subsections, we provide deterministic equivalent formulations for the optimization problem described in equations (6.76)-(6.77). To handle the uncertainty in the problem, we consider two moments-based uncertainty sets to define the sets D k , where k = 1, ..., K.

We begin by assuming that we have knowledge of the mean vector µ k and the covariance matrix Σ k associated with the random vector ζT k . Based on this information, we define the following uncertainty ser for each k = 0, 1, ..., K.

D 1 k (µ k , Σ k ) = F k E[ ζT k ] = µ k E[( ζT k -µ k )( ζT k -µ k ) T ] = Σ k
, where F k is a probability distribution of ζT k . In the case where we have the mean vector µ k and the covariance matrix Σ k of ζT k available, we can obtain a deterministic reformulation for the distributionally robust joint chance constraint (6.77) as follows [START_REF] Cheng | Distributionally robust stochastic knapsack problem[END_REF].

                 µ T k x + α z k 1 -α z k ||Σ 1 2 k x||≤ b k , k = 1, ..., K K k=1 z k = 1, z k ≥ 0, k = 1, .., K. (6.79)
The deterministic equivalent problem for (6.76)-(6.77) under these assumptions is given by.

min µ T 0 x, (6.80)
s.t. µ T k x + α z k 1 -α z k ||Σ 1 2 k x||≤ b k , k = 1, ..., K (6.81)
K k=1 z k = 1, x ≥ 0, (6.82)
z k ≥ 0, k = 1, .., K. (6.83) Lemma 60. The function z → α z 1-α z , with 0 < α < 1 is convex ∀z > 0.
Proof. Let z > 0 and 0 < α < 1, we have

α z 1-α z = exp 1 2 (zlog(α) -log(1 -α z ))
. We have z → α z is a convex function and the function z → log(1 -z) is non-increasing and concave, there follows that z → log(1-α z ) is concave. We have then that z →1 2 (zlog(α)-log(1 -α z )) is convex as an addition of two convex functions. Furthermore, z → e z is a non-increasing convex function we conclude then that z → exp 1 2 (zlog(α) -log(1 -α z ))

is convex. The conclusion follows.

Corollary 61. Problem (6.80)-(6.83) is biconvex on (x, z)

Now we consider the case where the mean of ζk lies in an ellipsoid of size γk1 ≥ 0 with center µ k and the covariance matrix of ζk lies in a positive semidefinite cone of center Σ k . We define for every k = 0, 1, ..., K,

D 2 k (µ k , Σ k ) = F k (E F k [ ζT k ] -µ k ) T Σ -1 k (E F k [ ζT k ] -µ k ) ≤ γ k1 COV F k ( ζT k ) ⪯ γ k2 Σ k
, where γ k2 ≥ 0 and COV F k is a covariance operator under probability distribution F k . The deterministic reformulation for the distributionally robust joint chance constraint (6.77) in this case where the mean lies in an ellipsoid and the covariance matrix lies in a positive semidefinite cone is provided in [START_REF] Peng | Games with distributionally robust joint chance constraints[END_REF].

                   µ T k x + α z k 1 -α z k √ γ k2 + √ γ k1 ||Σ 1 2 k x||≤ b k , k = 1, ..., K K k=1 z k = 1, z k ≥ 0, k = 1, .., K. (6.84)
The objective function can be formulated as [START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF].

min x∈ n + µ T 0 x + √ γ 01 ||Σ 1 2 0 x||. (6.85)
The constraints set (6.84) is biconvex and the objective function (6.85) is convex.

To study the optimality conditions of the robust joint chance-constrained problem. We give the equivalent deterministic problem for each uncertainty set in a general form as follows.

min f (x), (6.86)

s.t. g k (x, z) ≤ 0, k = 1, ..., K, (6.87)
h(z) ≤ 0, (6.88)

l(x) ≤ 0, (6.89)
where, f

(x) = µ T 0 x, if D k = D 1 k µ T 0 x + √ γ 01 ||Σ 1 2 0 x||, if D k = D 2 k , h(z) = ( K k=1 z k -1, 1- K k=1 z k , -z 1 , -z 2 , ..., z K ) T , l(x) = -x and g k (x, z) =    µ T k x + α z k 1-α z k ||Σ 1 2 k x||-b k , if D k = D 1 k µ T k x + ( α z k 1-α z k √ γ k2 + √ γ k1 )||Σ Definition 22.
Let U the feasible set of (6.86)-(6.89), U x = {z | g k (x, z) ≤ 0, h(z) ≤ 0, k = 1, ..., K} and U z = {x | g k (x, z) ≤ 0, l(x) ≤ 0, k = 1, ..., K}. (x * , z * ) is a partial optimum of (6.86)-(6.89) if f (x * ) ≤ f (x), ∀x ∈ U z * .

Lemma 62. Let U the feasible set of (6.86)-(6.89) and (x * , z * ) ∈ U. If there exists β (1) , β (2) , γ and λ such that (x * , z * ) verifies ∇ x f (x) + β (1) T ∇ x g(x, z) + λ T ∇ x l(x) = 0, (6.90) λ ≥ 0, λ T l(x) = 0, β (1) ≥ 0, β (1) T g(x, z) = 0, (6.91)

β (2) T ∇ z g(x, z) + γ T ∇ z h(z) = 0, (6.92)
β (2) ≥ 0, β (2) T g(x, z) = 0, γ ≥ 0, γ T h(z) = 0, (6.93)

where g(x, z) = (g 1 (x, z), ..., g K (x, z)) and (x) + = max(0, x). Then (x * , z * ) is a partial KKT point of (6.86)-(6.89).

Theorem 63. The partial KKT system (6.90)-(6.93) is equivalent to the following system ∇ x f (x) + ∇ x g(x, z) T (β (1) + g(x, z)) + + ∇ x l(x)(λ + l(x)) + = 0 (6.94) ∇ z g(x, z) T (β (2) + g(x, z)) + + ∇ z h(z) T (γ + h(z)) + = 0 (6.95) (β (1) + g(x, z)) + -β (1) = 0 Proof. By (β (1) + g(x, z)) + = β (1) and (λ + l(x)) + = λ, we have ∇ x f (x) + ∇ x g(x, z) T (β (1) + g(x, z)) + + ∇ x l(x)(λ + l(x)) + = 0 is equivalent to ∇ x f (x) + β (1) T ∇ x g(x, z) + λ T ∇ x l(x) = 0

We then obtain the equation (6.93) of the partial KKT system. Furthermore, observe that • (β (1) + g(x, z)) + -β (1) = 0 if and only if β (1) ≥ 0, g(x, z) ≤ 0 and β (1) T g(x, z) = 0,

• (λ + l(x)) + -λ = 0 if and only if λ ≥ 0, l(x) ≤ 0 and λ T l(x) = 0, which leads to the equation (6.91) of the partial KKT system. We obtain the remaining equations following the same lines. The converse part of the theorem is straightforward.

The following theorem gives the optimality conditions of problem (6.86)-(6.89) similarly to the previous Chapters. We use wavelet mutation operator from [START_REF] Ling | Hybrid particle swarm optimization with wavelet mutation and its industrial applications[END_REF] and performing for the i th particle if d < ζ as follows y i (j + 1) = y i (j) + µ(h i -y i (j)) , ρ > 0 y i (j) + µ(y i (j) -l i ) , ρ < 0 As the value of ρ approaches 1, the mutated element of the particle will move towards the maximum value of y i (j + 1), whereas approaching -1, the mutated element goes towards the minimum value of x i (j + 1). The magnitude of |ρ| determines the size of the search space for x i (j + 1), with larger values indicating a wider search space. Conversely, smaller values of |µ| lead to a smaller search space for fine-tuning. To achieve fine-tuning, the value of the dilation parameter a is adjusted based on the current number of iterations j relative to the total number of iterations T . Specifically, a is a function of j/T , we take a = e 10 j T . We generate ϕ randomly from [-2.5a, 2.5a]. The algorithm details are given in Algorithm 3 where y = (x, z, β, γ) Algorithm 3 The neurodynamic duplex -Let y 1 (0) and y 2 (0) be randomly generated in the feasible region.

-Let ỹ(0) = ŷ(0) = y(0) the initial best previous position and best position, respectively.

-Set the convergence error ϵ. while ||y(j + 1) -y(j)||≥ ϵ do Compute the equilibrium points ȳ1 (j) and ȳ2 (j) of RNN 1 and RNN 2 based on (6.28)-(4.96).

if f (x 1 (j)) < f (x(j)) then ỹ(j + 1) = ȳ1 (j) else ỹ(j + 1) = ỹ(j) end if if f (x 2 (j)) < f (x(j)) then ỹ(j + 1) = ȳ2 (j) else ỹ(j + 1) = ỹ(j) end if if f (x(j)) < f (x(j)) then ŷ(j + 1) = ỹ(j + 1) else ŷ(j + 1) = ŷ(j)

end if

Compute the value of y(j + 1) following (6.105)-(6.106).

if d < ζ then

Perform the wavelet mutation (6.108). end if j=j+1 end while 6.3.2 . Convergence Analysis Lemma 66. [START_REF] Uryasev | Stochastic Optimization: Algorithms and Applications[END_REF] Suppose that the objective function f is measurable, and the feasible region U is a measurable subset, and for any Borel subset B of U with positive Lebesgue measure we have converges to an equilibrium point, then the state of a neurodynamic model with two timescales, as described by equations (6.100)-(6.104) globally converges to a partial optimum of problem (6.86)-(6.89).

Proof. The Lagrangian function of problem (6.86)-(6.89) is given by L(x, z, β, λ, γ) = f (x) + β T g(x, z) + γ T l(x) + λ T h(z). As λ * T l(x) ≤ λ * T l(x * ) = 0, γ * T h(z) ≤ γ * T h(z * ) = 0 and β * T g(x * , z) ≤ β * T g(x * , z * ) = 0 from the KKT conditions, then f (x * ) ≤ f (x) and this for every z ∈ U *

x . By Definition 22, we have that x * is a partial optimum of (6.86)-(6.89).

Theorem 68. If the state of the two-timescale neurodynamic model (6.100)-(6.104) converges to a partial optimum and the initial states and time constants of the two neurodynamic models are sufficiently different. Then the duplex of two two-timescale neural networks in Figure 6.6 system is globally convergent to a global optimal solution of problem (6.86)-(6.89) with probability one.

Proof. By Theorem 58, the two-timescale neurodynamic models RNN 1 and RNN 2 are proven to converge to a partial optimum. From Algorithm 2, the solution sequence is generated as follows ŷ(j + 1) = ỹ(j + 1) if f (x(j)) < f (x(j)), ŷ(j + 1) = ŷ(j) else.

We observe that the generated solution sequence {f (ŷ(j))} ∞ j=1 is monotonically increasing.

Let M i,j be the supporting set of the initial state of RNN i at iteration j. where Φ is the set of the global optimal solutions of (6.86)-(6.89). The conclusion follows.

. Numerical experiments

To evaluate the performances of our approach, we consider a standard profit maximization problem. A manufacturing firm produces n products with N different machines. The times required to manufacture each unit are random variables. The mean vector µ j and the covariance matrix Σ j describing the uncertainty sets of the time vector t j = {t ij } 1≤i≤n , where t ij is the time required to manufacture one unit of each of product i using machine j and the daily capacity of each machine j given by b j are given. The objective of the study is to determine the daily number of units to be manufactured for each product without exceeding the available 138 machining times. We write our robust joint chance-constrained maximization problem as follows.

min sup x ≥ 0, (6.130) where vector c is a random variable and corresponds to the profit per unit for each product, t ij is the time required to manufacture one unit of product i using machine j, b j is the time capacity of machine j, p is a given probability level, D 0 is an uncertainty set for the distribution F 0 of c and D is an uncertainty set for the distribution F of the random variables.

F 0 ∈D 0 -E cT x ,
All the algorithms in this Section are implemented in Python. We run our algorithms on Intel(R) Core(TM) i7-10610U CPU @ 1.80GHz. The random instances are generated with numpy.random, and we solve the ODE systems with solve_ivp of scipy.integrate. The deterministic equivalent programs are solved with the package gekko and the gradients and partial derivatives are computed with autograd.grad and autograd.jacobian. For the following numerical experiments, the values of µ j and c the mean of c are uniformly generated in [2.0, 4.0], the components of the matrix Σ j are uniformly drawn in the interval x ≥ 0, z j ≥ 0, j = 1, ..., N, (

The neurodynamic duplex vs. convex approximations Cheng et al. [26] propose two convex approximations to solve problem (6.128)- (6.130). A linear approximation that gives an upper bound to the minimization problem and a tangent approximations that leads to a lower bound. In this first subsection, we compare the objective value obtained using the neurodynamic duplex with those obtained using the linear and the tangent approximations. We compute the gap between the two bounds and the global optimum given by the neurodynamic duplex by GAP = Bound lower, upper -ND Bound lower, upper ,

where Bound lower is the value of the lower bound, Bound upper is the value of the upper bound, and ND is the value obtained using the neurodynamic duplex. We recapitulate the obtained results in Table 6.6. Column one gives the value of the confidence parameter p. Column two gives the final value of the neurodynamic duplex. Columns three and four show the lower bound and its gap with the neurodynamic duplex, respectively. Finally, columns five and six present the upper bound and the gap with the neurodynamic approach. We observe that the final value obtained with the dynamical duplex remains between the two bounds for the different values of p with gaps less than 0.5%, demonstrating that the neurodynamic approach effectively converges to the global optimum.

Moreover, we remark that as p increases the value of the objective function increase which is coherent since lower values of p induce larger risk area.

The distributionally robust optimization approach vs. stochastic optimization approaches

To evaluate the robustness of the proposed duplex for the two uncertainty sets D 1 and D 2 , we additionally solve problem (6.128)-(6.130) when the random variables follow uniform and normal distributions and p = 0.95. We compare the solution of our proposed distributionally robust approach with the solution of the stochastic programming approach.

We generate 100 instances for (t ij ) 1≤I≤n,1≤j≤N using the mean vectors and the covariance matrix when the true distribution of the stochastic variables is one of the five following distributions: uniform distribution, normal distribution, log-normal distribution, logistic distribution and Gamma distribution. We calculate the number of times when the constraints were violated over the 100 generated scenarios for each stochastic and robust solutions. Table 6.8 recapitulates the obtained results, where column one gives the true distribution, columns two, three, four and five give the number of violated scenarios for the solution obtained using the uniform approach, the normal approach, the first robust approach and the second robust approach, respectively. The relative expected profit is computed relatively to the value achieved by the solution of the stochastic program with uniform distribution. We observe that the distributionally robust approaches are more conservative compared to the stochastic approaches. We invest between 4.3% and 12.2% of the expected profit in order to ensure the joint constraint. In fact, the average number of violated scenarios for the robust approaches are 0 while the numbers of violated scenarios for the stochastic solutions are significant, i.e., when Gamma is the true distribution of the random variables, the average number of the violated scenarios are 24 and 9 for the uniform and the normal solutions, respectively.

Stochastic solutions

Robust solutions Uniform Normal Dynamical recurrent neural networks have emerged as a powerful tool for tackling optimization problems, showcasing notable advancements in both theory and applications within optimization theory. However, despite these significant strides, several important challenges remain unresolved, creating an appealing avenue for further research and exploration in neurodynamic optimization. This thesis has been dedicated to investigating the practical application of recurrent dynamical neural networks in solving various geometric optimization problems. The dissertation's contributions encompass the development of novel theoretical frameworks, algorithmic techniques, and computational methodologies that effectively incorporate chance-constrained geometric optimization. Additionally, the research conducted in this thesis extends to practical applications across diverse domains, including telecommunications, engineering, transportation, and shape optimization. Through investigating these real-world problem scenarios, valuable insights and solutions have been derived, facilitating informed decision-making processes.

Chapter 4 of this thesis is dedicated to the exploration of geometric programs with joint chance constraints. We focus specifically on the case where the stochastic parameters are independently distributed according to a normal distribution. To address this problem, we begin by deriving a biconvex deterministic reformulation of the initial chanceconstrained problem through the application of a standard variable transformation technique. Next, we review existing state-of-the-art methods that rely on convex approximation. In contrast, we propose a novel neurodynamic approach that solves the stochastic program without the need for any approximations. We establish the stability and global convergence of the proposed neurodynamic approach using Lyapunov analysis. To validate the effectiveness of our approach, we conduct numerical experiments that demonstrate its robustness and superior performance compared to a sequential algorithm. Furthermore, we extend our research to consider the scenario of joint rectangular geometric chance-constrained programs. In this case, we obtain the biconvex equivalent formulation by employing the log transformation and leveraging the arithmetic-geometric mean inequality. In the final part of the chapter, we present a real-world application of our proposed method. Specifically, we apply it to solve a problem involving the maximization of signal-to-interference noise ratio for massive multiple-input, multiple-output (MIMO) systems.

In Chapter 5, we introduced the copula theory as an efficient tool to model chanceconstrained optimization problems with dependent stochastic variables. We first consider the case of linear programs. We propose a neurodynamic approach to solve the equivalent deterministic program based on the optimality conditions given by the partial KKT system. We apply the proposed algorithm to solve a profit maximization problem. Further, we extend our approach to solve geometric-dependent chance-constrained programs. The deterministic equivalent was obtained using an Archimedian copula combined with a log transformation. A generalized shape optimization problem was solved using the neurodynamical approach with different elliptical distributions. Chapter 6 of this thesis focuses on distributionally robust chance-constrained optimization. Initially, we explore the linear case and consider two different uncertainty sets based on moments. For each uncertainty set, we derive the corresponding deterministic equivalent formulation. To tackle the obtained equivalents, we propose a duplex consisting of two two-timescale recurrent neural networks. This novel approach demonstrates convergence almost surely to a global optimum, providing a robust solution method for distributionally robust chance-constrained optimization problems in the linear case. In addition to the linear case, we extend our investigation to geometric programs. We propose two uncertainty sets: one with known first-order moments and the other with a known first-order moment and nonnegative support. We further consider the scenarios of independent and dependent random variables, resulting in four different deterministic programs. For the first three programs, where the obtained equivalents are convex, we utilize a single timescale dynamical neural network to solve them efficiently. However, for the fourth program, biconvex, we adapt the neurodynamic duplex approach to address its complexity and find the optimal solution effectively. By incorporating distributionally robust techniques and leveraging the power of recurrent neural networks, this chapter contributes to advancing the field of chance-constrained optimization. The proposed methodologies offer robustness and efficiency in solving linear and geometric distributionally robust chance-constrained programs. The convergence properties and adaptability of the neurodynamic duplex highlight its effectiveness in finding global optima for complex optimization problems.

In conclusion, the outcomes obtained in this dissertation not only enrich the existing knowledge base but also offer fresh perspectives and practical tools and methodologies for addressing geometric chance-constrained problems in real-world settings. The research presented herein contributes to advancing the field, paving the way for further exploration and application of recurrent dynamical neural networks in optimization contexts.

. Perspectives

In this dissertation, we explore the research landscape of geometric chance-constrained problems and propose novel neurodynamic algorithms to tackle these optimization challenges. While our work contributes valuable insights and algorithmic solutions, there are still open issues and future directions that warrant further investigation in the application of dynamical neural networks to chance-constrained optimization.

One important area for future research is the development of more efficient and scalable neurodynamic algorithms for solving large-scale chance-constrained optimization problems. The size and complexity of real-world problems often pose computational challenges, and exploring techniques to enhance the computational efficiency of neural network-based approaches can greatly improve their applicability and practicality.

Furthermore, the exploration of alternative neural network architectures and learning paradigms could yield further improvements in addressing chance-constrained optimization problems. For instance, investigating recurrent neural networks with attention mechanisms or reinforcement learning-based approaches can offer new perspectives and potentially overcome limitations associated with traditional neurodynamic algorithms.

Moreover, the theoretical analysis of neurodynamic algorithms applied to chanceconstrained optimization problems is an important research direction. Examining the convergence properties, stability guarantees, and robustness of these algorithms under different problem settings and assumptions can provide a deeper understanding of their behavior and facilitate their reliable application in practical scenarios.

Additionally, exploring the combination of dynamical neural networks with other optimization techniques and frameworks, such as evolutionary algorithms, swarm intelligence, or metaheuristics, can lead to hybrid approaches that leverage the strengths of different methods and offer improved performance in solving chance-constrained optimization problems.

Finally, we must highlight that the efficiency and quality of the developed algorithms can be improved through further research and development. One approach is to implement ODE solvers based on artificial intelligence techniques. Neural networks or reinforcement learning can potentially enhance the speed and accuracy of solving the dynamical differential system. These AI techniques can learn patterns and optimize the solution process, leading to faster and more precise results.

Additionally, other optimization techniques, such as parallel computing, GPU acceleration, or distributed computing, can accelerate the algorithm's execution time further. These techniques leverage hardware advancements to process computations in parallel, reducing the overall time required for solving the system.

By exploring these avenues of research and development, the algorithm can be made more efficient, allowing for faster and more effective solutions to chance-constrained geometric optimization problems.

Figure 2 . 1 :

 21 Figure 2.1: Two-bar truss under load [100]

8 .

 8 We finally find y * = 42.426, h * = 30 in and d * = 2.475 in.

x 1 = 3 .

 13 5213, x2 = 2.9089, x3 = 7043, the dual optimal values are δ 1 = 0.75658, δ 2 = 0.17857, δ 3 = 0.06485 δ 4 = 0.35714, δ 5 = 8.427 × 10 -8 and δ 6 = 0.18515.

  ) dy dt = -y + (y + αc(x)) + . (3.8) The dynamical equation represented by equations (3.7)-(3.8) can be implemented using a one-layer recurrent neural network, as depicted in Figure 3.1.

Figure 3 . 1 :

 31 Figure 3.1: Block diagram of neural network (3.7)-(3.8) [142]

Figure 4 . 1 :Theorem 22 .

 4122 Figure 4.1: Feedback interconnection of neural network (4.17)-(4.21)

  )-(4.21) is stable in the Lyapunov sense and converges to (r * , x * , β * , λ * , γ * ), where (r * , x * ) is a KKT point of problem(4.11).

(4. 23 )

 23 We have that dV (y(t)) dt = ( dΦ dt ) T Φ+Φ T dΦ dt +(y -y * ) T dy dt and since dΦ dt = ∇Φ ∇y dy dt = ∇Φ(y)Φ(y), then dV (y(t)) dt = Φ T (∇Φ(y) T +∇Φ(y))Φ+(y-y * ) T Φ(y). By Lemma 23, we have Φ T (∇Φ(y))Φ ≤ 0 and Φ T (∇Φ(y) T )Φ ≤ 0. By Lemma 23 and Lemma 24, we observe that (y -y )T (Φ(y)) = (y -y )T (Φ(y) -Φ(y)) ≤ 0.

  )-(4.21) is convergent in the sense of Lyupanov to an equilibrium point ŷ = (r, x, β, λ, γ) where (r, x) is a KKT point of problem(4.11).

Table 4 . 1 :

 41 Comparison of our neural network vs the deterministic formulation for the 3D-transportation problem Regarding the numerical experiments, we have set the values as follows: ϵ = 0.1, 0.5, 0.01), α = β = γ = δ = 0.5.

× 100 ,

 100 as(Obj valueSA-Obj valueNN) Obj valueSA where Obj valueNN and Obj value SA represent the objective values obtained with the neural network and the sequential algorithm, respectively.

Figure 4 . 4 .

 44

Figure 4 . 3 : 2 :

 432 Figure 4.3: Number of VS for the deterministic formulation, the individual and joint chance constraints

Figure 4 . 4 :

 44 Figure 4.4: Number of VS of our neural network and the sequential algorithm

Lemma 27 .

 27 implementation of the neural network (4.60)-(4.66) is provided in Figure4.6.To study the stability and the convergence of the proposed neural network, we first show the equivalence between a KKT point (4.54)-(4.59) and an equilibrium point of (4.60)-(4.66). Let x and y two n-dimensional vectors, we have (x + y) + = x ⇔ x ≥ 0, y ≤ 0 and x T y = 0 .

  Let y * = (r * , z * , x * , µ * , λ * , γ * , ω * ) an equilibrium point of (4.60)-(4.66) and consider the Lyapunov function defined by V 1 (y) = ||Φ(y)|| 2 + 1 2 ||y -y * || 2 . We have that dV 1 (y) dt ≤ 0. In fact, dV 1 (y) dt = ( dΦ dt ) T Φ + Φ T dΦ dt + (y -y * ) T dy dt . Or since dΦ dt = ∇Φ(y)Φ(y), then dV 1 (y) dt = Φ T (∇Φ(y) T + ∇Φ(y))Φ + (y -y * ) T Φ(y). We use Theorem 29 and Lemma 31 to conclude.

For

  the numerical experiments, we set ϵ = 0.05, 1 A f loor ∼ N (1.0/20.0, 0.01), 1 A wall j ∼ N (1.0/60.0, 0.001), α wall = 0.5, β wall = 1.0, α f loor = 0.5 and β f loor = 1.0.

(4. 77 )

 77 P min ≤ p i ≤ P max , ∀i ∈ U.

Figure 4 . 6 :Figure 4 . 7 :Figure 4 . 8 :Figure 4 . 9 :

 46474849 Figure 4.6: The architecture of the neural network (4.60)-(4.66)

Figure 4 . 10 :

 410 Figure 4.10: Out of 100 scenarios, the constraints were violated 3 times .

Figure 4 . 11 :

 411 Figure 4.11: The Dynamical Neural Network vs. The Alternate Convex Search

Figure 4 . 12 :

 412 Figure 4.12: Signal to Interference plus Noise Ratio, illustration.

Figure 4 . 13 :Figure 4 . 14 :

 413414 Figure 4.13: Convergence of the neural network different starting points y 0

Figure 4 . 15 :

 415 Figure 4.15: Evolution of VS function to ϵ .

Definition 15 .

 15 A copula C : [0, 1] d -→ [0, 1] of dimension d is a joint cumulative distribution function with for which the marginals are uniformly distributed on the interval [0, 1]. Theorem 32. (Sklar's Theorem 1959) Let F a d-dimensional cumulative distribution func-tion with marginals F 1 , F 2 , .., F d , there exists a copula C such that

Definition 18 .Theorem 33 .

 1833 If lim u→+∞ ψ(u) = +∞, then C is called a strict Archimedean copula and ψ is called a strict generator. Let ψ : [0, 1] -→ R + a convex, strictly decreasing function such that lim u→+∞ ψ(u) = +∞ and

Figure 5 . 1 :

 51 Figure 5.1: Quantile functions of selected standard elliptical distributions

Figure 5 . 2 :Theorem 40 .

 5240 Figure 5.2: A simplified implementation of the neural network (5.8)-(5.13)

3 .× 100 ,

 3100 Column one gives the number of variables m. Columns two and three show the objective value and the number of VS obtained through our neural network, respectively. Columns four and five present the objective value and the number of VS obtained using the sequential algorithm. The sixth column gives the gap between the two objective values, calculated as follows: GAP = (Obj valueSA-Obj valueNN) Obj valueSA where Obj valueNN and Obj value SA represent the objective values obtained by the neural network and sequential algorithm, respectively.

( 5 .Figure 5 . 3 :

 553 Figure 5.3: A block diagram for the neural network (5.23)-(5.28)

  and

dE 1

 1 (z(t)) dt ≤ 0 then the neural network (5.23)-(5.28) is globally stable in the sense of Lyapunov.

8 :

 8 Results of problem (5.29)-(5.31) for different values of θ We observe that the neural network converges for the different values of θ and the number of violated scenarios doesn't exceed 6 which proves the robustness of the approach as shown by Figure 5.4.

Figure 5 . 4 :

 54 Figure 5.4: Number of VS for different values of θ

1 ,

 1 k = 1, ...., K, where c k i , i = 1, ..., I k , k = 0, ...., K are positive constants and the exponents a k ij , i = 1, ..., I k , j = 1, ..., M, k = 0, 1, ..., K are real constants.

( 6 . 6 )

 66 Theorem 50.[START_REF] Liu | Distributionally robust chance constrained geometric optimization[END_REF] If σ k i,l ≥ 0 for all i, l and k, problem (JCP log ind ) is a convex programming problem.

  ) λk ≤ πk , k = 1, ...., K, r j -βk i ≤ 0, i = 1, ..., I k , k = 1, ...., K.(6.22)

logNSy

  ) into a convex problem. Nevertheless, we can obtain the following biconvex equivalent by takingt j = log(r j ), λk = log(λ k ), βk = log(-β k ), π = log(π) k ≥ K -ϵ, 0 ≤ y k ≤ 1 , k = 1, ..., K, (6.30) y k exp(-λk ) + I k i=1 expλk + βk i + logµ k i ≤ 1, k = 1, ...., K, λk ≤ 0, k = 1, ...., K,(6.31)λk ≤ πk , k = 1, ...., K, r j -βk i ≤ 0, i = 1, ..., I k , k = 1, ...., K. nd ) can be written in the following general form min r f (z), (6.34) s.t. g(z) ≤ 0, where f and g are two convex functions.

( 6 . 36 )

 636 KKT) are satisfied. ∇f (z) + ∇g(z) T γ = 0 (6.35) γ ≥ 0, γ T g(z) = 0 To solve problem (6.34), we propose a dynamical recurrent neural network driven by the following ODE system κ dz dt = -∇f (z) + ∇g(z) T (γ + g(z)) + (6.37) κ dγ dt = -γ + (γ + g(z)) + (6.38)

(6. 45 ) 21 .

 4521 Definition Let (z, y) ∈ U, (z, y) is called a partial optimum of (6.42) if and only if

Figure 6

 6 Figure 6.6.

Figure 6 . 1 :

 61 Figure 6.1: A block diagram depicting a duplex neurodynamic system with a two-timescale configuration

Figure 6 . 2 : 1 :

 621 Figure 6.2: 3D-box shape [1]

Figure 6 . 3 :

 63 Figure 6.3: Transient behaviors of the state variables

Figure 6 . 4 :

 64 Figure 6.4: Convergence of the dynamical neural network (6.37)-(6.38) for different initial points for (6.70).

Figure 6 . 5 :

 65 Figure 6.5: Convergence of the power variables

(6. 96 )( 6 . 97 )( 6 . 98 )

 96697698 (β(2) + g(x, z)) + -β (2) = 0 (λ + l(x)) + -λ = 0 (γ + h(z)) + -γ = 0 (6.[START_REF] Miller | Chance constrained programming with joint constraints[END_REF]

Figure 6 . 6 :

 66 Figure 6.6: A block diagram depicting a duplex neurodynamic system with a two-timescale configuration

(6. 108 )

 108 where h i and l i are the upper and the lower bound for y i , respectively. ζ > 0 is a given threshold and ρ is defined using a wavelet function

∞ k=1 ( 1 -Theorem 67 .

 k=1167 P k (B)) = 0. Let {y(k)} ∞ k=1 be a sequence generated by a stochastic optimization algorithm. If {f (y(k))} ∞ k=1 is a nonincreasing sequence, then it converges with probability one to the set of global optimal solutions. If the state of the following neurodynamic model with a single timescale,κ dx dt = -x + x -(∇ x f (x) + β T ∇ x g(x, z)) + , + (z -(β T ∇ z g(x, z) + γ T ∇ z h(z))) + ,(6.111)κ dβ dt = -β + (β + g(x, z)) + ,(6.112)κ dλ dt = -λ + (λ + l(x)) + . (6.113) κ dγ dt = -γ + (γ + h(z)) + . (6.114)

(6. 115 )

 115 As an equilibrium point of (6.110)-(6.114) corresponds to a KKT point (x * , z * , β * , λ * , γ * ) of (6.86)-(6.89)[START_REF] Xia | A recurrent neural network for nonlinear convex optimization subject to nonlinear inequality constraints[END_REF] verifying∇ x L(x * , z * , β * , λ * , γ * ) = 0, (6.116) ∇ z L(x * , z * , β * , λ * , γ * ) = 0, (6.117) ∇ β L(x * , z * , β * , λ * , γ * ) = 0, (6.118) ∇ γ L(x * , z * , β * , λ * , γ * ) = 0, (6.119) ∇ λ L(x * , z * , β * , λ * , γ * ) = 0, (6.120) β * T g(x * , z * ) = 0, β * ≥ 0, (6.121) γ * T l(x * ) = 0, γ * ≥ 0. (6.122) λ * T h(z * ) = 0, λ * ≥ 0. (6.123)We fix x * and take z ∈ U *x , problem (6.86)-(6.89) becomes convex, we have then L(x * , z * , β * , λ * , γ * ) ≤ L(x * , z, β * , λ * , γ * ),(6.124) which leads to f (x * ) + β * T g(x * , z * ) + λ * T h(z * ) ≤ f (x * ) + β * T g(x * , z) + λ * T h(z).

  As indicated by equation(6.108), the mutation operation ensures that the initial states of the RNNs are forced to be in the feasible region U. Hence, for every iteration index J ≥ 1, the supporting sets fulfill the following conditionWe have then v(U) = v(M) > 0. By Lemma 66, we have lim j->∞ P(ŷ(j) ∈ Φ) = 1(6.127)

t

  ij x i ≤ b j , j = 1, ..., N ≥ p,(6.129)

2 jxp z j 1 -p z j √ γ k2 + √ γ k1 ||Σ 1 2 j

 212 [1.0, 3.0] and we generate the values of b j uniformly in [50.0, 60.0], γ k1 = 5 and γ k2 = 5. The resulting deterministic equivalent problems of (6.128)-(6.130), where the uncertainty sets are D 1 and D 2 are given respectively by min cT x, (6.131)s.t. µ T j x + p z j 1 -p z j ||Σ1 x||≤ b j , j = 1, ..., N, ≥ 0, z j ≥ 0, j = 1, ..., N, x||≤ b j , j = 1, ..., N,

D 1 8 : 141 7

 18141 (µ, Σ) D 2 (µ, Σ) Number of violated scenarios for the stochastic and the robust solutions

Theorem 8 .

 8 Let y * = (x * , u * , v * ) T be the equilibrium of the neural network (3.16)-(3.17). Then, x * is a KKT point of problem (CNLP). If ∇f (x) is positive definite and ∇g k (x), k = 1, ..., m is positive semi-definite, or if ∇f (x) is positive semi-definite and ∇g k (x), k = 1, ..., m is positive definite, then the neural network (3.16)-(3.17) is globally stable in the Lyapunov sense. Additionally, it

	Lemma 9. The equilibrium point of (3.16)-(3.17) is unique.
	Lemma 10. For any initial point y(t 0 ) = (x(t 0 ), u(t 0 ), v(t 0 )) T , there exists a unique contin-
	uous solution y(t) = (x(t), u(t), v(t)) T for the system described by equations (3.16)-(3.17).
	Theorem 11.

Table 4 . 3 :

 43 Comparison between the deterministic approach and the individual and joint chance constraint neural networks

	Data Deterministic Approach	Individual chance constraints	Joint chance constraints
		Obj value VS	Obj value VS GAP	Obj value VS GAP
	(2,4)	0.020	63	0.026	17 23	0.028	10 30
	(3,8)	0.062	88	0.080	24 29	0.088	11 41
	(4,14) 0.111	97	0.130	38 28	0.144	5	41
	(5,22) 0.209	100	0.254	54 21	0.273	4	30
	(6,32) 0.355	99	0.400	45 12	0.421	13 18
	(7,44) 0.520	100	0.607	59 16	0.648	12 24

Table 4 .

 4 4: Neural network vs. the sequential algorithm algorithm, column six shows the corresponding number of VS, and column seven shows the gap between the two objective values. The gaps range from 8.9% up to 17.8%. Notice that the number of VS is comparable for the two approaches.

  .5. Column one gives the number of variables m. Columns two, three, and four give the objective value of problem (4.69), the number of VS, and the corresponding CPU time in seconds, respectively. Columns five, six, and seven give the objective value of problem (4.68), the number of VS, and the corresponding CPU time, respectively. We observe that the objective values of the two problems are relatively close. Nevertheless, the problem (4.68) covers better the risk region. In fact, we remark that the number of violated scenarios for m = 5 for problem (4.68) is equal to 3, whereas the number of violated scenarios for problem (4.69) is equal to 6 as shown in Figures 4.9 and 4.10.

Table 4 .

 4

		5: Results of the generalized transportation problem for different values
	of m
	m	Individual constraints	Joint constraints
		Obj value VS CPU Time	Obj value VS CPU Time
	3	0.039	3	35.27	0.042	0	25.82
	5	0.117	6	61.60	0.120	3	92.20
	7	0.230	4	86.36	0.236	2	80.32
	10	0.440	3	203.84	0.456	1	169.57
	15	0.909	5	444.26	0.907	5	531.83
	20	1.384	4	1111.72	1.402	2	874.75
	Table 4.6: Results of the generalized transportation problem for different values
	of ϵ
	ϵ	Individual constraints	Joint constraints
		Obj value VS CPU Time	Obj value VS CPU Time
	0.1	0.115	8	39.34	0.121	1	66.67
	0.15	0.113	7	39.44	0.115	5	38.02
	0.2	0.112	28 74.80	0.114	15 44.32
	0.3	0.109	31 49.87	0.111	24 49.84
	0.4	0.107	48 44.21	0.109	39 52

Table 4 .

 4 7: Individual constraints vs. Joint constraints for different values of K

	(4.97)

.

[START_REF] Luedtke | A branch-and-cut decomposition algorithm for solving chanceconstrained mathematical programs with finite support[END_REF]

)

dζ dt = (ζ + h(x)) + -ζ.

Table 4 .

 4 8: Neural network vs. the sequential algorithm for different values of K

	K	Sequential Algorithm	VS	Neural Network	VS	GAP
	2	5.40	12	5.10	6	5.88
	3	25.77	10	25.61	8	0.62
	5	28.97	11	28.88	9	0.31
	7	68.79	10	68.56	8	0.33
	10	70.81	21	69.68	14	1.62
	15	84.43	7	84.39	6	0.04
	20	117.37	13	117.33	10	0.03

Table 4 . 9

 49

	ϵ	Sequential Algorithm	VS	Neural Network	VS	GAP
	0.05	30.23	4	29.89	2	1.13
	0.1	29.47	15	29.07	9	1.37
	0.15	28.96	22	28.53	11	1.50
	0.2	28.56	32	28.10	19	1.63
	0.3	27.87	54	27.40	26	1.71
	0.4	27.30	63	26.81	34	1.82

: Neural network vs. the sequential algorithm for different values of ϵ

Table 5 .

 5 3 indicates that our neural network outperforms the sequential algorithm, as the upper bounds obtained are better and fewer scenarios are violated.

	m Neural Network VS Sequential Algorithm VS GAP
	3		0.283	5	0.334		8 15.3
	7		0.399	0	0.444		9	10.1
	10		0.601	1	0.667		9	9.9
	15		0.253	3	0.277		6	8.6
	20		0.663	4	0.685		5	3.2
	Table 5.3: Neural network vs. the sequential algorithm for normal distribution
	m		θ = 1		θ = 2			θ = 10
		(independent constraints)
		Obj value VS		Obj value VS	Obj value VS
	10	2.78	11		1.97	14	1.06	21
	15	13.91	11		10.46	22	6.19	31
	20	1.07	8		0.85	9	0.60	10
	30	3.13	1		2.50	7	1.75	9

We solve

(5.14)

for different values of m and different values of θ; the dependence parameter of Gumbel-Hougaard Copula. Tables 5.4, 5.5, 5.6 show the obtained results

Table 5 .

 5

4: Results for different values of m for normal distribution

Table 5 .

 5

		independent constraints)
		Obj value VS	Obj value VS	Obj value VS
	10	2.30	24	1.34	30	0.86	36
	15	4.22	19	2.47	26	1.53	31
	20	5.19	10	3.79	15	2.33	18
	30	11.49	2	6.83	5	4.21	8
	Table 5.5: Results for different values of m for Laplace distribution
	m		θ = 1	θ = 2		θ = 10
		(independent constraints)
		Obj value VS	Obj value VS	Obj value VS
	10	0.75	1	0.28	4	0.21	7
	15	2.21	4	1.35	5	0.75	6
	20	3.30	5	2.37	7	1.38	9
	30	4.98	1	3.56	5	2.05	9

6: Results for different values of m for Logistic distribution matrix Σ k .

  1 -65x 2 -62x 3 , P(t 11 x 1 + t 12 x 2 + t 13 x 3 ≤ 680, t 21 x 1 + t 22 x 2 + t 23 x 3 ≤ 660, t 31 x 1 + t 32 x 2 + t 33 x 3 ≤ 710) ≥ 0.95,

	Machine			Time per unit (min)			Machine
		Product 1	Product 2	Product 3	capacity
		mean	standard	mean	standard	mean	standard	(min/day)
			deviation		deviation		deviation
	Machine 1	t 11 = 12 σ 11 = 2.3	t 12 = 12 σ 12 = 2.0	t 13 = 12 σ 13 = 2.0	C 1 = 680
	Machine 2	t 21 = 12 σ 21 = 2.9	t 22 = 12 σ 22 = 2.3	t 23 = 13 σ 23 = 2.0	C 2 = 660
	Machine 3	t 31 = 14 σ 31 = 2.9	t 32 = 14 σ 32 = 2.9	t 33 = 12 σ 33 = 2.9	C 3 = 720
				Table 5.7: Data
	The resulting deterministic equivalent problem of (5.29)-(5.31) is given as follows
		min -60x 1 -65x 2 -62x 3 ,
		s.t. t11 x 1 + t12 x 2 + t13 x 3 + ϕ -1 ((0.95) y 1 1 θ ) σ 2 11 x 2 1 + σ 2 12 x 2 2 + σ 2 13 x 2 3 ≤ 680,
								(5.29)
			s.t. (5.30)
				x 1 , x 2 , x 3 ≥ 0.				(5.31)

Table 5 .

 5 10: Results for different data sizes with α = 0.10

	α		θ = 1				θ = 2	θ = 10
		(independent constraints)
		Obj value VS CPU Time	Obj value VS CPU Time	Obj value VS CPU Time
	0.10	-2803.75 1	9.50		-2856.78 1	9.46	-2897.36 4	2.60
	0.15	-2863.54 2	9.89		-2926.52	5	9.19	-2976.50 8	2.99
	0.20	-2911.68	2	9.42		-2983.73 9	9.70	-3041.72	11 3.22
	0.25	-2953.43 9	10.86		-3034.27 15 10.33	-3100.25	24 3.83
				Table 5.9: Results for different values of α
	n k			θ = 2				θ = 10
		Tangent		VS Neural	VS GAP		Tangent	VS Neural	VS GAP
		approximation		network				approximation	network
	5 3 -2422.18		0	-2420.41 0	0.07%	-2425.46	1	-2425.16 1	0.01%
	7 5 -2524.58		5	-2523.07 5	0.05%	-2545.96	10 -2545.78 10 0.007%
	10 5 -2638.22		3	-2634.13 2	0.15%	-2675.62	9	-2674.92 9	0.02%
	15 10 -2839.43		4	-2837.83 4	0.05%	-2867.28	7	-2867.05 7	0.08%

  .[START_REF] Duffin | Cost minimization problems treated by geometric means[END_REF]) is locally Lipschitz continuous, given that ∇f , ∇g and (γ + g) + and are locally Lipschitz continuous. By applying the local existence theorem of ordinary differential equations, we can conclude that there exists a unique continuous solution trajectory (z(t), γ(t)) for (6.37)-(6.38).

	Theorem 54. The neural network proposed in equations (6.37)-(6.38) exhibits global sta-
	bility in the Lyapunov sense. Furthermore, the dynamical network globally converges to
	a KKT point denoted (z

* , γ * ) where z * is the optimal solution of the problem (6.34).

Table 6 .

 6

	m	Independent case CPU Time	Dependent case CPU Time
	3	1.03	1.05	1.30	1.39
	5	2.09	5.11	2.15	5.20
	10	14.79	4.83	5.04	15.10
	15	7.76	47.80	7.99	58.04
	20	10.68	97.72	10.87	100.91

2: Results for different values of m

Table 6 .

 6 3: (6.70) for different values of m for D = D 3more robust, i.e., the number of VS when the true distribution is the Logistic distribution is equal to 23 and 19 for the nonnegative support solutions and is equal to 0 for the robust solutions.

		m	Independent case	Dependent case
			Obj value CPU Time VS	Obj value CPU Time VS
		3	0.204	2.28	3	0.491	10.12	0
		5	1.03	6.25	2	1.82	98.68	0
		10	6.99	15.26	2	9.79	86.35	0
		15	18.43	23.84	3	23.45	201.13	0
		20	32.09	94.76	5	38.71	744.26	0
		30	42.37	100.23	3	51.56	1155.42	0
					Normal solutions		Robust solutions
					Independent Dependent	Independent Dependent
		Objective Value		0.86	0.99		2.43	4.14
	Number	Uniform distribution	22	15		0	0
	of	Normal distribution		18	11		1	0
	violated	Log-normal distribution	7	4		2	1
	scenarios Logistic distribution		23	19		0	0
		Gamma distribution		16	12		2	2

Table 6 .

 6 4: Number of violated scenarios for the stochastic and the robust solutions

Table 6 . 5

 65

	K	Individual constraints	Joint constraints
		Obj value VS	Obj value VS
	5	27.27	5	29.07	0
	10	47.36	4	50.23	0
	15	66.03	5	68.76	1
	20	123.48	3	127.43	0

: Results for different values of K

Table 6 .

 6 6: Results for different values of p for D 1 (µ, Σ)

	p		Neurodynamic duplex	Tangent approximation	Linear approximation
			Obj value	Obj value GAP	Obj value GAP
	0.95		-36.25	-36.41	0.43%	-36.20	-0.13%
	0.9		-40.48	-40.51	0.07%	-40.46	-0.04%
	0.8		-45.30	-45.41	0.24%	-45.22	-0.17%
	0.7		-47.31	-47.38	0.14%	-47.28	-0.06%
	0.6		-48.09	-48.13	0.08%	-48.07	-0.06%
	n N	Neurodynamic duplex	Tangent approximation	Linear approximation
			Obj value	Obj value GAP	Obj value GAP
	7	4	-22.86	-22.97	0.47%	-22.77	-0.39%
	10 5	-22.51	-22.65	0.61%	-22.44	-0.31%
	15 10	-21.36	-21.61	1.15%	-20.82	-2.59%
	20 15	-21.28	-21.78	2.29%	-20.93	-1.67%
	25 20	-19.79	-20.78	4.67%	-19.01	-4.10%

Table 6 .

 6 7: Results for different values of n and N for D 1 (µ, Σ)

Cette thèse vise à explorer l'application des approches neurodynamiques dans le domaine de l'optimisation géométrique sous contraintes de probabilité. La thèse se concentre sur le développement d'algorithmes efficaces et de méthodes de solution pour résoudre des problèmes d'optimisation complexes avec des contraintes géométriques dans des conditions incertaines. La thèse se compose de sept chapitres, comprenant une revue de littérature complète, un contexte théorique et diverses études de cas.Le chapitre 2 propose une revue de littérature sur l'optimisation géométrique. Le chapitre discute des approches de solution et des algorithmes pour résoudre des programmes géométriques, y compris les méthodes duales et primales. Divers algorithmes et techniques ont été passés en revue, tels que la condensation, les méthodes à points intérieurs et les algorithmes primaux-duaux. Le chapitre explore également l'application de la programmation géométrique dans différents domaines, notamment l'ingénierie mécanique, l'ingénierie chimique, le contrôle de puissance et la finance. De plus, des exten-

Lemma 57.[START_REF] Uryasev | Stochastic Optimization: Algorithms and Applications[END_REF] Suppose that the objective function f is measurable, and the feasible region U is a measurable subset, and for any Borel subset B of U with positive Lebesgue

k x||-b k , if D k = D 2 k .

Proof. Let (z, y, ω) is an equilibrium point of (6.50)- (6.52). We have then dz dt = 0 ⇔ -∇f (z) + ∇ z g(z, y) T (ω + g(x, z)) + = 0, (6.53) dy dt = 0 ⇔ -∇ y g(z, y) T (ω + g(z, y)) + = 0, (6.54) dω dt = 0 ⇔ -ω + (λ + g(z, y)) + = 0.

(6.55)

We obtain then system (6.47)- (6.49). By Theorem 55, the conclusion follows. The converse part of the Theorem is traightforward.

The process begins by initializing the state variables of the neurodynamic models. Subsequently, each model undergoes a precise local search based on its dynamics to optimize its performance. Once all neurodynamic models have converged to their equilibria, the initial states of the recurrent neural networks are optimized using the particle swarm optimization (PSO) updating rule. In this context, we represent the position of the i th particle as Λ i = (Λ i1 , ..., Λ in ) T , and its velocity as v i = (v i1 , ..., v in ) T . The inertia weight w ∈ [0, 1] determines the extent to which the particle retains its previous velocity. The best previous position that yielded the maximum fitness value for the i th particle is denoted as Λi = ( Λi1 , ..., Λin ) T , and the best position in the entire swarm that yielded the maximum fitness value is represented by Λ = ( Λ1 , ..., Λn ) T . The initial state of each neurodynamic model is updated using the PSO updating rule, as described in reference [START_REF] Clerc | The particle swarm -explosion, stability, and convergence in a multidimensional complex space[END_REF].

where the iterative index is represented by j, while the two weighting parameters are denoted as c 1 and c 2 and r 1 and r 2 represent two random values from the interval [0, 1].

To achieve global optimization, the diversity of initial neuronal states is crucial. One approach to enhance this diversity is by introducing a mutation operator, which generates a random Λ i (j + 1). This random generation of Λ i (j + 1) helps increase the variation among the initial neuronal states. To measure the diversity of these states, we employ the following function

We utilize the wavelet mutation operator proposed in [START_REF] Ling | Hybrid particle swarm optimization with wavelet mutation and its industrial applications[END_REF], which is performed for the i-th particle if d < ζ. The mutation operation is carried out as follows

where h i and l i are the upper and the lower bound for Λ i , respectively. ζ > 0 is a given threshold and µ is defined using a walvet function (2) then (x * , z * ) is a KKT point of (6.86)-(6.89).

. A neurodynamic duplex

Based on the projection equations (6.47)-(6.49), we propose a duplex of two two-timescale recurrent neural network models for solving (6.86)-(6.89)

)

)

(6.103)

where (x, z, β, γ, λ) are now time-dependent variables and κ 1 and κ 2 are two time scaling constants with κ 1 ̸ = κ 2 . We consider duplex of two two-time-scale recurrent neural network (6.100)-(6.104) for solving (6.86)-(6.89) one with κ 1 > κ 2 and the second with κ 1 < κ 2 as shown in Figure 6.6. The zoom on RNN1 shows the circuit implementation of a single two-timescale recurrent neural network (6.100)-(

Theorem 65. (x * , z * , β * , γ * , λ * ) is an equilibrium point of (6.100)-(6.104) if and only if (x * , z * ) is a KKT point of (6.86)-(6.89) and β * , γ * and λ * are the associated Lagrange variables.

Proof. (x * , z * , β * , γ * , λ * ) is an equilibrium point of (6.100)-(6.104) if and only if dx dt = 0, dz dt = 0, dβ dt = 0, dλ dt = 0 and dγ dt = 0, we then obtain equations (6.94)- (6.99). By Theorems 63 and 64.

First, the state variables of the neurodynamic models are initialized. Then, each model undergoes a precise local search based on its dynamics for the optimization process. Once all neurodynamic models have converged to their equilibria, the initial states of the recurrent neural networks are optimized using the updating rule of particle swarm optimization (PSO). We denote y i = (y i1 , ..., y in ) T the position of the i th particle and v i = (v i1 , ..., v in ) T its velocity. The inertia weight w ∈ [0, 1] determines the degree to which the particle's previous velocity is retained. The best previous position yielding the maximum fitness value for the i th particle is denoted as ỹi = (ỹ i1 , ..., ỹin ) T , and the best position yielding the maximum fitness value in the swarm is represented by ŷ = (ŷ 1 , ..., ŷn ) T .

The initial state of each neurodynamic model is updated using the PSO updating rule given by [START_REF] Clerc | The particle swarm -explosion, stability, and convergence in a multidimensional complex space[END_REF], i.e, v i (j + 1) = wv i (j) + c 1 r 1 (ỹ i -y i (j)) + c 2 r 2 (ŷ i -y i (j)), (6.105) y i (j + 1) = y i (j) + v i (j + 1). (6.106)