Article AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 Contents lists available at http://qu.edu.iq Al-Qadisiyah Journal for Engineering Sciences Journal homepage: http://qu.edu.iq/journaleng/index.php/JQES * Corresponding author. E-mail address: jboscosamon@gmail.com (Jean Samon) https://doi.org/10.30772/qjes.v15i3.835 2411-7773/© 2022 University of Al-Qadisiyah. All rights reserved. This work is licensed under a Creative Commons Attribution 4.0 International License. Analysis of coupling methods: a review Jean Bosco Samon*, Jemimah Maurelle Ndomou, and Damasse Harold Fotsa Laboratory of Mechanics, Materials and Photonics, Department of Mechanical Engineering, ENSAI-University of Ngaoundere, Bp 455, Cameroon A R T I C L E I N F O Article history: Received 15 Jone 2022 Received in revised form 20 August 2022 Accepted 25 September 2022 Keywords: Multiphysics Multiphysics system Mechatronics Mechatronics system Coupling methods A B S T R A C T The multiphysics field is a branch of physics whose objective is to the couple at least two physical systems. Each is governed by its own principles of evolution or equilibrium such as balance laws or constitutive laws. Many engineering problems can only be described correctly by coupling fields of physics that have historically been developed and taught separately. These problems require on the one hand a good understanding of each physical domain, but above all an analysis of the coupling mechanisms of these physical domains, in order to propose a relevant model capable of solving the problem. A challenge in the multiphysics (mechatronics) field is the construction of coupled multiphysics models from experimental observations, as well as the analysis of their mathematical properties. The mathematical analysis of the coupled model must be able to show the well-posedness of the problem at the defined boundary and initial values. For this reason, we have identified several coupling methods: Newton, Gauss-Seidel, JNFK, and direct and explicit coupling. From these methods, it appears that the Newton method is suitable for the coupling of the different disciplines of Mechatronics. A summary table shows the comparative advantages of each method. © 2022 University of Al-Qadisiyah. All rights reserved. 1. Introduction The compartmentalization of trades and physical disciplines, combined with the weakness of computer resources with the emergence of numerical calculation has led to the creation of specialized simulation tools. However, the phenomena studied often mix several disciplines. The need for Multiphysics simulations thus appeared very early on [20]. The field of Multiphysics is a branch of physics whose objective is to couple two or more physical systems, each governed by its own principles of evolution or equilibrium, such as balance laws or constitutive laws [1]. Many engineering problems can only be described correctly by coupling fields of physics that have historically been developed and taught separately [21]. This is called coupling. There are two traditional responses to this need. The first is to incorporate simplified models (which can be fed by reference calculations) from other physics disciplines into a specialized code. The second is the sequential execution of several codes exchanging information. Both approaches have in common the implementation of penalty levers to guarantee the conservatism of the physical response obtained by the modelling. The increase in computing capacity and the desire to improve the modelling associated with the complexity of the technologies studied are increasing the interest in Multiphysics simulations. These objectives are recurrent, but it is not at all certain that they will be achieved if attention is not paid to the way the coupling is carried out. Having separate simulations communicate may introduce limitations in stability, accuracy and robustness that are stronger than those of the components. In addition, the exchange of information may itself be more costly than the execution of the coupled codes. These challenges have led the scientific community to develop efficient and robust coupling techniques. This article is part of this http://qu.edu.iq/ https://doi.org/10.30772/qjes.v15i3.835 http://creativecommons.org/licenses/by/4.0/ JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 193 process. After defining the coupling system, we will present two families of associated coupling methods. We have deduced a method suitable for mechatronic systems. 2. Some definitions of the concepts According to [2] coupling is an approach that can be used as a unifying basis for modeling the different domains of mechatronics. We can therefore say that coupling is the unambiguous connection between two or more inputs of the same system. It is then for us to define a coupled system. According to [23] a system is a set of entities that interact and evolve over time. The entities can be continuous or discrete fields. Their evolution is governed by ODEs1 or PDEs2. He thus defines a coupled system as a system made up of several interacting subsystems. They can be differentiated by the model (conservation law, behavioral, etc.) or the numerical methods used. These subsystems interact through interfaces. Geometrically, an interface can correspond to a boundary or to a part or the whole spatial domain. As a sketch, we have fig. 1 below. Figure 1 - Arrangement of the subsystems [2] Figure 2 - Family of couplings [3] This figure shows us that the subsystems can be adjacent (fluid-structure, solid body dynamics), or overlap each other (overlapped). Partial overlaps are also possible (e.g. electromagnetic-mechanical). In their paper entitled "Stabilization of staggered solution procedures for fluid-structure interaction analysis", [3] discussed the difficulty of stability in the coupling method. It is this difficulty that has challenged the researchers to the point of leading them to a solution to the problem, namely the different coupling algorithms. The latter are grouped into two families qualified respectively as weak coupling and strong coupling which we shall begin by presenting. Let us note that the type of coupling between two subsystems (each physical for example) can be qualified as "Strong" or "Weak" as shown in figure 2. Concerning the "strong coupling": there is a significant influence of the two subsystems on each other. For "weak coupling": only one subsystem has a significant influence on the other subsystem. The fundamental difference between these two families comes from the way of imposing equilibrium between the different physical systems, i.e. the action-reaction principle. It can also be compared in terms of their speed of convergence. This is measured by the coupling force. The coupling "strength" is determined by the equations of the mathematical model. The value of the coupling "strength" between the subsystems is given by a scaled writing of the model equations which allows some scaled numbers to appear. 2.1 Weak coupling family This is the first type of coupling to be used in design offices. It is very present in industrial software. This type has the advantage of concentrating the calculation effort on the zones that require it. However, it is mainly reserved for cases where the local detail has little or no influence on the rest of the structure on a global scale. For this reason, a mathematical modeling of weak coupling will be presented in the following. Consider a dynamic system whose state is described by the vector u(t)={u_1 (t),u_2 (t)}and governed by the following differential system: { u̇ = ℱ(u, t); u(t = 0) = u0. (1) The coupling between the two subsystems described by their respective entities u1(t) et u2(t) or components ℱ1 𝑒𝑡 ℱ2 are both functions of u1(t) et u2(t)The coupling between the two subsystems, described by their respective entities or components, is said to be weak if only subsystem 1 has a significant influence on the second one. ℱ1(u1, u2, t) ≈ ℱ1(u1, t). (2) That is: { u̇1 u̇2 } = { ℱ1 (u1, t) ℱ2(u1, u2, t) } (3) In the weak coupling model, each of the partial differential equations is solved separately. In the weak coupling model, each of the partial differential equations is solved separately, which requires the transfer of the results of one partial differential equation to the other. In such a case, each system may involve a different mathematical type of variable (real or complex), different time constants or pulses, and different study areas. 2.2 Strong coupling family In the case where the influence of the local model extends beyond its definition domain, it is easy to understand that weak coupling methods will produce a more or less erroneous solution depending on the situation [4-5]. Therefore, the only way to build an exact coupled solution is to bring the information from the solution of the local problem to the global level. This is the principle of strong coupling methods. In this type of strong coupling, we have direct and iterative coupling methods. 2.2.1 Direct coupling methods 194 JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 These methods take their name from the fact that the coupling between two models is considered by means of a direct calculation, i.e., without using iterative processes. They can be classified into two main categories: structural reanalysis methods and static condensation structural reanalysis methods. The principle of structural reanalysis is to calculate the response of a structure after a (often localized) change from the principle of structural reanalysis is to calculate the response of a structure after a (often localized) modification from the knowledge of the response before modification. Most of the time, the procedure consists in evaluating the response of the initial system under different elementary stresses (e.g., m) representative of the modification made to the structure. This is then done by correcting the initial stiffness operator by a matrix of rank m. This correction can then be carried out in a non-intrusive manner with respect to the operator through the Sherman-Morrison and Woodbury formulae Woodbury formulae [8]. Let {Ui}i=1...m be a family of vectors constructed from the responses under elementary stresses. If U = [U1, ..., Um] then the response of the modified structure can be calculated in a non-intrusive way (i.e., without requiring the assembly and factorization of the modified system): (K + UUT )-1 = K-1 - K-1 U (UT K-1 U + I) -1 UT K-1 (4) Structural reanalysis methods are therefore based on an a posteriori modification of the solution, in contrast to static condensation methods. We take the example (Fig.3) of a structure occupying the domain Ω that is partitioned into a global part ΩG and a local part ΩL Figure 1 - Partition of the study area (Park and Felippa, 2000). It is assumed that a localized defect in ΩL is taken into account. Even if the solution in ΩG can be affected by the alteration of the model in the local area, the finite element representation of the model defined in ΩG always remains the same. If the model defined at ΩL were to evolve or be modified many times (for example during a parametric study), a complete recalculation of the whole structure would require a significant and unnecessary additional cost. Static condensation sub structuring provides a solution to this problem by allowing the construction of a model representing the behavior of the structure in ΩG involving only the degrees of freedom of the interface Γ [6-7]. We consider here the case of a displacement connection under the assumption that the mechanical problem is formulated in primal variable. The global model KG UG + C UL = FG can be condensed so that only the interface displacement has to be considered. The superscripts i indicate the internal degrees of freedom of the models, and the superscripts b indicate the degrees of freedom of the interface. The Schur complement of KG, noted can thus be used to substitute the global model in equation (5). This gives the following equation 6: In this case, static condensation is used as a decoupling method to connect a fixed global model to a detailed and possibly modifiable local model. Sometimes, this method is also used to provide an indelible finite element model in the form of a "black box" to protect an industrial secret for example. Indeed, the only knowledge of the Schur complement allows an exact consideration of the model, without having access to the details of the latter (geometry, loading cases). In the context of a generic and succinct presentation of the most used coupling methods to date, we refer here to static condensation. However, it should be noted that the process of sub structuring finite element models is more general and can also be used in dynamics and modal calculation. Also, nothing prohibits the joint use of static condensation with other coupling methods, as is the case in Hirai who combines static condensation with structural reanalysis [9]. 2.2.2 Iterative coupling methods As we have just seen, direct coupling methods imply a certain level of interoperability of the solvers used by each of the coupled models. Interoperability of the solvers used by each of the coupled models. The other major family of coupling methods commonly used are the iterative methods. Successive resolutions are performed on the different parts of the model that are to be coupled. Among the common iterative coupling methods are those based on sub structuring [22]. As presented earlier, the structure can be partitioned into substructures. The idea is then to solve the interface problem, but this time in an iterative way (for example via a conjugate gradient method). The other major class of iterative coupling methods concerns the so-called Schwarz methods [8]. The alternating Schwarz method, historically the oldest, proposes to split the study domain into several domains. The problem is then solved alternately in each of the sub-domains. If we assume a division of the domain into two parts Ω1 and Ω2, then the solution of the problem according to the alternating Schwarz method will be done as follows (Equation 7): JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 195 Alternating Schwarz: This method thus assumes that the problem is solved on each sub-domain successively. With the advent of computers with distributed architecture, the parallel Schwarz method (Equation 8) appeared. In this case, at each iteration, the problems per subdomain are solved in parallel, in order to take advantage of the potential of such computational means. Schwarz parallel In the case where we consider a division of the structure into two sub- domains, these two methods are equivalent since the sequences of solutions thus formed are the same in Ω1 or Ω2. We also speak of a multiplicative Schwarz method for the alternate method and an additive Schwarz method for the parallel method. The rewriting of these methods in the form of a fixed point reveals a sum of operators defined on Ω1 and Ω2 respectively for the parallel method and a product for the alternate method. Finally, we will mention here the "finite element patch" type methods [17- 18], generally used for solving multiscale problems. The idea is to consider the solution U as the sum of a global contribution UG and a local contribution UL provided by the patch so that U = UG + UL. As with the Schwarz methods, different strategies are then possible for the imposition of boundary conditions at each iteration. The main interest of such a method is to be able to calculate a local correction in a flexible way, i.e., to make the definition of the local patch independent of the characteristics of the global model. Contrary to direct methods, iterative methods require only weak interoperability between the different models (and thus solvers) used. This makes them very good candidates for the development of non- intrusive methods. 2.2.3 Connections between different models Whatever the coupling method used, it is necessary to establish a connection between the different models involved. Two types of connection can be distinguished: Interface connections, which aim to establish relationships between the quantities of interest at the interface. In this category, we will find weak connection methods such as the "Mortar" method [10-12]. The principle is to impose the equality of quantities (displacement, effort) in the weak sense, i.e., through a duality relation. Such a connection can be used within strong coupling methods to lead to hybrid finite element methods. Another example is the Nitsch method [13-15] which allows a Dirichlet edge condition to be imposed in a weak form, by adding an extra work term in the vibrational formulation of the problem. Volume connections through an energy average. In this case, the connection area between the different models is no longer an interface but a covering area: this is the Harlequin method [16]. While a volume connection may offer more flexibility in defining the connection area, it is also much more intrusive than an interface connection method. This is because a volume energy connection involves the simultaneous consideration of the quantities defined in each of the models. An interface connection, on the other hand, only requires independent unidirectional exchanges between the different models. Despite all these efforts, it sometimes appears that multi-scale/multi-model calculation methods cannot be used in practice. Indeed, in an industrial context, most numerical simulations are performed using commercial software that has been developed and certified for predefined applications. Thus, it is not always convenient to use such software to implement a given multi-scale/multi- model calculation. Moreover, the recent development of supercomputers now makes it possible to carry out very ambitious simulations. Current practices therefore increasingly consist in coupling different models rather than integrating all the required specificities into a single model. To overcome these difficulties, many authors have developed methods. This is what will be presented in the rest of this article. 3. Stationary coupling methods Suppose we have two codes 𝑓1 and 𝑓2 working on the variables 𝑥2 € R d2 and 𝑥1 R d1 . 𝑓1 gives us x1 knowing x2 , which we will note 𝑥1 = 𝑓1 (𝑥2) and vice versa, 𝑓2 gives us 𝑥2 as a function of 𝑥1 , i.e. 𝑥2 = 𝑓2 (𝑥1) . Carrying out the coupling, therefore, consists of looking for the global solution which verifies equation 9. Often the coupling will be formulated with a dependence of the functions 𝑓𝑖 in 𝑥𝑖 : There may be two reasons for this. This may be a choice not to try to converge each of the codes precisely until the common solution is reached. This is close to an intricate solution method. It also happens that it makes no sense to converge a code on xi without simultaneously updatingxj. By default, we will assume that the coupling is of the form (10). We will also sometimes rewrite the system in the following form: It is possible to a couple more than two equations with the techniques we present here. As the generalization is simple, to simplify the writings, we will limit ourselves to couplings between two disciplines. 3.1 Gauss-seidel method a. Presentation of the method The Gauss-Seidel method (often referred to as the Picard method or 196 JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 iterations in the literature), which implements a flexible coupling, is very common and easy to implement. It is applied to AX=B (linear) systems with the matrix A having to be diagonal dominant. In this method, it is simply a matter of iterating on the resolution of the codes until convergence. And one does not stop the iterations from precision. Nevertheless, given that the main concern is to see if the sequence of iterations converges and if so to what order. We will start by presenting the sequence of iterations, then we will find an equivalence with the fixed point equation and then we will check the convergence. Let us note with an exponent n the iteration index. From initial data (x1 (0) , x2 (0) ), we iterate by posing: { x1 (n+1) = f1(x1 (n) , x2 (n) ) x2 (n+1) = f2(x1 (n) , x2 (n) ). (12) ε1 (n+1) = ‖x1 (n+1) − x1 (n) ‖ et ε2 (n+2) = ‖x2 (n+1) − x2 (n) ‖ We generally stop when the increments become small enough. We use this criterion for want of a better one, and we see that it does not give any information on the distance to the solution. It is easy to see that the method if it converges, can only converge to a solution to the problem. However, convergence is not systematic. b. Equivalence with the fixed point The so-called "fixed point" algorithm is used to solve an equation of the type x = f(x). It consists in iterating as follows: x(n+1) = f(x(n)). (13) The Gauss-Seidel method can be seen as a fixed point by posing: x(n) = ( x1 (n) x2 (n) ) et f(x(n)) = ( f1(x1, x2) f2(x1, x2) ). (14) If one of our functions, say f1does not depend on the variable it is used to calculate, i.e. x1 (as is often the case), we notice that each iteration is equivalent to x2 (n+1) = f2( f1(x2 (n) ) , x2 (n) ) . This is another way of getting back to the fixed point, which has the advantage of reducing the dimension of the space in which we work. The following discussion of the Gauss- Seidel method uses this equivalence. We will therefore be interested in the properties of the fixed point. c. General condition for convergence Convergence condition: Let Ebe a complete metric space and f an application from E into E. If f is k contract i (i. e. ∀(a, b) ∈ E2 , d(f(a), f(b)) kd(a, b) with k < 1, d is the distance between two points of E.), then there exists a unique fixed point xsol of f and any sequence satisfying x(n+1) = f(x(n)) converges to xsol. d. Practical conditions for convergence The last remark leads us to look at the neighborhood of the solution to determine whether the algorithm can converge. We use the Taylor- Lagrange inequality for a function f in Rm in Rm of sufficient regularity. We note Ja the Jacobian of at a. If there is a real M such that ‖Ja+h‖ ≤ M ∀t ∈ ]0,1[ then: ‖f(a + h) − f(a)‖ ≤ M ‖h‖ . (15) We note xsol the fixed-point solution (we assume the existence and uniqueness of the solution here). It is assumed that there exists a real M such that ‖Jx‖ ≤ M ∀x ∈ B (xsol , r). B (xsol , r) is the ball of center xsol and radius r, any real. Then we have more generally: ∀a, b ∈ B(xsol, r ), ‖f(b) − f(a)‖ ≤ M ‖b − a‖ (16) Thus, if M < 1, this shows that f is contracting on B(xsol, r )and that the algorithm converges if x(0) ∈ B(xsol, r ) from the above. It is therefore sufficient that ‖Jx‖ ≤ M < 1 around the solution. If we use the norm L 2we can write the following sufficient condition (ρ(A) denotes the spectral radius of the matrix A, i.e. ρ(A) = maxi|λi| where the λi are the eigenvalues, possibly complex, of A). ‖Jx‖2 = √(ρ(Jx t Jx) ≤ M < 1 (17) This equation is valid in the vicinity of the solution, and the initialization is chosen in this vicinity. As before, the norm of the Jacobian allows us to reduce the speed of convergence of the algorithm. We can look for a less demanding condition by using its properties: ∀A, for any matrix norm we have: ρ(A) ≤ ‖A‖ . ∀A, ∀ε there exists a matrix norm such that ρ(A) ≤ ‖A‖ ≤ ρ(A) + ε. Thus if we have, in the solution xsol, ρ(Jxsol ) < 1, for a certain norm we have ‖Jxsol ‖ ≤ M < 1. By continuity of the norm, we then have with the same norm ‖Jx‖ ≤ M′ < 1 on a neighborhood of xsol . We thus have a sufficient condition of local convergence (thanks to the equivalence of norms in finite dimension. It suffices to show the convergence of the sequence for any norm), less demanding than the previous one since the first proposition of the lemma tells us that ρ(A) ≤ √ρ(At A) ∀A. We will try to show that ρ(Jxsol ) < 1 is in fact a necessary condition for "good We will try to show that ρ( ) < 1 is in fact a necessary condition of "good behavior of the algorithm" in a sense that we will specify. We will then say, by abuse of language, that it is a necessary condition of we will then say, by abuse of language, that it is a necessary condition of convergence. e. Variant with relaxation A variant of the fixed-point method is to introduce a relaxation α > 0. This is called over-relaxation when α > 1. The iterations are then written. 𝑥(𝑛+1) = (1 - α) 𝑥 (𝑛) + αf (𝑥 (𝑛)) (18) Let us see what this gives us about the criterion and the speed of convergence. We see that by posing g = (1 - α) Id + αf with Id identity matrix, we return to the classical fixed point with 𝑥(𝑛+1) = g(𝑥(𝑛)). The previous results therefore apply to g. Let us first consider the one- dimensional case. JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 197 One-dimensionnel case : We have 𝑔′ = (1-α)+α𝑓 ′. The fixed point converges at its maximum speed when the derivative of the function is zero in the solution. This leads us to an optimal relaxation: α = 1 1− 𝑓′(𝑥𝑠𝑜𝑙) . This value will generally be unknown since it is based on a property of (f) in the solution 𝑥𝑠𝑜𝑙 which is unknown a priori. We will therefore look more generally at the behaviour of the algorithm as a function of α. The convergence criterion |𝑔′(𝑥𝑠𝑜𝑙 ) |< 1 gives us 𝑓 ′(𝑥𝑠𝑜𝑙 ) ∈ ]1 - 2 𝛼 ; 1[. The speed of convergence depends on the maximum value of |𝑔′ | around the solution x. Compare 𝑔′ to 𝑓 ′. 𝑔′ - 𝑓 ′ = (1 - α)(1 - 𝑓 ′). Since 𝑓 ′ < 1, we have 𝑔′ > 𝑓 ′ α < 1. Now compare 𝑔′ à -𝑓 ′. 𝑔′+ 𝑓 ′ = (1 - α) + (α + 1) 𝑓 ′. Thus, 𝑔′(x) > -𝑓 ′(x) , α < 1+ 𝑓′(x) 1− 𝑓′(x) . Thus, we have: • If 𝑓 ′(𝑥𝑠𝑜𝑙 ) < -1, i.e. if the algorithm, without relaxation, does not converge with an oscillating behavior, a relaxation α ∈ ]0, 2 1− 𝑓′(𝑥𝑠𝑜𝑙) [ allows it to converge. The optimum is α = 1 1− 𝑓′(𝑥𝑠𝑜𝑙) . • If 𝑓 ′(𝑥𝑠𝑜𝑙 ) ∈ ] - 1, 0[, one can improve the speed of convergence with α ∈ ] 1 + 𝑓′(𝑥𝑠𝑜𝑙) 1− 𝑓′(𝑥𝑠𝑜𝑙) ,1[. The optimum is at α = 1 1− 𝑓′(𝑥𝑠𝑜𝑙) . • If 𝑓 ′(𝑥𝑠𝑜𝑙 ) ∈ 0, 1[, one can improve the speed of convergence with α ∈ ]1, 1 + 𝑓′(𝑥𝑠𝑜𝑙) 1− 𝑓′(𝑥𝑠𝑜𝑙) [. This is called over-relaxation. The optimum is at α = 1 1− 𝑓′(𝑥𝑠𝑜𝑙) . • If 𝑓 ′(𝑥𝑠𝑜𝑙 ) > 1, i.e. if the algorithm, without relaxation, has a non- convergent behavior without oscillation, it will keep this behavior ∀α > 0. 𝑓 ′(𝑥𝑠𝑜𝑙 ) Since the value of α is not known a priori, the difficulty is to determine α. More than one-dimensional case: To improve the convergence of g we will try to minimize the spectral radius of its Jacobian at x, the solution to the problem. We have the following general property, which is easy to check, where σ(A) denotes the spectrum of A, i.e. the set of its eigenvalues: ∀A, ∀P polynomial P (σ(A)) = σ(P (A)). Thus, the previous analysis, done in one dimension, applies to each of the eigenvalues of the Jacobian of g. Let us denote 𝜆𝑖 its eigenvalues, and 𝜆𝑖 (α) those obtained with relaxation. We have: 𝜆𝑖 (α) = 1 + α(𝜆𝑖 - 1). We seek to minimize 𝑚𝑎𝑥𝑖 |𝜆𝑖 (𝛼)|. As before, 𝜆𝑖 > 1 ⇔ 𝜆𝑖 (α) > 1 ∀α > 0. Therefore, relaxation does not allow the algorithm to converge if ∃𝜆𝑖 > 1, and we, therefore, consider that 𝜆𝑖 < 1 ∀i. Thus, we have 𝜆 𝑖 ′ (α) = (𝜆𝑖 - 1) < 0. All functions 𝜆𝑖 (α) have the same monotonicity and are all 1 in 0. Moreover, it is easy to verify that if 𝜆𝑖 ≠ 𝜆𝑗 then 𝜆𝑖 (α)≠ 𝜆𝑗 (α) ∀α ≠ 0. Let 𝜆𝐼 be the eigenvalue which realizes 𝑚𝑎𝑥𝑖 | 𝜆𝑖 |. It is easy to see that as long as 𝜆𝐼 (α) realizes 𝑚𝑎𝑥𝑖 |𝜆𝑖 (𝛼)|we can still decrease the spectral radius. And since, |𝜆𝑖 (𝛼)| = |𝜆𝑗 (𝛼)| ⟹ 𝜆𝑖 (α) = - 𝜆𝑗 (α), we see that when ∃ i ≠ 𝐼 𝑡𝑒𝑙 𝑞𝑢𝑒 𝜆𝑖 ≠ 𝜆𝐼 𝑒𝑡 | 𝜆𝐼 (𝛼)| = |𝜆𝑖 (𝛼)| . Since the monotonies are identical, we can no longer minimize both |𝜆𝐼 (𝛼)| and |𝜆𝑖 (𝛼)| . The optimum is thus reached. We notice that, since we have 𝜆𝑖 (α) < 1 ∀α > 0, the convergence of the algorithm is guaranteed at the optimum. 𝜆𝑖 (α)= -𝜆𝑗 (α) gives us: α = 1 1− 𝜆𝑖 +𝜆𝑗 2 . Moreover, since the 𝜆𝑖 (α) cannot intersect, and 𝜆𝐼 corresponds to either the minimum or the maximum of the eigenvalues, the first 𝜆𝑖 ≠ 𝜆𝐼 which realizes | 𝜆𝐼 (𝛼)| = |𝜆𝑖 (𝛼)|is either the maximum (if 𝜆𝐼 is the minimum) or the minimum of the eigenvalues. The optimal relaxation, which we have seen ensures convergence, is therefore: α = 1 1− 𝑚𝑎𝑥𝑖𝜆𝑖 + 𝑚𝑖𝑛𝑖𝜆𝑖 2 . As in the one-dimensional case, the optimal relaxation we have identified requires knowledge of properties of (f) in the solution of the problem. It will therefore generally be necessary to use approximate methods. A very simple empirical procedure, which works whatever the dimension of the problem, often gives good results. It consists in dividing the relaxation by two each time the convergence criterion increases (be careful with its definition). It is possible, but not essential, to recalculate the iterations where the convergence criterion has increased. f. Jacobi iterations The Jacobi method is similar to the Gauss-Seidel method. The difference is that the result of the previous calculation is not used to improve the result of the next calculation. This is equivalent to iterating as follows (n is the iteration number): { x1 (n+1) = f1(x1 (n) , x2 (n) ) x2 (n+1) = f2(x1 (n) , x2 (n) ). (19) We can see that it is also a fixed-point method. The convergence is slower than for Gauss-Seidel, but the algorithm is parallelizable. The Gauss-Seidel method is generally preferred to the Jacobi method. Thus we can say that the Gauss-Seidel method is easy to implement and often converges. In particular, the convergence takes place whatever the initialization on the coupling is sufficiently weak. However, the method cannot converge if the couplings involved are too strong. Adding a relaxation can be a solution in this case, but the difficulty is then to parameterize it. Moreover, the linear convergence is all the faster that the functions studied are weakly coupled. 3.2 Newton’s method a. Presentation of the method This is another common method, which has many variants. In its original formalism, it is an intricate (and therefore firm) coupling method. Let us put the system to be solved as follows (F1 and F2 can be taken from the formulation (20) for example): F (x) = ( F1(x1, x2) F2(x1, x2) ) = 0, x = (x1, x2) . (20) It is recalled that, in general, xi and Fi(xi, xj) are vectors of R di . Note that d= d1 + d2 the number of components of x and define the Jacobian Jx of the system at x as follows (the Fi are here and for the following the components of F, and thus 198 JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 Jx=( ∂F1 ∂x1 (x) ⋯ ∂F1 ∂xd (x) ⋮ ⋱ ⋮ ∂Fd ∂x1 (x) ⋯ ∂Fd ∂xd (x) ) (21) Like the Gauss-Seidel method, the Newton method is iterative. Let us note here again with an exponent n the iteration index. From an initial data x(0) we iterate as follows: We solve in δx the linear system: Jx(n) δx = −F(x (n). And then we pose x(n+1) = x(n) + δx . We generally stop when ‖δx ‖ or ‖F(x(n))‖ becomes small enough. As before, we see that the method can only converge to a solution of the problem. b. Study of the method The method is based on a Taylor expansion to order 1 around the solution xsol : F(xsol) = 0 ≈ F(x (n)) + Jx(n) (xsol − x (n)). (22) Assuming F is sufficiently regular, let us use the Taylor-Lagrange inequality to quantify the error. As before, let us note D2Fx the linear application of Rd in L(Rd), such that (D 2Fxab)k =∑ ∂2Fk(x) ∂xi ∂xj aibji,j . If there is a real M such that: ‖D2Fx(n)+t(xsol−x(n))‖ ≤ M ∀t ∈ ]0; 1[ then: ‖−𝐹(𝑥(𝑛)) + 𝐽𝑥(𝑛) (𝑥𝑠𝑜𝑙 − 𝑥 (𝑛))‖ ≤ M‖xsol−x (n)‖ 2 2 (23) Since F(𝑥𝑠𝑜𝑙) = 0. We introduce the vector 𝜂𝑛of norm less than 1, such that: 𝐹(𝑥(𝑛)) = 𝐽𝑥(𝑛) (𝑥𝑠𝑜𝑙 − 𝑥 (𝑛)) + 𝜂𝑛 M‖xsol−x (n)‖ 2 2 . (24) We can therefore write the norm of the error as follows, assuming the Jacobian is invertible: ‖ℰ (𝑛+1)‖ = ‖𝑥 (𝑛+1) − 𝑥𝑠𝑜𝑙 ‖ = ‖𝑥 (𝑛) − 𝐽 𝑥(𝑛) −1 𝐹(𝑥(𝑛)) − 𝑥𝑠𝑜𝑙 ‖ (25) ‖ℰ (𝑛+1)‖ = ‖𝑥(𝑛) − 𝐽 𝑥(𝑛) −1 ( 𝐽𝑥(𝑛) (𝑥𝑠𝑜𝑙 − 𝑥 (𝑛)) + 𝜂𝑛 M‖xsol−x (n)‖ 2 2 ) − 𝑥𝑠𝑜𝑙 ‖ (26) ‖ℰ (𝑛+1)‖ = ‖𝐽 𝑥(𝑛) −1 𝜂𝑛‖ M‖ℰ (𝑛)‖ 2 2 ≤ ‖𝐽 𝑥(𝑛) −1 ‖ M‖ℰ (𝑛)‖ 2 2 . (27) Let K = max𝑚𝑎𝑥𝑥∈𝐵(𝑥𝑠𝑜𝑙,‖𝑥𝑠𝑜𝑙−𝑥(0)‖( 1 2 ‖𝐽𝑥 −1‖‖D2Fx‖). 𝐵(𝑥𝑠𝑜𝑙 , ‖𝑥𝑠𝑜𝑙 − 𝑥 (0)‖ is the ball of center 𝑥𝑠𝑜𝑙 and radius ‖𝑥𝑠𝑜𝑙 − 𝑥 (0)‖ convergence occurs if the following criterion is met: 𝐾‖𝑥(0) − 𝑥𝑠𝑜𝑙 ‖< 1. (28) Furthermore, it can be shown by recurrence that we have: ‖ℰ (𝑛)‖ ≤ (𝐾‖𝑥(0)−𝑥𝑠𝑜𝑙‖) 2𝑛 𝐾 . (29) This shows us that we have a quadratic convergence, thus faster than the Gauss-Seidel one, if the initialization is good enough. Newton’s algorithm can only be defined correctly if we can define the Jacobian of F at all points. However, it happens that the phenomena we are studying introduce breaks in the slope. We are thinking of phase changes. In this case, the previous equations assure us that the convergence of the algorithm always takes place, for a sufficiently good initialization, whether we choose the partial derivatives on the right or on the left at the breaks in the slope. c. Variant with relaxation In the same way that we introduced a relaxation for Gauss-Seidel, it is possible to define one for Newton, which we will note α ∈]0; 1] and which occurs as follows (δx is always calculated in the same way): x(n+1)=x(n)+αδx . (30) It is easy to see that when α ≠ 1 one loses the quadratic convergence speed of Newton's algorithm. Relaxation cannot therefore be used to accelerate convergence, but it can allow convergence in certain situations where the initial algorithm would not have done so. The Jacobian is then mainly used to determine the direction in which the x(n) and α determines the size of the "step" taken. In principle, we can see that we are approaching a gradient descent type of algorithm, used. Let us insist on the fact that α = 1 will always be optimal when approaching the solution for the Newton algorithm, which was not the case for Gauss-Seidel. When using a relaxation with Newton's algorithm, one must therefore have the possibility, at each iteration, to increase α. As with Gauss-Seidel, there is a simple empirical procedure for using relaxation with the Newton method. This involves adding an internal iteration level. Each Newton iteration starts with α = 1. This iteration is then recalculated by dividing the relaxation by two each time, as long as the resulting residual is not smaller than that of the previous iteration. d. Newton's algorithm does not correct It is simply a Newtonian algorithm in which it is assumed that the solution of the system Jx(n) δx = −F(x (n)) is not exact. Starting with an initial data x(0) we iterate as follows: We determine δx which verifies: ‖Jx(n) δx + F(x(n))‖ ≤ 𝜂(n)‖F(x(n))‖. And we pose x(n+1) = x(n)+ α(n)δx. 𝜂 (n) ∈ [0; 1] fixes the required precision on the resolution of the linear system and α(n)> 0 corresponds to a relaxation as already introduced. As the resolution of the linear system is assumed to be inaccurate, the algorithm does not benefit from the quadratic convergence speed of the Newton and the introduction of the relaxation is less troublesome. 3.3 JFNK , Jacobian-free Newton-Krylov, method a. Presentation of the method In particular, it is the basis of the MOOSE coupling platform [19]. The JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 199 JFNK is an implementation of Newton's method based on a Krylov method for solving linear systems, which is applied at each iteration of Newton to the system 𝐽𝑥(𝑛) δx = -F (x (n) ). Krylov's methods have the advantage of involving the system matrix only through matrix-vector products, which we are able to approximate by simple code calls. The JFNK method thus allows to define a sequential (codes are kept) but firm coupling technique. b. Principle of Krylov methods for solving Ax = b Let mbe the dimension of the vector b. Krylov's methods consist in finding an approximation of the solution on a vector space Kn of reduced dimension n ≤ m, which is constructed as follows: Kn = Vect(r0; Ar0,..., A n-1 r0), where r0 = b - Ax0, with x0 an initial estimate of the solution of the linear system. In our case x0 is usually taken to be zero and we simply have r0 = b. To simplify the writings, we will work in the following directly with the variable b. Note that (b, Ab, ..., An-1 b) does not necessarily form a family of free vectors, but that if this is not the case, the solution x can be generated by the family (b, Ab, ..., An-2 b). c. Use and interest of the JFNK method in the context of coupling In a problem involving a Jacobian, one can use the following Taylor expansion to approximate the matrix-vector product: 𝐽𝑢 𝜈 = 𝐹(𝑢+𝜀𝜈)− 𝐹(𝑢) 2𝜀 + 𝑂(𝜀). (31) An evaluation of F corresponds to a call to the codes we couple. The method can therefore be used without the need to compute or store the Jacobian of our coupled problem but requires frequent calls to the codes. We thus need one call of the codes per GMRES iteration, in addition to the necessary call per Newton iteration. We see that our approximation of the product 𝐽𝑢 𝜈 product introduces an error in O (𝜀). This error can be reduced by using a centered difference, but then we need two calls of the codes per iteration of GMRES. 𝐽𝑢 𝜈 = 𝐹(𝑢+𝜀𝜈)− 𝐹(𝑢−𝜀𝜈) 2𝜀 + 𝑂(𝜀 2). (32) Let us insist on the fact that the JFNK is only a way to implement Newton's algorithm by circumventing the difficulties linked to the construction of the Jacobian. Newton's algorithm by circumventing the difficulties associated with the construction of the Jacobian. The discussion about Newton's method applies here. The JFNK introduces new parameters which are the 𝜀 used to calculate the product 𝐽𝜈 and the convergence criterion of the GMRES method. The former must be chosen small enough to make the Taylor expansion as fair as possible, but large enough to avoid numerical noise. The GMRES convergence criterion is a compromise between accuracy and calculation time. 4. Transient coupling methods 4.1. Method directly derived from stationary techniques In the previous methods, the time variable does not intervene. We can then ask ourselves how to achieve a coupling in kinetics. To explain it, let us start by giving a general form of the coupling, as we did previously: { ∂x1 ∂t = G1 (x1 , x2 ) ∂x2 ∂t = G2 (x1 , x2 ) (33) Except for very special systems, after choosing a numerical scheme, we obtain resolutions of the type: { 𝑥1(𝑡 + ∆𝑡) = 𝑥1(𝑡) + ℎ1(𝑥2(𝑡 + ∆𝑡), ∆𝑡) 𝑥2(𝑡 + ∆𝑡) = 𝑥2(𝑡) + ℎ2(𝑥1(𝑡 + ∆𝑡), ∆𝑡) (34) Where the hi are continuous functions such that lim ∆𝑡→0 ℎ𝑖 (𝑥, ∆𝑡) = 0 ∀𝑥. By posing: { 𝐹1 (𝑢1, 𝑢2) = 𝑢1 − 𝑥1(𝑡) − ℎ1(𝑢2, ∆𝑡) 𝐹2 (𝑢1, 𝑢2) = 𝑢2 − 𝑥2(𝑡) − ℎ2(𝑢1, ∆𝑡) (35) We return to the form already studied: { 𝐹1 (𝑢1, 𝑢2) = 0 𝐹2 (𝑢1, 𝑢2) = 0 The methods already seen are therefore directly applicable, but we have an additional parameter that allows us to influence the strength of the coupling thus defined at each timestep: the time step itself. In particular, when using a Gauss-Seidel method, the convergence criterion (∃𝑀 < 1, ‖𝐽‖ ≤ 𝑀) is automatically verified for ∆t sufficiently small. Unfortunately, the time step necessary for Gauss-Seidel is sometimes too small to be exploitable (calculation of the transient too long), and we then turn to the other methods already presented. Let us note that there is a purely technical difficulty, but sometimes difficult to circumvent, for the implementation of evolving coupled schemes: the methods presented are iterative, it is thus necessary to be able to recalculate the same time step several times, without modifying the variables at the beginning of the time step (the 𝑥𝑖 (𝑡)in our notation). However, the codes sometimes forget 𝑥𝑖 (𝑡) during the calculation of 𝑥𝑖 (𝑡 + ∆𝑡). 4.2. Explicit coupling To save calculation time when the couplings are weak, but also because of the technical difficulty we have just mentioned, it is often desired to calculate a transient coupled without iteration at each time step. To do this, we will look at how an error evolves during the calculation of a transient, whether or not the phenomenon being modelled arises from a coupling. Let us write a general form for an evolution problem: 𝜕𝑥 𝜕𝑡 = 𝐹(𝑥). It is assumed that a classical explicit scheme is used, which gives us: x (t + ∆t) = x(t) + ∆t F (x(t)). (36) Let us note 𝑥𝑒 the exact solution of the problem and suppose that at the 200 JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 beginning of the time step time step x(t) = 𝑥𝑒 (t) + 𝜀(t). We will also note 𝜀′ the error introduced by the scheme itself (which is of order 1 in time) on the time step we are calculating: 𝑥𝑒 (t) + ∆t F (𝑥𝑒 (t)) = 𝑥𝑒 (t + ∆t) + 𝜀′ (t + ∆t). (37) With this we have : x(t + ∆t) = 𝑥𝑒 (t) + 𝜀(t) + ∆t F (𝑥𝑒 (t) + 𝜀(t)) = 𝑥𝑒 (t) + 𝜀(t) + ∆t (F(𝑥𝑒 ( t)) + F(𝑥𝑒 ( t) + 𝜀 (t)) - F(𝑥𝑒 ( t))) = 𝑥𝑒 (t + ∆t) + 𝜀'(t + ∆t) + 𝜀 (t) + ∆t (F (𝑥𝑒 ( t) + 𝜀 (t)) - F (𝑥𝑒 ( t))) Eq (38) Thus we can write 𝜀(t + ∆t) = f(𝜀 (t)) with : f(u) = u + ∆t (F (𝑥𝑒 ( t) + u) - F (𝑥𝑒 ( t)) )+ 𝜀'(t + ∆t) (39) As we move forward in time, the error committed evolves almost as if we were doing fixed point iterations with the function f thus defined. However, the equivalence is not perfect since 𝑥𝑒 and 𝜀′ evolve with time. To understand what is going on, let us suppose that these quantities are constant (which amounts to saying that we verify the stability of a stationary state), and that we are in one dimension. In this case, convergence occurs if |𝑓′(𝑢)| < 1, 𝑠𝑜𝑖𝑡 𝑠𝑖 𝐹′(𝑢) 𝜖 ] −2 ∆𝑡 , 0[ . The upper bound, F'(u) < 0, means that the system under study must itself be (statically) stable. As for the lower bound, we notice that by decreasing the time step it can be arbitrarily small. Thus, with a sufficiently small-time step, the error committed by an explicit coupling scheme can remain bounded. Nevertheless, the error committed is in O(∆t). The explicit scheme is therefore accurate to first order in time. The formalism used is valid when we evolve a coupled problem in time, and we do not converge the coupling at each time step. Thus, the way in which the coupling is performed can reduce the accuracy of the solution to first order in time, even though the solvers themselves have better accuracy. This point is particularly important. By saving the work of thinking about the coupling scheme, one risks losing the benefit of all the work invested in the codes themselves. 5. Synthesis The Gauss-Seidel method is very easy to implement and is the first method to be proposed in the genesis of coupling. This does not mean that the Newton algorithm is superior in all respects to the Gauss-Seidel algorithm. The main advantage of Newton's method is its quadratic rather than linear convergence speed. It is, therefore, to be preferred if a good coupling is desired. If one is simply interested in the solution of the coupled problem, the Gauss-Seidel algorithm is generally sufficient. The use of a relaxation, which can be determined empirically, is sufficient to make it very robust. Only those cases, which are actually quite rare, where the algorithm diverges without oscillating cannot be stabilized with a relaxation. Newton's algorithm, on the other hand, can be used theoretically on any case, but on condition that the initialization is good enough. Even if the use of a relaxation is likely to greatly improve its behavior, the Newton algorithm is sometimes less robust than the Gauss-Seidel algorithm, because it is more sensitive to noise. The comparative table summarizes our analysis by presenting the advantages and disadvantages of each method presented above. It shows that Newton's method is the most appropriate for the coupling of disciplines in a mechatronic system. Table1 - Table of methods with their advantages and disadvantages Methods Advantages Disadvantages Mechatronic systems Direct methods Suitable for medium and large sized systems. For large systems, it requires a lot of memory size. Not adequate because there is no stability. Gauss Seidel methods Suitable for linear systems and linear convergence. Choice of the stopping condition. Not adequate because mechatronic systems are linear. Newton's methods Appropriate for nonlinear systems and quadratically convergent. Choice of initial data. Adequate. JFNK method Suitable for nonlinear systems, More memorable, includes more parameters which makes calculations easier. Choice of initial data. Adequate but with too many parameters. Explicit pairing method Time saving. Reduced accuracy due to partitioned convergence. Not Adequate because poor precision emerges. 6. Comparative table At the end of this coupling analysis, we have noticed that all systems, whatever the size, can be coupled by coupling methods. Depending on the problem to be solved we have to choose the adequate method. Direct methods if the system is of average size but also of large size while envisaging much memory size. Hence the interest in using iterative methods, especially Gauss Seidel and Jacobi. Nevertheless, they are used for linear systems. However, mechatronic systems are non-linear. This justifies the choice of Newton's method. Therefore, for non-linear systems, we have to use the Newtonian method while taking care to choose the initial data carefully. 7. Conclusion At the end of this coupling analysis, we have noticed that all systems, whatever the size, can be coupled by coupling methods. Depending on the problem to be solved we have to choose the adequate method. That is to say direct methods if the system is of average size but also of large size while envisaging much memory size. Hence the interest in using iterative methods, especially Gauss Seidel and Jacobi. Nevertheless, they are used for linear systems. However, mechatronic systems are non-linear. This justifies the choice of Newton's method. Therefore, for non-linear systems, we have to use the Newtonian method while taking care to choose the initial data carefully. Acknowledgements JEAN SAMON ET AL. /AL-QADISIYAH JOURNAL FOR ENGINEERING SCIENCES 15 (2022) 192–201 201 As corresponding author of this paper, I am the scientific coordinator of the Laboratory of Mechanics, Materials, and Photonics. It is in this laboratory that this research was conducted. I would like to thank all my colleagues for their frank collaborations. References [1] David et al 2013 Special Issue on Multiphysics simulations: Challenges and opportunities, Volume 27 Number 1 February 2013 The International Journal of High Performance Computing Applications 27(1) 4–83. The. DOI: 10.1177/1094342012468181 hpc.sagepub.com [2] Miladi C., R. Plateauxb , J.-Y. Choleyb , C. Karraa , A. Riviereb M. Haddar : New topological approach for the modelling of mecatronic systems: application for piezoelectric structures, European Journal of Computational Mechanics, 2013 Vol. 22, Nos. 2–4, 209–227, http://dx.doi.org/10.1080/17797179.2013.820896 [3] K. C. Park and Carlos A. Felippa : A Variational Principle for the Formulation of Partitioned Structural Systems, Department of Aerospace Engineering Sciences and Center for Aerospace Structures University of Colorado, Campus Box 429 Boulder, CO 80309, USA, 1999. [4] Cresta P: Décomposition de domaine et stratégies de relocalisation non-linéaire pour la simulation de grandes structures raidies avec flambage local, thèse 2008 https://tel.archives-ouvertes.fr/tel-00363656/document [5] Lionel G. : Approche globale / locale non-intrusive : application aux structures avec plasticité localisée. Mécanique [physics.med-ph]. École normale supérieure de Cachan - ENS Cachan, 2009. Français. https://tel.archives- ouvertes.fr/tel-00449687/document [6] E. Wyarta, M. Duflota, D. Coulonb, P. Martinya, T. Pardoenc , J.-F. Remacled, F. Lania : Substructuring FE–XFE approaches applied to three-dimensional crack propagation, Journal of Computational and Applied Mathematics. doi:10.1016/j.cam.2006.03.066 [7] E.Wyart D.Coulon T.Pardoen F.Remacle F.Lani : Application of the substructured finite element/extended finite element method (S-FE/XFE) to the analysis of cracks in aircraft thin walled structures, Engineering Fracture Mechanics, Volume 76, Issue 1, January 2009, Pages 44-58 https://doi.org/10.1016/j.engfracmech.2008.04.025 [8] Mehmet A. Akgün, John H. Garcelon, Raphael T. Haftka : Fast exact linear and non-linear structural reanalysis and the Sherman–Morrison–Woodbury formulas, 26 January 2001 https://doi.org/10.1002/nme.87 [9] Mickaël D. Passieux J. Michel S. Stéphane G : Non-intrusive Coupling: Recent Advances and Scalable Nonlinear Domain Decomposition , Archives of Computational Methods in Engineering. 2014. DOI 10.1007/s11831-014-9132- x [10] Bernardi C. Nicolas F. Robert G. Owens : An error indicator for mortar element solutions to the Stokes problem, IMA Journal of Numerical Analysis ( Volume: 21, Issue: 4, Oct. 2001), Page(s): 857 – 886. DOI: 10.1093/imanum/21.4.857 [11] Belgacem F : The Mortar finite element method with Lagrange multipliers Laboratoire de Mathématiques pour l'Industrie et la Physique, Unité Mixte de Recherche CNRS–UPS–INSAT–UT1 (UMR 5640), Université Paul Sabatier, 118 route de Narbonne, F-31062 Toulouse Cedex 04, France , FR, volume 84, pages173–197 (1999) [12] Bernardi C, Yvon M. Francesca R. : Basics and some applications of the mortar element method, Volume28, Issue2, 08 May 2014, Pages 97-123. https://doi.org/10.1002/gamm.201490020 [13] Anita H. Peter H. : An unfitted finite element method, based on Nitsche’s method, for elliptic interface problems, Computer Methods in Applied Mechanics and Engineering Volume 191, Issues 47–48, 22 November 2002, Pages 5537-5552,https://doi.org/10.1016/S0045-7825(02)00524-8 [14] Becker R. Peter H. Rolf S. : A finite element method for domaindecomposition with non-matching grids, ESAIM: Mathematical Modelling and Numerical Analysis , Volume 37 , Issue 2 , March 2003 , pp. 209 – 225, DOI: https://doi.org/10.1051/m2an:2003023 [15] Nguyen V. Pierre K. Marco B. Stéphane P. Elvio B. : Nitsche’s method for two and three dimensional NURBS patch coupling, Computation Mechanics, Volume 53, pages1163–1182 (2014) https://link.springer.com/article/10.1007/s00466-013-0955-3 [16] Paul T. Hachmi B. Nadia E. Tinsley O. Serge P : On the application of the Arlequin method to the coupling of particle and continuum models, Computational Mechanics volume 42, pages511–530 (2008) https://link.springer.com/article/10.1007/s00466-008-0291-1 [17] Picasso M. Jacques R. Vittoria R. : Multiscale algorithm with patches of finite elements, Research Article, Communications in Numerical Methods in Engineering, 31 May 2007 https://doi.org/10.1002/cnm.1019 [18] Alexei L. and Pironneau O : Numerical zoom for advection diffusion problems with localized multiscales, Mathematics, Applied, October 2010 https://doi.org/10.1002/num.20642 https://onlinelibrary.wiley.com/journal/10982426 [19] Derek G. Chris N. Glen H. Damien L : MOOSE: A parallel computational framework for coupled systems of nonlinear equations, Nuclear Engineering and Design, Volume 239, Issue 10, October 2009, Pages 1768-1778 https://doi.org/10.1016/j.nucengdes.2009.05.021 [20] Cyril P. : Couplages multi-physiques : évaluation des impacts méthodologiques lors de simulations de couplages neutronique/thermique/mécanique. https://www.theses.fr/2016SACLX007 [21] Shreyansh J, Victoire M. , H. N. S. Narayana, Simon B, Joseph D, Victor C, Tianchi C., L. Heuzé, Philippe M, René-Marc M, J. Kabla, Teck L & Benoit L. : The role of single-cell mechanical behaviour and polarity in driving collective cell migration. https://www.nature.com/articles/s41567-020-0875-z.pdf [22] Mao J : Implementation of the analytic nodal method in the ndf code and applications to candu reactors, Ecole polytechnique de Montréal, 2000 https://publications.polymtl.ca/8798/ [23] Thomas Heuzé : Couplages Multiphysiques, 2020 https://hal.archives- ouvertes.fr/hal-02448292 http://dx.doi.org/10.1080/17797179.2013.820896 https://tel.archives-ouvertes.fr/tel-00363656/document https://tel.archives-ouvertes.fr/tel-00449687/document https://tel.archives-ouvertes.fr/tel-00449687/document https://www.sciencedirect.com/science/article/abs/pii/S0013794408001203#! https://www.sciencedirect.com/science/article/abs/pii/S0013794408001203#! https://www.sciencedirect.com/science/article/abs/pii/S0013794408001203#! https://www.sciencedirect.com/science/article/abs/pii/S0013794408001203#! https://www.sciencedirect.com/science/article/abs/pii/S0013794408001203#! https://www.sciencedirect.com/journal/engineering-fracture-mechanics https://www.sciencedirect.com/journal/engineering-fracture-mechanics https://www.sciencedirect.com/journal/engineering-fracture-mechanics/vol/76/issue/1 https://doi.org/10.1016/j.engfracmech.2008.04.025 https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Akg%C3%BCn%2C+Mehmet+A https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Garcelon%2C+John+H https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Haftka%2C+Raphael+T https://doi.org/10.1002/nme.87 https://link.springer.com/journal/11831 https://link.springer.com/journal/11831 https://ieeexplore.ieee.org/xpl/RecentIssue.jsp?punumber=8016799 https://ieeexplore.ieee.org/xpl/tocresult.jsp?isnumber=8144758 https://doi.org/10.1093/imanum/21.4.857 https://link.springer.com/article/10.1007/s002110050468#auth-Faker_Ben-Belgacem https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Bernardi%2C+Christine https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Maday%2C+Yvon https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Rapetti%2C+Francesca https://onlinelibrary.wiley.com/toc/15222608/2005/28/2 https://doi.org/10.1002/gamm.201490020 https://www.sciencedirect.com/science/article/abs/pii/S0045782502005248#! https://www.sciencedirect.com/science/article/abs/pii/S0045782502005248#! https://www.sciencedirect.com/journal/computer-methods-in-applied-mechanics-and-engineering https://www.sciencedirect.com/journal/computer-methods-in-applied-mechanics-and-engineering https://www.sciencedirect.com/journal/computer-methods-in-applied-mechanics-and-engineering/vol/191/issue/47 https://doi.org/10.1016/S0045-7825(02)00524-8 https://www.cambridge.org/core/search?filters%5BauthorTerms%5D=Roland%20Becker&eventCode=SE-AU https://www.cambridge.org/core/search?filters%5BauthorTerms%5D=Peter%20Hansbo&eventCode=SE-AU https://www.cambridge.org/core/search?filters%5BauthorTerms%5D=Rolf%20Stenberg&eventCode=SE-AU https://www.cambridge.org/core/journals/esaim-mathematical-modelling-and-numerical-analysis https://www.cambridge.org/core/journals/esaim-mathematical-modelling-and-numerical-analysis https://www.cambridge.org/core/journals/esaim-mathematical-modelling-and-numerical-analysis/volume/4AD21D52450796C5948832D41421AB15 https://www.cambridge.org/core/journals/esaim-mathematical-modelling-and-numerical-analysis/issue/184274FA06AEAF55E5323CCB98B9969E https://doi.org/10.1051/m2an:2003023 https://link.springer.com/article/10.1007/s00466-013-0955-3#auth-Vinh_Phu-Nguyen https://link.springer.com/article/10.1007/s00466-013-0955-3#auth-Pierre-Kerfriden https://link.springer.com/article/10.1007/s00466-013-0955-3#auth-Marco-Brino https://link.springer.com/article/10.1007/s00466-013-0955-3#auth-St_phane_P__A_-Bordas https://link.springer.com/article/10.1007/s00466-013-0955-3#auth-Elvio-Bonisoli https://link.springer.com/article/10.1007/s00466-013-0955-3 https://link.springer.com/article/10.1007/s00466-008-0291-1#auth-Nadia-Elkhodja https://link.springer.com/article/10.1007/s00466-008-0291-1#auth-J__Tinsley-Oden https://link.springer.com/article/10.1007/s00466-008-0291-1#auth-Serge-Prudhomme https://link.springer.com/journal/466 https://link.springer.com/article/10.1007/s00466-008-0291-1 file:///C:/Users/Damasse/Desktop/Articles%20FOTSA/Article%20Ndomou/soumettre/Picasso file:///C:/Users/Damasse/Desktop/Articles%20FOTSA/Article%20Ndomou/soumettre/Jacques%20R https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Rezzonico%2C+Vittoria https://onlinelibrary.wiley.com/journal/10990887 https://onlinelibrary.wiley.com/journal/10990887 https://doi.org/10.1002/cnm.1019 https://onlinelibrary.wiley.com/action/doSearch?ContribAuthorRaw=Lozinski%2C+Alexei file:///C:/Users/Damasse/Desktop/Articles%20FOTSA/Article%20Ndomou/soumettre/%20Pironneau https://doi.org/10.1002/num.20642 https://onlinelibrary.wiley.com/journal/10982426 https://www.sciencedirect.com/science/article/abs/pii/S0029549309002635#! https://www.sciencedirect.com/science/article/abs/pii/S0029549309002635#! https://www.sciencedirect.com/science/article/abs/pii/S0029549309002635#! https://www.sciencedirect.com/science/article/abs/pii/S0029549309002635#! https://www.sciencedirect.com/journal/nuclear-engineering-and-design https://www.sciencedirect.com/journal/nuclear-engineering-and-design https://www.sciencedirect.com/journal/nuclear-engineering-and-design/vol/239/issue/10 https://doi.org/10.1016/j.nucengdes.2009.05.021 https://www.theses.fr/191823120 https://www.theses.fr/2016SACLX007 https://www.nature.com/articles/s41567-020-0875-z#auth-Shreyansh-Jain https://www.nature.com/articles/s41567-020-0875-z#auth-Victoire_M__L_-Cachoux https://www.nature.com/articles/s41567-020-0875-z#auth-Gautham_H__N__S_-Narayana https://www.nature.com/articles/s41567-020-0875-z#auth-Simon-Beco https://www.nature.com/articles/s41567-020-0875-z#auth-Joseph-D_Alessandro https://www.nature.com/articles/s41567-020-0875-z#auth-Victor-Cellerin https://www.nature.com/articles/s41567-020-0875-z#auth-Victor-Cellerin https://www.nature.com/articles/s41567-020-0875-z#auth-Tianchi-Chen https://www.nature.com/articles/s41567-020-0875-z#auth-M_lina_L_-Heuz_ https://www.nature.com/articles/s41567-020-0875-z#auth-Philippe-Marcq https://www.nature.com/articles/s41567-020-0875-z#auth-Ren__Marc-M_ge https://www.nature.com/articles/s41567-020-0875-z#auth-Alexandre_J_-Kabla https://www.nature.com/articles/s41567-020-0875-z#auth-Chwee_Teck-Lim https://www.nature.com/articles/s41567-020-0875-z#auth-Chwee_Teck-Lim https://www.nature.com/articles/s41567-020-0875-z#auth-Benoit-Ladoux https://www.nature.com/articles/s41567-020-0875-z.pdf https://publications.polymtl.ca/8798/