Allen Starke
email: allen1.starke@ufl.edu

Janise Mcnair
email: mcnair@ece.ufl.edu

Rodrigo Trevizan
email: rodtrevizan@ufl.edu

Arturo Bretas
email: arturo@ece.ufl.edu

Joshua Peeples
email: jpeeples@ufl.edu

Alina Zare
email: azare@ece.ufl.edu

Toward Resilient Smart Grid Communications Using Distributed SDN with ML-Based Anomaly Detection

Keywords: software defined networks, anomaly detection, machine learning, security, resilience

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

The next-generation power grid, named the Smart Grid, has drawn the attention of academia, industry and government agencies due to the great impact of such systems on the distribution of power within and between various regions. These next generation systems integrate control, communication and computation to achieve stability, efficiency and robustness of the physical processes. While a great amount of research has addressed these objectives, science and technology related to secure SG communications is still relatively immature. Additionally, many critical cyber-physical infrastructures are increasing dependency of control of physical processes on communication networks, thus becoming exposed to various cyber-threats or faulty operation. An example is the network of smart grid subsystems shown in Figure 1. The smart grid subsystems are local agents, composed of Distributed Energy Storage Systems (DESS), such as flywheels and grid-connected batteries, a Synchronous Generator, a Phasor Measurement Unit (PMU) and a Distributed State Estimator (DSE). The PMUs in the smart grid form a communication network, where smart grid data is regularly exchanged with other PMUs to coordinate and analyze energy performance measures. In the attack model considered, a PMU subsystem may be attacked in order to disable it or to inject faulty data; or a PMU element may fail, generating faulty data as a result of the failure. If a PMU subsystem is attacked or fails, the communication of invalid data may introduce substantial errors when exchanged with the other connected PMU subsystems, creating an avalanche effect. We propose to secure the PMU subsystem network by introducing a hybrid distributed and decentralized (HDD) software-defined network (SDN). The HDD-SDN architecture will leverage data and statistics from PMU communications to gain situational awareness and will provide machine intelligence to detect and protect against abnormal network behavior within the distributed PMU communication network. SDN is a networking paradigm in which the forwarding hardware is decoupled from the control decisions. The network intelligence is logically centralized in software-based controllers (the control plane), and the network devices become simple packet forwarding devices (the data plane) that can be programmed via an open interface [START_REF] Monsanto | Composing software defined networks[END_REF]. The controller is the only source responsible for determining routing paths, developing policies, partitioning the network, as well as other network administrative functionality. The traditional SDN operation is a centralized, global controller. While the centralized approach does strengthen the capability of the controller to manage the entire network, it is well known that a centralized approach also creates a vulnerability for a single point of failure. This is a significant barrier to using SDN for large-scale cyber-physical networks. New approaches must ensure that SDN controllers are fault tolerant on a larger scale and retain the advantages of the centralized perspective (global view) even while the implementation takes a distributed approach.

Related Work and Contributions of This Work

To this point, there has been limited research on SDN for monitoring of smart grid communications. Previously, SDN for Smart Grid has been proposed using centralized, non-real-time network monitoring and control. For example, the work of [START_REF] Lin | Self-healing attack-resilient pmu network for power system operation[END_REF], presents a self-healing PMU smart grid network using SDN. When a cyber-attack takes place on a PMU then that node is isolated and reconfigured to a previous stable state. These works do not consider cyber-attacks on the centralized SDN architecture itself nor the restoration and reconfiguration of the SDN controllers in a distributed SDN design. Very recently, distributed SDN (D-SDN) solutions have been proposed for specific categories of Smart Grid, including (PEVs) [START_REF] Chen | Sdn-based framework for the pev integrated smart grid[END_REF]. A few works proposed a distributed architecture for cybersecurity [START_REF] Pisharody | Brew: A security policy analysis framework for distributed sdn-based cloud environments[END_REF].The work in [START_REF] Nkosi | Classification of sdn distributed controller approaches: A brief overview[END_REF], provides detailed insights on various approaches for developing a distributed SDN architecture.

Our proposed HDD-SDN is a real-time, distributed approach that uses current information available at the PMUs and can respond in real-time to failures and attacks. HDD-SDN analyzes traffic flow as well as smart grid measurement data and use machine learning techniques to identify, fix, and then attempt to prevent similar cyber-attacks, while continuing to detect anomalous behavior in the future. Our work will leverage the physically distributed controller approach described in [START_REF] Nkosi | Classification of sdn distributed controller approaches: A brief overview[END_REF] and employs machine intelligence-based anomalous flow detection to increase resilience of smart grid systems as well as the underlying SDN communication network.

The rest of the paper is divided as follows. Section 2 introduces the hybrid distributed and decentralized controller software-defined network architecture for smart grid systems. Section 3 describes the execution of multiple anomaly detection algorithms to be able to maintain reliability while reducing latency and CPU/memory use. Section 4 provides the experimental setup, a combination of Smart Grid test bed and SDN network simulation. Section 4 also provides results on the machine learning clustering algorithms that were evaluated for anomaly detection in network traffic and generated PMU data. Finally, Section 5 provides conclusions and future work.

Network Architecture

As mentioned previously, Figure 1 shows an example application for the HDD-SDN architecture, namely the control and monitoring of PMU communications in a smart grid system. The problem of protecting and controlling the power grid is reduced into simpler, more tractable engineering problems by subdividing the power system into small regions or zones. For example, a fault in a transmission line or a fault in the control of the power output of a generator are problems solved locally by monitoring the variables measured by local sensors. The proliferation of the fault comes when data exchange begins between the local PMU and other PMUs.

The proposed HDD-SDN architecture is shown in Figure 2. The network consists of a collection of sub-regions, which can represent the PMU sub-systems within the smart grid architecture. In the hierarchical architecture, multiple controllers share their sub-region network topology and state information with their neighbors in the distributed system. This is accomplished using the distributed store function of the open-source network operating system ONOS [START_REF]ONOS-Open Network Operating System[END_REF].

The distributed store provides a global view of the entire networking system, and ensures sub-systems are performing as efficiently as possible. Each controller is responsible for managing the nodes under its sub-region/domain, and for updating important information from their sub-region to the distributed system. In this sense, the controllers are connected in a mesh using a specific TCP port for interactions and using keep-alive messages to monitor each controller's status. Consistency levels among the users can be "strongly consistent", which in ONOS indicates frequent updates of network topology state to the distributed stores using the RAFT protocol, or "eventually consistent", which implies less frequent updates using the Anti-Entropy protocol [START_REF] Muqaddas | Inter-controller traffic to support consistency in onos clusters[END_REF]. In the case of this project, we use an "eventually consistent" model for updating the distributed store of anomalous behavior taking place in a region. More information about the types of stores used in ONOS can be found here: [START_REF] Muqaddas | Inter-controller traffic to support consistency in onos clusters[END_REF]. Multiple controller support, introduced in OpenFlow 1.2 [START_REF] Kopeikin | Multi-uav network control through dynamic task allocation: Ensuring data-rate and bit-error-rate support[END_REF], allows a switch to connect to multiple controllers simultaneously. Each controller has either a MASTER role, EQUAL role, or a SLAVE role to each switch in the network. Next, we describe the decentralized approach that provides the SDN and Smart Grid network with additional benefits compared to a centralized strategy.

Decentralized SDN Architecture

The sub-region operation for the proposed HDD-SDN architecture is shown in Figure 3. Child controllers are placed in each sub-region to monitor data and control planes. This provides a redundant, decentralized SDN architecture for each sub-region. The primary tasks of the child controller are to intercept: (1) network data packets relayed in the data plane between PMU devices, and (2) network control packets transmitted from the parent controller to the PMU's switching devices. The child controller cannot make changes to the network it is monitoring. It is only allowed to communicate with the distributed store (hierarchical Fig. 4. Decentralized Sub-Regional Architecture Details distributed controller cluster) if its parent controller demonstrates anomalous behavior to the network. If this condition is detected, the child controller can flag the parent's anomalous behavior and relay the indication to its neighboring parent controllers. A neighbor parent controller can then take control of the affected sub-region, and the faulty controller can be reconfigured.

The proposed approach also allows the parent controller to offload security features, such as firewalls, deep packet inspection (DPI), or anomaly detection techniques to the child controller, since the logical decentralized child controller can be aware of all statistics in the network. Delegating jobs to separate controllers releases the burden placed on one controller, and allows room for numerous possibilities for managing or reconfiguring the network. Further details about the parent and child controller roles and operation are shown in Figure 4 and described below. We depend on network statistics from the entire system to identify anomalous behavior, including (1) anomalous behavior in the network traffic and payload content, i.e. control data sent between the SDN controller and PMU devices, network data with current and voltages reading exchanged between PMU devices and a state estimator, (2) anomalous behavior in topology changes, i.e. additional node changes within the network, such as new addition of external nodes and unwanted changes within network configurations; and (3) a variety of network performance changes, i.e broken links, lost connections, and performance decreases.

Parent Controller

Efficient resource allocation and management takes place on the parent controller to provide congestion control, load balancing, and traffic engineering for the smart grid and to distribute the computation load among the PMUs. The parent controller is also responsible for reconfiguring the network when the child controller demonstrates anomalous behavior in the system. The parent controller interacts with the devices in the region using the Openflow protocol for updat- ing flows, while simple network management protocol (SNMP) and NETCONF are used for gaining statistical information and reconfiguring devices when necessary. Each parent controller of each sub-region allocates enough memory to store the network states of their neighboring regions, allowing for each controller to monitor the network performance of the neighboring regions.

Child Controller

As stated previously, the child controller has no ability to change topological information of the region. Its sole purpose is monitoring the status and performance of the parent controller and the devices that form the sub-regions network using management protocols like SNMP, and NETCONF for retrieving information of the system. As indicated in Figure 2, the child controller contains a few security applications specifically for determining anomalies of the packets sent in the control or data plane. The anomaly detection (AD) module uses a combination of machine learning algorithms for clustering different types of network traffic and payload data, and an identification dictionary (IDD) module. The decentralized architecture grants use of multiple anomaly detection algorithms in parallel, without greatly increasing the CPU and memory usage, as it would if management and anomaly detection were executed on a centralized controller.

Feature Scaling and Preprocessing

The application flow of the child controller is described in Figure 5. Preprocessing features is a standard step for implementing machine learning algorithms for anomaly detection [START_REF] Limthong | Real-time computer network anomaly detection using machine learning techniques[END_REF]. However, while most techniques attempt to normalized a wide range of different features into a standard range of 0 to 1 [START_REF] Limthong | Real-time computer network anomaly detection using machine learning techniques[END_REF], in this approach, we increased the scaling factor by 1000 while incrementing through the columns. This granted a substantial change in the calculation of centers when certain features are present and allowed the clustering algorithm to better distinguish between different cluster groups.

Anomaly Detection in the HDD-SDN Child Controller

Various machine learning techniques can be implemented in parallel to detect anomalous behavior. A commonly used clustering algorithm, k-means, is implemented, along side a sparse clustering extension called incremental k-means [START_REF] Yadav | Incremental k-means clustering algorithms: A review[END_REF]. The k-means algorithm clusters data points by alternatively assigning data points to clusters and updating cluster representatives. Data points are assigned to the cluster in which they have the minimum Euclidean distance to cluster centers computed using Equation 1:

d(p n , c k ) = D d=1 (p n -c k ) 2 1/2 (1)
where p n = [p 1 , p 2 , . . . , p D ] ∈ R D is a D dimensional vector containing the features associated with the n th data point (in our case, a network packet) and

c k = [c 1 , c 2 , . . . , d D ] ∈ R D is a D dimensional
vector of the k th cluster representative. Cluster representatives are using the vector mean of the data points assigned to that cluster. Since the k-means algorithm requires repeated pair-wise distance computations between each data point and each cluster representative, the computational load is significant for large data sets. Incremental k-means is a k-means approximation that is applicable to large scale sparse dynamic datasets (such as network traffic). During each iteration, the standard k-means algorithm uses the entire dataset for recalculating cluster centers. In contrast, Incremental k-means updates the previous centers with only newly input data [START_REF] Yadav | Incremental k-means clustering algorithms: A review[END_REF].

The k-means approach and its extensions have an obvious limitation in that the number of clusters must be known before executing the algorithm. However, the number of clusters are rarely known in advance, particularly for dynamic datasets containing anomalies. Cluster validity metrics are one mechanism to address this issue. Cluster validity metrics provide a quantitative measure of clustering effectiveness for a particular data set. Thus, a data set can be repeatedly clustered with a different number of clusters and, then, each result evaluated using a validity metric to determine the appropriate number of clusters. Validity metrics appropriate for the k-means algorithm include Dunn's index and the Davies-Bouldin index [START_REF] Dunn | Well-seperated clusters and optimal fuzzy partitions[END_REF][START_REF] David | A cluster seperation measure[END_REF]. A draw-back of this potential approach is the need to repeatedly cluster the data using a different number of clusters each time. Alternatively, the number of clusters can be adapted in real time with incoming data by generating new clusters when data points appear that are far from all current cluster representatives (similar to approaches used to generate new clusters in Dirichlet Process Mixture Modeling) [START_REF] Neal | Markov chain sampling methods for dirichlet process mixture models[END_REF]. These new clusters could then be evaluated and tagged as corresponding to anomalous behaviors when appropriate. When anomalous behavior is identified, short-term mitigation techniques can take place for a quick response to the problem. For example, a firewall can be activated on an infected IP address.

Identification and Detection Dictionary in the HDD-SDN Child Controller Figure 4 also shows the Identification and Detection Dictionary (IDD) process in the Child Controller. The IDD uses a simple logical algorithm for comparing extracted features from the anomalous traffic to features associated with different types of cyber-attacks as identified during a training phase. Different types of cyber-attacks affect different features of network traffic. If the anomalous behavior is completely different from what is stored and cannot be identified, then the system has learned a new type of cyber-attack or a new type of fault in the system. This information will be recorded and classified as a new anomaly. Then it is sent to the distributed store for the other regions to be able to take action in mitigating that type of attack or fault, or tagging reported data with the current attack or fault. After an anomaly is identified, the ONOS system can now produce long-term changes in the network configuration to prevent future attacks of this type,such as changing flows, reconfiguring the attacked device, rerouting around the offending sub-region, etc.

The extracted features from the anomalous traffic in a sub-region sometimes cannot be directly used in another sub-region. Different regions may have different network topology and configurations. An abnormal traffic pattern in one region may be normal in another region. We envision the child controller to be able to develop mitigation schemes for their respective sub-domains in which they reside. For instance if anomalous behavior is detected in sub-region 1, the child controller (monitor controller) in that region will develop a mitigation scheme for its sub-region only. This information is transmitted to it's respective parent controller to take action. The child controller will then relay information about the type of attac, what network features were affected, feature threshold values to be used by other child controllers in separate sub-regions. It is up to the other sub-regions to monitor their network features and determine if they are being changed in a similar manner. If their thresholds aren't met then no change will take place for the other regions. If their thresholds are met then the child controllers of those regions will develop their own mitigation scheme to correct the anomalous behavior.

Experimental Results

The proposed HDD-SDN environment was simulated in Mininet using the wireless devices environment. The SDN controller used was the open-source network operating system ONOS [START_REF]ONOS-Open Network Operating System[END_REF]. ONOS provides a large API for users to develop their own networking applications to meet requirements of custom networking scenarios.

Simulation Setup

Using Mininet, we developed networking scenarios and topologies similar to realistic smart grid environments. In addition, we deployed an instant virtual network on a stand-alone computer and were able to expand this network by allowing the connection of multiple external nodes and other computation resources, including other PCs, mobile devices, VMs, etc. In this environment the state estimator was the client node connected to multiple PMUs setup as servers or client nodes.

We are developing our testbed to integrate a smart grid network emulator called OPAL-RT to generate the smart grid information. Virtual measurements of the smart grid emulator are taken from OPAL-RT, streamed from the server to the clients, and analyzed by the HDD-SDN application in ONOS. Network state information (i.e. number of devices, number of links, flow information, network traffic, etc.) collected by the external SDN parent and child controllers were stored in a database (InfluxDB, AWS, etc.) to be analyzed by an instance AD module operating in the child controller. Any necessary network reconfigurations were sent back to the parent SDN controller for the proper actions to take place. For experimentation purposes, the cluster algorithms were trained using the KDD CUP 1999 dataset that contains a standard set of data, which includes a wide variety of classified intrusions emulated in a military network environment [START_REF]KDD Cup 1999 Data[END_REF]. With this data we extracted key features of network traffic identical to features which can be obtained using ONOS controller during network monitoring to compare anomaly detection techniques. We simulated measurements of the recorded and exchanged measurements by PMUs using MATLAB.

Performance Evaluation

For this experiment, we executed both k-means and incremental k-means for clustering and, subsequently, anomaly detection on network traffic data and PMU measurement data. In order to evaluate how well the network traffic was classified, the data was separated into training and testing sets. The training set is used to generate cluster centers for the k-means and incremental k-means algorithms. Once the cluster centers are calculated, a label (normal or anomaly) is manually assigned to each cluster based on the IDD approach described above. The testing set is then added to the dataset and assigned the label of the cluster center closest to that data point. The size of the training set was n = 194,565 (70 percent) and the size of the testing set was n = 83385 (30 percent). The results of the classification of the network traffic are shown in Figure 6. We recorded true positive (TP), true negative (TN), false positive (FP), and false negative (FN) results. Using these values we calculated true positive rate, T P R (sensitivity), true negative rate, T N R (specificity), positive predictive value P P V (precision), negative predictive value, N P R, false positive rate F P R (fallout), false negative rate, F N R, false discovery rate F DR, and overall accuracy, ACC, with the inclusion of execution time, CPU and memory usage as follows:

T P R = T P/(T N + F P )

T N R = T N/(T N + F P ) P P V = T P/(T P + F P )

N P V = T N/(T N + F N ) F P R = F P/(F P + T N ) F N R = F N/(T P + F N ) F DR = F P/(T P + F P ) ACC = (T P + T N )/(T P + F P + F N + T N ) (2)
As shown in Figure 6, we observed that k-means produces a high accuracy (96.9 percent) of identifying anomalies within network traffic. Incremental kmeans produces a slightly lower accuracy (95.6 percent), but demonstrated an increased speed with respect to the k-means by approximately 22 seconds. For the CPU and memory usage, k-means used more resources than incremental kmeans. This is likely due to the fact that kmeans is a batch algorithm, using the entire dataset during every iteration to recalculate cluster centers. On the other hand, it was observed that the accuracy of anomaly detection for the OPAL-RT generated PMU data is lower than the accuracies provided when clustering network traffic. This is possibly because features of the PMU data are not as readily distinguishable from each other as expected, causing a higher false positive rating than expected. The proposed anomaly detection technique using clustering still preserves approximately 90 percent accuracy, mostly due to the true negative value exceeding others. It can be seen that CPU and memory usage are much lower for this data type, due to use of fewer features when clustering. This is a positive indicator for scaling the system to larger PMU networks.

For performance overhead, we implemented a single ONOS controller with network management applications running on a device with the following specs: Ubuntu OS 14.04 (64-bit), AMD A6-6310 APU processor, 3.3 GiB memory, Java version 1.7.0 to manage a simple network consisting of 5 switches and 10 hosts. CPU usage for this system ranged from 1.3-17.3% of the device. We predict adding an additional controller on the same device to manage the same network will double this usage. With the ONOS system we can deploy separate controllers as VMs operating on different servers if need be. The research in [START_REF] Muqaddas | Inter-controller traffic to support consistency in onos clusters[END_REF] highlights the amount of increased bandwidth between decentralized controllers communication as the number of nodes controlled is increased.

Conclusion and Future Work

In conclusion, the benefits of implementing an adaptive distributed and decentralized SDN in place of the common networking or in place of traditional SDN strategy has been discussed. It was shown that implementing the HDD-SDN architecture can provide safe fail-over using redundant systems and additional resilience in the presence of faulty or attacked data or communication nodes in the smart grid system. The paper also evaluated the use of a combination of machine learning clustering algorithms for parallel processing of anomaly detection and discussed potential approaches for automated determination of the number of clusters needed. K-means produced a high accuracy for identifying anomalies within normal traffic and incremental k-means produced a slightly lower accuracy with increased speed and fewer CPU and memory resources, indicating a possibility for scaling the system to much larger networks. In future work, the process for preprocessing of the features before executing the machine learning will be examined. In addition, we continue to develop the SDN architecture integration with the OPAL-RT Smart Grid test bed.

Fig. 1 .

 1 Fig. 1. Distributed PMUs in the Smart Grid. The distributed PMUs form a communication network, where smart grid data is regularly exchanged with other PMUs to coordinate and analyze energy performance measures. If a PMU subsystem is attacked or fails, the communication of invalid data may introduce substantial errors when exchanged with the other connected PMU subsystems.

Fig. 2 .Fig. 3 .

 23 Fig. 2. Distributed SDN Architecture and Store

Fig. 5 .

 5 Fig. 5. Example Flow Diagram of Applications in Child Controller

Fig. 6 .

 6 Fig. 6. Kmeans and Incremental Kmeans Confusion Matrix for Network Traffic

Fig. 7 .

 7 Fig. 7.Table Showing Accuracy of Anomaly Detection Algorithms for Network Traffic and PMU Data

  Table Showing Accuracy of Anomaly Detection Algorithms for Network Traffic and PMU Data

Acknowledgment

The authors would like to thank the Harris Corporation Excellence in Research program for providing funding for this research.


Pierre Dissaux

Frank Singhoff

Stood and Cheddar: AADL as a Pivot Language for Analysing Performances of Real Time Architectures

Keywords: AADL, Real-Time, Performance analysis, Stood, Cheddar, Design Patterns

Coupling of modelling and analysis tools requires that both ends strictly comply with the same semantic definition of the exchanged model. This is particularly important for real-time systems and software architectures. Such a guaranty can be brought by the common use of the Architecture Analysis and Design Language (AADL) all along the tool-chain. This paper discusses modelling and analysis options of various real-time architectural patterns expressed in AADL though an experiment with Stood and Cheddar tools.

Introduction

The SAE Architecture Analysis and Design Language (AADL) is a textual and graphical language support for model-based engineering of embedded real time systems that has been approved and published as SAE Standard AS-5506 [1]. AADL is used to design and analyse software and hardware architecture of embedded real-time systems. Many tools provide support for the modelling and the analysis of AADL models. Ocarina implements Ada and C code generators for distributed systems [START_REF] Hugues | Rapid Prototyping of Distributed Real-Time Embedded Systems Using the AADL and Ocarina[END_REF]. TOPCASED, OSATE and Stood provide AADL modelling features [START_REF]TOPCASED web site[END_REF][START_REF] Dissaux | Using AADL for mission critical software development[END_REF][START_REF] Sei | OSATE : an extensible Source AADL tool environment[END_REF]. The Fremont toolset and Cheddar implement AADL performance analysis methods [START_REF] Sokolsky | Schedulability analysis of AADL models[END_REF][START_REF] Singhoff | Scheduling and Memory requirement analysis with AADL[END_REF][START_REF] Singhoff | AADL Modeling and Analysis of Hierarchical Schedulers[END_REF]]. An updated list of supporting tools can be found on the official AADL web site: http://www.aadl.info. This article deals with the interoperability between AADL tools: we show how AADL can be used as a pivot language between a modelling tool (Stood) and a performance analysis tool (Cheddar). Stood is a software design tool that provides an extended support of the AADL modelling language in addition to its compliancy with the HOOD methodology. With Stood, it is possible to manage a complete software project by building libraries of reusable components, reversing legacy code and specifying the real time application as well as its execution platform. Most of the modelling activities can be performed graphically and the corresponding AADL code is automatically generated by the tool.

The Cheddar framework is a set of Ada packages which aims at performing performance analysis of real time architectures. It includes analytical scheduling methods and most of classical real time scheduling algorithms. The Cheddar framework also offers a domain specific language (and its interpreter, compiler, …) for the design and the analysis of schedulers which are not already implemented into the framework.

In this article, in order to illustrate the interoperability between Stood and Cheddar, we propose a set of AADL design patterns to model usual real time synchronization paradigms [START_REF] Tanenbaum | Modern Operating Systems[END_REF]. This paper is organized as follows: In section 2, we present performance analysis methods that are expected to be applied on AADL design patterns. These AADL design patterns are then described in section 3. Finally, we conclude and describe ongoing works in section 4.

Real-Time performance analysis with AADL

2.1 Real time scheduling analysis and Queueing system analysis From an AADL model, we can perform performance analysis based on real time scheduling theory and queueing system theory. Real time scheduling theory helps the system designer to analyse the timing behaviour of a set of tasks with scheduling algorithms or with algebraic methods usually called feasibility tests. For example, with the well known Liu and Layland real time task model [START_REF] Liu | Scheduling algorithms for multiprogramming in a hard real-time environnment[END_REF], each task periodically performs a treatment. This "periodic" task is defined by three parameters: its deadline (Di), its period (Pi) and its capacity (Ci). Pi is a fixed delay between two release times of the task i. Each time the task i is released, it has to do a job whose execution time is bounded by Ci units of time. This job has to be ended before Di units of time after the task wake up time. Some algebraic methods can provide a proof that an architecture will meet its periodic task performance requirements. Scheduling algorithms allow the designer to compute scheduling simulations of the architecture to analyse. Usually, simulations can not lead to a proof. However, in some cases (with deterministic schedulers and with periodic tasks for example), scheduling simulation may lead to a schedulability proof if the designer is able to compute a scheduling during the base period [START_REF] Leung | On the Complexity of Fixed-Priority Scheduling of Periodic, Real-Time Tasks[END_REF]. Different kinds of feasibility tests exist: tests based on processor utilization factor and tests based on worst case task response time which are designed to check task deadlines [START_REF] George | Preemptive and Non-Preemptive Real Time Uni-Processor Scheduling[END_REF]; tests based on buffer utilization factor which are designed to check buffer overflow [START_REF] Singhoff | Scheduling and Memory requirement analysis with AADL[END_REF]. For example, the worst case response time feasibility test consist in comparing the worst case response time of each task with its deadline. Joseph, Pandia, Audsley et al. [START_REF] George | Preemptive and Non-Preemptive Real Time Uni-Processor Scheduling[END_REF] have proposed a way to compute the worst case response time of a task with pre-emptive fixed priority scheduling by:

[1]
Where Ri is the worst case response time of the task i. This feasibility test can be easily extended to take into account task waiting time on shared resources, jitter on task release time, task precedency relationship …

Queueing system theory may be used to perform analysis of real time architecture when waiting times exist. Queueing system theory allows to study performance of a system composed of servers, clients and storage places [START_REF] Kleinrock | Computer Applications[END_REF]: people waiting in a room for a doctor, network switch routing data, ... If clients becomes active while a server is busy, their requests are stored in a queue. By defining the average rate of clients request arrivals and the average rate of requests that the server can handle, a queueing system model allows to predict the average system occupation factor L, the average customer waiting time W, and the probability Pn of having n clients in the queue.

Investigated performance criteria

An AADL model is a set of hardware and software components. An AADL operational system is a set of process components encompassing thread and data components that are bound to an execution platform composed of processor, memory and bus components. Component relationships are modelled by ports and data access connections. Three different kinds of ports exist: data ports, event ports and event data ports. Data ports represent connection points for transfer of data values. Event ports represent connection points for transfer of control through raised events that can trigger thread dispatch or mode transition. Event data ports represent connection points for transfer of events with data, i.e., messages that may be queued. Data access represent asynchronous read or write operations on a shared data component.

From an AADL model, several performance criteria can be computed with the algebraic methods proposed by the real time scheduling theory and by the queueing systems theory. Some examples of these performance criteria are:

A. The worst case task response times; B. The bounds on the thread waiting time due to data access; C. The deadlocks and priority inversions due to data access; D. The numbers of messages in queued message communication links. This criterion allows memory footprint analysis; E. The numbers of context switches, preemption; F. The processor utilization factor; G. … In the sequel, we only focus on the 4 th first criteria.

Examples of AADL design patterns: from modelling to analysis

In the next sections we present four design patterns that can be used to express usual inter-thread communications.

For each pattern, an applicative test case is described under the form of an AADL model which has been formatted in purpose to highlight some of the possible performance analysis presented in the previous section. These criteria can not be investigated independently: in the case of the Blackboard design pattern for instance, one must compute first B before computing A.

Modelling pattern Analysis criteria

j i hp j j i i i C P R C R . ) ( ∑ ∈         + =
its output data ports at complete time. This design pattern does not require the use of a shared data component. In this simple case, the execution platform consists in one processor running a scheduler that is compliant with fixed priority scheduling such as Rate Monotonic [START_REF] Liu | Scheduling algorithms for multiprogramming in a hard real-time environnment[END_REF].

Example: Corresponding modelling case study consists in three periodic threads linked together by two data port connections. The AADL graphical representation of this case study and fragments of its textual specification are given below: Analysis: This design pattern leads to a very static scheduling which is difficult to change but which is also very easy to analyse. Indeed, in this design pattern, threads are independent, which allows the use of numerous simple feasibility tests such as the processor utilization factor test or the worst case response time. These simple feasibility tests are available for most usual real time schedulers such a EDF, LLF and fixed priority schedulers. Verifying this design pattern only consists in checking that thread deadlines will be met. Since no shared data component is used in this design pattern to implement data exchange, no thread blocking time on shared resource has to be computed.

Mutex protected communication pattern.

Main drawback of the previous pattern is its lack of flexibility at run time: each thread will always execute, read and write data at pre-defined times, even if useless. In order to introduce more flexibility at that level, asynchronous inter-thread communications must be considered. An example of such a run-time environment is given by the Ravenscar profile. Description: Ravenscar is a part of Ada 2005 [5,[START_REF] Burns | Concurrent and Real Time programming in Ada[END_REF]. Ravenscar is a set of Ada program restrictions usually enforced at compilation time, which allows a more flexible design of the architecture and which guaranties that this architecture remains compliant with real time scheduling theory/queueing systems theory analysis assumptions.

Ravenscar is an Ada subset from which one can write applications composed of a set of threads. In Ravenscar, threads access shared data components asynchronously according to priority inheritance protocols. Ravenscar assumes that threads are scheduled with a fixed priority scheduler and that data are accessed with the ICPP protocol [START_REF] Burns | Concurrent and Real Time programming in Ada[END_REF]. However, it may be possible to apply Ravenscar to some dynamic schedulers such as EDF if a data access protocol exists for such schedulers (eg. PLCP protocol for EDF scheduler [5]). Ravenscar is then a subset of concurrency features which can be defined in many real time executive such as POSIX 1003.1b [START_REF] Gallmeister | POSIX 4 : Programming for the Real World[END_REF], ARINC 653 [START_REF][END_REF] and Java-RT [START_REF]Specification: JSR001 Realtime Specification for Java[END_REF].

Exemple: To illustrate asynchronous inter-thread communications, the chosen case study shows an implementation of the classical P and V procedures of a Mutex. A complete description of this protocol in AADL requires the use of the Behavior Annex [START_REF]AADL Behavior Annex[END_REF] to provide a detailed specification of the internal realisation of the P and V subprograms. Analysis: With this design pattern, additional analysis such as deadlock detection, can be performed. As threads are not independent anymore, response time analysis becomes more complicated to investigate. Verifying thread deadline with this design pattern requires to first compute the thread waiting time before computing the worst case thread response time. This waiting time is computed according to the duration of the critical sections and the shared resource access protocol (eg. Priority Ceiling Protocol [START_REF] Sha | Priority Inheritance Protocols : An Approach to Real Time Synchronization[END_REF]).

SUBPROGRAM IMPLEMENTATION

Blackboard communication pattern

Description: Ravenscar allows a thread to allocate/release several shared resources (eg. AADL data). Real time scheduling theory usually models such a shared resource as a semaphore, to represent, for example, a critical section. In classical operating system, it exists many synchronization design patterns such as critical section, barrier, Readers-Writer, private semaphore, and various Producer-Consummer synchronization design patterns [START_REF] Tanenbaum | Modern Operating Systems[END_REF]. Programming languages also propose specific synchronization design patterns such as Java synchronized methods [START_REF] Zakhour | The Java Tutorial: A Short Course on the Basics, 4th Edition (The Java Series)[END_REF] or Ada protected types [5]. The Blackboard communication is one of them. The Blackboard design pattern implements a Readers-Writer synchronization protocol. At a given time, only one writer can get the access to the Blackboard in order to update the stored data, as opposed to the readers which are allowed to read the data simultaneously. The usual implementation of this protocol implies that readers and writers do not perform the same semaphore access.

Example: The following case study shows the AADL implementation of a multi-threads application with asynchronous communications through a Blackboard using the Readers-Writer protocol. Analysis: With this kind of design pattern, we can expect to compute the same performance criteria that we have shown for the Mutex design pattern. However, with complex synchronization design patterns, before computing worst case thread response time, we have first to analyse shared resource blocking time. For example, with the Blackboard design pattern, a thread may get the access to a different number of semaphores depending on its type (eg. reader or writer thread). It means that we have to first evaluate the semaphore access of each thread in order to compute shared data waiting time, for a given critical section and for a given priority inheritance protocol.

Queued buffer communication pattern

Description: In the Blackboard design pattern, at any time, only the last written message is made available to the threads. Some real time executives (eg. ARINC 653) provide communication features which allow to store all written messages in a memory unit. AADL also propose such a feature with event data ports or shared data components. An event data port may contain several messages. We assume that buffer messages are handled with a FIFO protocol.

Example: The case of the Producer-Consumer can be used to illustrate the queued buffer communication pattern. Like the Blackboard, the buffer component can be described in AADL by its provided services and its internal variables. Analysis: This design pattern needs the same schedulability analysis as the previous one. An other typical analysis for such design pattern is also to estimate the memory footprint of the buffer in order to ensure no loss of data when the rate of the producer thread is temporarily greater than the one of the consumer thread. In [START_REF] Legrand | Contribution à l'ordonnancement des systèmes temps réel comprenant des tampons[END_REF][START_REF] Legrand | Performance Analysis of Buffers Shared by Independent Periodic Tasks[END_REF], Legrand et al. have proposed a set of feasibility tests based on queueing system. These feasibility tests were adapted to AADL in [START_REF] Singhoff | Scheduling and Memory requirement analysis with AADL[END_REF]. It was shown how perform memory footprint analysis with AADL models containing event data ports. For example, if both producers and consumers are periodic thread, a worst case number of messages in a buffer can be computed with this feasibility test: L=2.n if threads are harmonic [START_REF] Liu | Scheduling algorithms for multiprogramming in a hard real-time environnment[END_REF] or

DATA IMPLEMENTATION

L=2.n+1 otherwise

Where L is the maximum number of messages and n is the number of producers. In the case of sporadic threads, the same worst case analysis can be performed, but [START_REF] Legrand | Performance Analysis of Buffers Shared by Independent Periodic Tasks[END_REF] also have proposed feasibility tests for average analysis.

Conclusions and perspectives

This article deals with the interoperability between a modelling tool called Stood and an analysis tool called Cheddar. Coupling of modelling and analysis tools requires that both ends strictly comply with the same semantic definition of the exchanged model. For such a purpose, we have chosen AADL as a pivot language between Stood and Cheddar. AADL is a textual and graphical language support for model-based engineering of embedded real time systems. AADL is a flexible language which allows the modelling of both synchronous and asynchronous systems. Then, in order to illustrate interoperability between Stood and Cheddar, we have proposed a set of AADL design patterns. For each design pattern, we have listed a set of performance criteria that can be checked. Each design pattern is also illustrated by an example. This article presents a first step to achieve interoperability between Stood and Cheddar. The next steps must refine the set of AADL design pattern presented above. Other synchronization patterns must also be investigated (eg. Private semaphore, Ada protected type, …).

Figure 1 :

 1 Figure 1: synchronous data-flows pattern Real-time attributes of each thread can be described individually by a set of pre-defined properties: THREAD IMPLEMENTATION T1.others PROPERTIES Dispatch_Protocol => Periodic; Period => 29ms; Deadline => 29ms; Compute_Execution_Time => 7ms..7ms; END T1.others;

Figure 2 :

 2 Figure 2: analysis performed by Cheddar We can note however, that this kind of analysis still depends on the validity of execution time estimates.Figure 2 depicts an analysis of this test case by Cheddar. The top part of the window shows the thread scheduling. The worst case response time computed from this scheduling and with feasibility tests are shown in the bottom part of the window.

Figure 3 :

 3 Figure 3: mutex protected communication pattern The use of the AADL Behavior Annex is required again to express the internal functional structure of the applicative threads. Such a level of details is needed to describe the atomic actions and the critical sections. THREAD IMPLEMENTATION Thread_A.others PROPERTIES Dispatch_Protocol => Periodic; Period => 10 ms; ANNEX Behavior_Specification {** states s0: initial state; s1, s2, s3, s4: state; s5: complete state; transitions acquire_M1: s0-[]->s1 {P!(Mutex_1);}; acquire_M2: s1-[]->s2 {P!(Mutex_2);}; critical_section: s2-[]->s3 {…}; release_M1: s3-[]->s4 {V!(Mutex_1);}; release_M2: s4-[]->s5 {V!(Mutex_2);}; **}; END Thread_A.others; THREAD IMPLEMENTATION Thread_B.others … ANNEX Behavior_Specification {** … transitions acquire_M2: s0-[]->s1 {P!(Mutex_2);}; acquire_M1: s1-[]->s2 {P!(Mutex_1);}; … **}; END Thread_B.others;

Figure 4 :

 4 Figure 4: Blackboard communication patternThe communication object is an instance of an AADL data component that can be described by its provided services, implemented by subprogram features, and its internal structure consisting in a set of state variables.

Figure 5 :

 5 Figure 5: Queued buffer communication pattern

DATA IMPLEMENTATION T_BlackBoard.o SUBCOMPONENTS Contents: DATA T_Item; Readers: DATA Behavior::Integer; Is_Idle: DATA Behavior::Boolean; Is_Reading: DATA Behavior::Boolean; Is_Writing: DATA Behavior::Boolean; END T_BlackBoard.o;


A Bersani

F Dell'isola

P Seppecher

Lagrange multipliers in infinite dimensional spaces, examples of application

Keywords:

The Lagrange multipliers method is used in Mathematical Analysis, in Mechanics, in Economics and in several other fields, to deal with the search of the global maximum or minimum of a function, in presence of a constraint. The usual technique, applied to the case of finite-dimensional systems, transforms the constrained optimization problem into an unconstrained one, by means of the introduction of one or more multipliers and of a suitable Lagrangian function,to be optimized. In Mechanics, several optimization problems can be applied to infinite-dimensional systems. Lagrange multipliers method can be applied also to these cases.

Introduction

In this entry we show that the theorem of Lagrange multipliers in infinite dimensional systems [START_REF] Dell'isola | Lagrange Multipliers in Infinite-Dimensional Systems, Methods of[END_REF] can be a very powerful tool for dealing with constrained problems also in infinite dimensional spaces. This tool is powerful but must be used carefully. As penalization is often invoked as an intuitive and numerically efficient approach to constrained problems, we also show that the penalization approach may present the same drawbacks as a rough application of Lagrange multipliers method does. To that aim we describe two examples issued from continuum mechanics [START_REF] Glüge | Continuum Mechanics Basics, Introduction and Notations[END_REF].

The first one goes back to Lagrange himself [START_REF] Lagrange | [END_REF]: an incompressible fluid, or an incompressible linear elastic material is a material whose displacement field is constrained to be divergence-free. We show how easily the Lagrange multipliers method gives the system of balance equations and introduces the pressure as a supplementary unknown of the problem. The approach through penalization gives the same balance equations and the displacement field becomes divergence-free at the limit only. It is not very difficult to understand why the pressure which is not independent from the displacement until the limit, becomes an independent quantity at the limit. Hence both approaches are efficient for understanding what an incompressible material is and how it behaves.

The second example requires us to be more cautious. It consists in considering that second gradient materials are nothing else than micromorphic materials subjected to the constraint that the micro-deformation coincides with the gradient of the displacement field. We show that a careful application of the Lagrange multipliers method leads to the correct system of equilibrium equations and we warn against its possible erroneous applications. We show that the penalization method, when over-interpreted, may lead to the same errors: one cannot infer from the admissible boundary conditions for a micromorphic material what can be the boundary conditions for the constrained model, that is for the second gradient model. The discrepancy between results obtained by studying directly the critical points of a second gradient energy and those obtained by a rough limit of micromorphic models has led to some confusion which still survives.

Incompressible materials

Consider an elastic material contained in a bounded Lipschitz domain Ω. Its elastic energy density depends on the first gradient of the displacement field u. The global energy has the form1

L(u) := Ω (u(x), ∇u(x)) dH 3
The material is said to be linearly incompressible2 if u satisfies almost everywhere in Ω the condition div(u) = 0.

Assuming a quadratic growth at infinity for L with respect to its second argument, the natural functional framework is the Sobolev space H 1 (Ω, R 3 ). In the sequel L 2 and H 1 stand for the Sobolev spaces L 2 (Ω, R 3 ) and H 1 (Ω, R 3 ).

The goal is to write the partial differential equation that any smooth critical point u of the global energy must solve. We thus assume that both u and are of class C 2 .

Application of Lagrange Multiplier Theorem

It is well known that the divergence operator is surjective from H 1 onto L 2 . Indeed, for any f ∈ L 2 , it is enough to remind that the Dirichlet problem ∆v = f on Ω with v = 0 on the boundary admits a solution v and thus that f = div(∇v). Duality in the Hilbert space L 2 is trivial: Lagrange multiplier theorem states the existence of a L 2 function -p such that any critical point u of L under the constraint [START_REF] Dell'isola | Lagrange Multipliers in Infinite-Dimensional Systems, Methods of[END_REF] satisfies

Ω (∂ 1 • h + ∂ 2 • ∇h -p • div(h)) dH 3 = 0. ( 2
)
or

Ω (∂ 1 • h + (∂ 2 -p Id) • ∇h) dH 3 = 0
for any smooth test function h. Here the partial derivatives of are taken at (u(x), ∇u(x)).

Hence, in the sense of distributions on Ω,

∂ 1 H 3 |Ω -div((∂ 2 -p Id) H 3 |Ω ) = 0. Now let us use the divergence theorem 3 : div(ϕ H 3 |Ω ) = div(ϕ) H 3 |Ω -(n • ϕ) H 2 |∂Ω (3)
where n stands for the outward normal to ∂Ω and H 2 |∂Ω stands for the restriction to the boundary ∂Ω of the two-dimensional Hausdorff measure. We get

(∂ 1 -div((∂ 2 -p Id)) H 3 |Ω + n • (∂ 2 -p Id) H 2 |∂Ω = 0.
In mechanics the quantities -∂ 1 and ∂ 2 -p Id are respectively denoted f and σ and interpreted as the bulk external force and the stress tensor. Using this notation, previous equation reads

-(f + div(σ)) H 3 |Ω + n • σ H 2 |∂Ω = 0.
As the measures H 3 |Ω , H 2 |∂Ω are orthogonal (or mutually singular) [START_REF] Rudin | Real and complex analysis[END_REF], this equation splits in

f + div(σ) = 0 H 3 a.e. in Ω, n • σ = 0 H 2 a.e. on ∂Ω.
We recover the standard equilibrium equations where the stress tensor σ = ∂ 2 -p Id involves explicitly the Lagrange multiplier -p. Its opposite p is interpreted as an unknown pressure which has to be determined together with the equilibrium displacement u by using the equilibrium equations together with the constraint div(u) = 0.

3 Indeed, for any C 1 function ϕ and any smooth test function ψ,

div(ϕ H 3 |Ω ), ψ = -ϕ H 3 |Ω , ∇ψ = - Ω ϕ • ∇ψ dH 3 = Ω div(ϕ) • ψ dH 3 - ∂Ω (n • ϕ)ψ dH 2 .
3

Approach through penalization

An intuitive way for dealing with incompressibility is to penalize the fact that div(u) does not vanish. Looking for minimizers of the total energy, one decides to look for a sequence of approximate minimizers by considering, instead of the original potential subjected to the constraint, the original Lagrangian to which one adds a penalization term: for instance

L ε (u) = Ω ( (u(x), ∇u(x)) + ε -1 (div(u)) 2 ) dH 3 . (4)
This corresponds to a slightly compressible elastic material. Denoting -p ε := 2ε -1 div(u), any critical u satisfies, for any smooth test function h,

Ω (∂ 1 • h + ∂ 2 • ∇h -p ε div(h)) dH 3 .
This equation is identical to [START_REF] Glüge | Continuum Mechanics Basics, Introduction and Notations[END_REF]. It then leads to the same system of equilibrium equations. The difference is that now the stress tensor is directly related to u through the so-called constitutive equation σ = ∂ 2 l + 2ε -1 div(u) Id. When passing to the limit ε → 0, div(u) tends to zero, the system of equilibrium equations is preserved but the constitutive equation is partially lost as ε -1 div(u) is an undetermined form.

Second gradient as constrained generalized model

Consider an elastic material contained in a bounded simply connected domain Ω with piecewise C 1 boundary. In the sequel L 2 and H s stand for the Sobolev spaces L 2 (Ω, R 3 ) and H s (Ω, R 3 ). The elastic material is said to be a second gradient material (or equivalently a strain-gradient one) if its elastic energy density depends on the gradient and the second gradient of the displacement field u. The global energy has the form Ω l(u(x), ∇u(x), ∇∇u(x)) dH 3 Such a model and particularly the associated boundary conditions present some difficulties of interpretation (see for instance [START_REF] Dell'isola | How contact interactions may depend on the shape of Cauchy cuts in Nth gradient continua: approach "à la D'Alembert[END_REF][START_REF] Isola | Cauchy tetrahedron argument applied to higher contact interactions[END_REF][START_REF] Corte | The postulations à la D'Alembert and à la Cauchy for higher gradient continuum theories are equivalent: a review of existing results[END_REF]). That is why many researchers prefer to see this model as a special case of a generalized continuum [START_REF] Forest | First vs. second gradient of strain theory for capillarity effects in an elastic fluid at small length scales[END_REF].

An elastic material is said to be a micromorphic generalized continuum if its elastic energy density depends on an extra matrix valued kinematic quantity v and on the first gradients of u and v. The global elastic energy has the form

L(u, v) := Ω (u(x), v(x), ∇u(x), ∇v(x)) dH 3
It is clear that, when such a model is constrained in such a way that v has to coincide with the first gradient of u, that is if

C(u, v) := v -∇u = 0,
then one recovers the previous second gradient model.

This indirect approach does coincide with the direct approach under the express condition that the Lagrange Multiplier Theorem is precisely respected. A rough application of Lagrange Multiplier Method may lead and have led to misunderstanding boundary conditions.

In the sequel we assume that the considered energies have quadratic growth at infinity so the Sobolev spaces L 2 (Ω, R 3 ), H s (Ω, R 3 ) or H s (Ω, R 3×3 ) are the natural framework. We simply denote L 2 and H s these spaces as the tensorial nature of the involved functions is clear from the context.

We are thus looking for smooth critical points of the functional L(u, v) defined on H 1 × H 1 under the constraint C(u, v) = 0. It is obvious that this problem can be simplified by setting4 l(u, η, ξ) := (u, η, η, ξ) and by looking for the critical points of

L(u, v) := Ω l(u(x), ∇u(x), ∇v(x)) dH 3
under the same constraint.

Application of Lagrange Multiplier Theorem

As our functionals are defined on the Hilbert space H 1 × H 1 , assumed to be differentiable, and as the constraint C(u, v) = 0 is linear, it is enough for applying Lagrange Multiplier Theorem [START_REF] Dell'isola | Lagrange Multipliers in Infinite-Dimensional Systems, Methods of[END_REF] to check that the range F of C is a Banach space and to identify its dual.

Actually F coincides with the Hilbert space H curl which has been extensively studied [START_REF] Dautray | Mathematical analysis and numerical methods for science and technology: volume 3 spectral theory and applications[END_REF][START_REF] Schweizer | On Friedrichs inequality, Helmholtz decomposition, vector potentials, and the div-curl lemma[END_REF]. This space is defined as the set of all functions in L 2 (Ω) whose curl (in the sense of distributions on Ω) belongs to L 2 (Ω). As curl(∇u) = 0, curl(C(u, v)) = curl(v) ∈ L 2 (Ω) and so the inclusion F ⊂ H curl is obvious. For proving the converse inclusion, let us consider w ∈ H curl and let us use Helmholtz decomposition of curl(w) ∈ L 2 (see theorem 4.2 of [START_REF] Schweizer | On Friedrichs inequality, Helmholtz decomposition, vector potentials, and the div-curl lemma[END_REF]) and write curl(w) = ∇φ + curl(v) with φ ∈ H 1 0 (Ω) and v ∈ H 1 (Ω). Obviously φ ∈ H 1 0 (Ω) satisfies ∆φ = 0 and thus, owing to the maximum principle, it vanishes. So curl(w -v) = 0 in the sense of distributions which is equivalent to the orthogonality of w -v with the curl of every smooth test function that is with every divergence-free smooth function. Hence (see [START_REF] Dautray | Mathematical analysis and numerical methods for science and technology: volume 3 spectral theory and applications[END_REF]

page 217) (w -v) ∈ L 2 (Ω) is a gradient: w-v = ∇u with u ∈ H 1 (Ω). We have w = v+∇u = C(-u, v) ∈ F .
The simplest way for identifying the dual F * consists in taking advantage of the Hilbert structure of H curl and applying Riez Theorem: to any element Λ of F * is associated a function f Λ in H curl such that, for any w ∈ H curl , Λ, w = Ω (f Λ • w + curl(f Λ ) • curl(w)) dH 3 . As smooth functions on Ω are dense (Theorem 2 p. 204 of [START_REF] Dautray | Mathematical analysis and numerical methods for science and technology: volume 3 spectral theory and applications[END_REF]) and continuously embedded in H curl , Λ is completely determined by its values for such smooth functions w. For such functions, one can integrate by parts in the sense of distributions on Ω and get Λ, w = f Λ +curl(curl(f Λ )), w . Any element of F * appears to be a distribution of order one on Ω.

So the theorem of Lagrange multipliers implies the existence, at any critical point (u, v), of a distribution5 Λ on Ω such that, for any pair of smooth test functions (h, k),

Ω (∂ 1 l • h(x) + ∂ 2 l • ∇h(x) + ∂ 3 l • ∇k(x)) dH 3 + Λ, (k -∇h) = 0 (5)
where the partial derivatives of l are taken at (u(x), ∇u(x), ∇v(x)). Using H 3 |Ω := 1 Ω H 3 , the restriction to Ω of the three-dimensional Hausdorff measure, this condition can be rewritten in terms of distributions on Ω

∂ 1 l H 3 |Ω , h + ∂ 2 l H 3 |Ω , ∇h + ∂ 3 l H 3 |Ω , ∇k + Λ, (k -∇h) = 0 or (∂ 1 l H 3 |Ω -div(∂ 2 l H 3 |Ω ) -div(Λ)), h + ( div(∂ 3 l H 3 |Ω ) + Λ), k = 0.
Hence, still in the sense of distributions on Ω,

∂ 1 l H 3 |Ω -div(∂ 2 l H 3 |Ω ) -div(Λ) = 0, and div(∂ 3 l H 3 |Ω ) + Λ = 0. Eliminating the Lagrange multiplier Λ is straightforward. We get ∂ 1 l H 3 |Ω -div(∂ 2 l H 3 |Ω ) + div(div(∂ 3 l H 3 |Ω
)) = 0. Now we make use twice of the divergence theorem [START_REF] Lagrange | [END_REF]. We obtain successively

(∂ 1 l -div(∂ 2 l)) H 3 |Ω + n • ∂ 2 l H 2 |∂Ω + div(div(∂ 3 l) H 3 |Ω ) -div(n • ∂ 3 l H 2 |∂Ω ) = 0, (∂ 1 l-div(∂ 2 l)+div(div(∂ 3 l))) H 3 |Ω +n•(∂ 2 l-div(∂ 3 l)) H 2 |∂Ω -div(n•∂ 3 l H 2 |∂Ω ) = 0.
We apply now a more general version of Stokes theorem, valid on a sub-manifold with boundary6 :

div(ϕ H 2 |∂Ω ) = div // (ϕ // ) H 2 |∂Ω -τ • ϕ H 1 |∂∂Ω -(n • ϕ) D n|∂Ω div(ϕ H 2 |∂Ω ), ψ = -ϕ H 2 |∂Ω , ∇ψ = - ∂Ω ϕ • ∇ψ dH 2 = - ∂Ω (n • ϕ)(n • ∇ψ) dH 2 - ∂Ω ϕ // • ∇ // ψ dH 2 = - ∂Ω (n • ϕ)(n • ∇ψ) dH 2 + ∂Ω div // (ϕ // )ψ dH 2 - ∂∂Ω (τ • ϕ // )ψ dH 1 .
where τ stands for the vector tangent to the sub-manifold and normal to its boundary, ϕ // for the projection of ϕ onto the tangent space, div // for the tangential divergence operator and D n|∂Ω for the normal derivative distribution. We obtain, by using this theorem on all faces (i.e. C 1 parts) of the boundary 7 ,

(∂ 1 l -div(∂ 2 l) + div(div(∂ 3 l))) H 3 |Ω -(n • (n • ∂ 3 l)) D n|∂Ω + n • (∂ 2 l -div(∂ 3 l) -div // ((n • ∂ 3 l) // )) H 2 |∂Ω -[[τ • (n • ∂ 3 l)]] H 1 |∂∂Ω = 0.
As the distributions [START_REF] Schwartz | [END_REF], we deduce the four criticality conditions,

H 3 |Ω , H 2 |∂Ω , H 1 |∂∂Ω , D n|∂Ω are orthogonal
(∂ 1 l -div(∂ 2 l) + div(div(∂ 3 l))) = 0 H 3 a.e. in Ω, (6)
n • (∂ 2 l -div(∂ 3 l) -div // ((n • ∂ 3 l) // )) = 0 H 2 a.e. on ∂Ω, (7)
[[τ • (n • ∂ 3 l)]] = 0 H 1 a.e. on ∂∂Ω, (8)
n • (n • ∂ 3 l) = 0 H 2 a.e. on ∂Ω. (9)
Of course these conditions can be made non homogeneous by adding some boundary terms in the Lagrangian. If we remind that constraint C(u, v) = 0 enables us to eliminate the variable v, these conditions correspond to the classical conditions (respectively: body force, surface boundary force, edge force, surface boundary double force) obtained by studying directly the second gradient Lagrangian [START_REF] Mindlin | Micro-structure in linear elasticity[END_REF][START_REF] Germain | La méthode des puissances virtuelles en mécanique des milieux continus[END_REF].

Note that the expression of the Lagrange multiplier is explicit. We have Λ = -div(∂ 3 l H 3 |Ω ) which becomes, by using (3),

Λ = -div(∂ 3 l) H 3 |Ω + (n • ∂ 3 l) H 2 |∂Ω .
It is now clear that Λ involves both a volume part and a boundary part. Assuming that Λ is a function generally leads to wrong results.

Note also that a seemingly natural idea could be to write the constraint under the form Ω v -∇u 2 dH 3 = 0 and to infer from the fact that this constraint takes values in R the existence of a scalar Lagrange multiplier. A straightforward computation would again lead to equation ( 5) but with Λ = 2λ(v -∇u). This is in complete contradiction with the previous remark. The reader must be aware that such a way of treating the constraint is erroneous as the function (u, v) → Ω v -∇u 2 dH 3 = 0 is not a submersion (see entry [START_REF] Dell'isola | Lagrange Multipliers in Infinite-Dimensional Systems, Methods of[END_REF]). 7 Each edge in ∂∂Ω is bordered by two faces. Thus the application of Stokes theorem on these two faces leads to two contributions on the edge involving the normals n 1 and n 2 of these faces and the vectors τ 1 and τ 2 tangent to the faces and normal to the edge. Notation

[[τ • (n • ∂ 3 l)]] stands for τ 1 • (n 1 • ∂ 3 l) + τ 2 • (n 2 • ∂ 3 l).

Approach through penalization

An intuitive way for dealing with constrained problems is to penalize the fact that the constraint is not verified. Looking for minimizers of the total energy, one decides to look for a sequence of approximate minimizers by considering, instead of the original potential subjected to the constraint, the original Lagrangian to which one adds a penalization term: for instance

L ε (u, v) = L(u, v) + ε -1 Ω C(u, v) 2 dH 3 . (10)
Denoting

µ ε := 2ε -1 C(u, v), any critical pair (u, v) satisfies Ω (∂ 1 l • h + ∂ 2 l • k + ∂ 3 l • ∇h + ∂ 4 l • ∇k + µ ε • (k -∇h)) dH 3 = 0 ( 11
)
which is exactly the same equation as equation ( 5) where the measure Λ ε := µ ε H 3 |Ω replaces the distribution Λ. The criticality conditions still read

∂ 1 H 3 |Ω -div(∂ 3 H 3 |Ω ) + div(µ ε H 3 |Ω ) = 0, ∂ 2 H 3 |Ω -div(∂ 4 H 3 |Ω ) + µ ε H 3 |Ω = 0, or equivalently ∂ 1 H 3 |Ω -div((∂ 2 + ∂ 3 ) H 3 |Ω ) + div(div(∂ 4 H 3 |Ω )) = 0, (12)
∂ 2 H 3 |Ω -div(∂ 4 H 3 |Ω ) + µ ε H 3 |Ω = 0, (13)
Using (3), equation ( 13) can be rewritten

(∂ 2 + div(∂ 4 ) + µ ε ) H 3 |Ω -(n • ∂ 4 ) H 2 |∂Ω = 0 (
Even if equation ( 12) has been already found and treated in the previous section, let us do it again taking into account (16). Using (3), taking into account (16) and using again (3), we get

0 = ∂ 1 H 3 |Ω -div((∂ 2 + ∂ 3 ) H 3 |Ω ) + div(div(∂ 4 ) H 3 |Ω -n • ∂ 4 H 2 |∂Ω ) = = ∂ 1 H 3 |Ω -div((∂ 2 + ∂ 3 -div(∂ 4 )) H 3 |Ω ) = = (∂ 1 -div(∂ 2 + ∂ 3 -div(∂ 4 ))) H 3 |Ω -n • (∂ 2 + ∂ 3 -div(∂ 4 )) H 2 |∂Ω .
which leads to

Recalling the relationship linking l and , bulk equilibrium equation ( 17) is identical to (6) while boundary force balance (18) is identical to [START_REF] Schwartz | [END_REF]. However equation ( 16) is much stronger than [START_REF] Dell'isola | How contact interactions may depend on the shape of Cauchy cuts in Nth gradient continua: approach "à la D'Alembert[END_REF] and there is no more edge equation corresponding to [START_REF] Forest | First vs. second gradient of strain theory for capillarity effects in an elastic fluid at small length scales[END_REF]. This discrepancy makes dangerous to interpret second gradient boundary conditions starting from a micromorphic model. The point is that one does not know how to pass to the limit as ε tends to zero. The considered penalization is a singular perturbation: we know from our previous study that equation ( 16), valid for any ε > 0, is no more valid in the limit. This comes from the fact that the function Λ ε converges to a measure partially concentrated on the boundary. Hence the way we decomposed distribution ( 14) is not valid in the limit.

  14)which leads to∂ 2 + div(∂ 4 ) + µ ε = 0 H 3 a.e. in Ω,(15)n • ∂ 4 = 0 H 2 a.e. on ∂Ω.

∂ 1 -

 1 div(∂ 2 + ∂ 3 -div(∂ 4 )) = 0 H 3 a.e. in Ω, (17) n • (∂ 2 + ∂ 3 -div(∂ 4 )) = 0 H 2 a.e. on ∂Ω.

Here H 3 stands for the Lebesgue measure on R 3 .

Rubber is an example of an incompressible elastic material. Incompressibility is also of high importance in hydrodynamics: in that case u must be interpreted as the velocity field and L(u) as the associated dissipation.

Note that ∂ 1 l = ∂ 1 , ∂ 2 l = ∂ 2 + ∂ 3 and ∂ 3 l = ∂ 4where ∂ i f stands for the partial derivative of a function f with respect to its i-th variable. The densities l and are assumed to be of class C 1 .

A frequent error consists in introducing the Lagrange multiplier Λ as a function. We will soon see that a part of Λ is a measure concentrated on the boundary of the domain. Note that the constraint v = ∇u has some implications on the boundary: owing to this constraint ∇u has, like v, a trace on the boundary of Ω which is not natural in H 1 (Ω).

Indeed, for any C 1 function ϕ and any smooth test function ψ,