Alexandra Carpen-Amarie

Djawida Dib
email: djawida.dib@inria.fr

Anne-Cécile Orgerie

Guillaume Pierre
email: guillaume.pierre@irisa.fr

Towards Energy-Aware IaaS-PaaS Co-design

Keywords: Cloud computing, energy-awareness, IaaS, PaaS

The wide adoption of the cloud computing paradigm plays a crucial role in the ever-increasing demand for energy-efficient data centers. Driven by this requirement, cloud providers resort to a variety of techniques to improve energy usage at each level of the cloud computing stack. However, prior studies mostly consider resource-level energy optimizations in IaaS clouds, overlooking the workload-related information locked at higher levels, such as PaaS clouds. In this position paper, we argue that cross-layer cooperation in clouds is a key to achieving an optimized resource management, both performance and energy-wise. To this end, we claim there is a need for a cooperation API between IaaS and PaaS clouds, enabling each layer to share specific information and to trigger correlated decisions. We identify the drawbacks raised by such co-design objectives and discuss opportunities for energy usage optimizations. Moreover, we outline the design of a set of extension modules for Libcloud to serve as building blocks for cross-layer information sharing and cooperation.

INTRODUCTION

Energy efficiency is a major concern in large-scale cloud infrastructures. In 2010, data centers worldwide reportedly consumed 203.4 to 271.8 TWh, which accounted for 1.1% to 1.5% of the global electricity use. This number is expected to grow even further [START_REF] Koomey | Growth in data center electricity use 2005 to 2010[END_REF]. To help make the planet greener while reducing their energy-related costs, cloud operators use a variety of techniques such as workload consolidation and dynamic virtual machine (VM) resizing. These strategies are however fundamentally limited by the lack of knowledge by the data center infrastructure about the applications that generate the workload. If the data center infrastructure understood better the full software stack it supports, it would be able to use more aggressive consolidation policies and derive further energy savings.

Platform-as-a-Service (PaaS) environments deliver users complete runtime environments for building, deploying and hosting software applications, while dispensing them from the complexity of managing the underlying software and hardware required for their applications. PaaS systems rely on computational resources obtained on-demand from Infrastructure-as-a-Service (IaaS) clouds, by acting as clients that can create/remove, suspend/resume, shutdown/restart, or migrate VMs.

This position paper claims that significant energy gains could be obtained by creating a cooperation API between the IaaS layer (in charge of handling elementary computing resources) and the PaaS layer (in charge of hosting applications). We discuss two complementary approaches for establishing such cooperation:

Cross-layer information sharing. Each of the two layers contains information which is potentially relevant for the other in order to reduce the energy footprint of a given application. For example, an energy-aware IaaS may have an estimate of the energy footprint of individual virtual machine types, which may be useful for the PaaS to preferentially select energy-efficient VM instance types. Conversely, a PaaS layer may build short-term traffic predictions which constitute useful hints for IaaS-level consolidation algorithms.

Cross-layer coordination. Both layers should coordinate their reconfiguration actions to reach a state where they help each other in achieving their goals rather than potentially take mutually detrimental decisions. For example, if the IaaS layer decides to migrate a specific VM, then the PaaS could in many cases temporarily reduce the load of the concerned VM to facilitate its migration. Conversely, if the PaaS layer gives early warnings to the IaaS about future requests for creating/destroying resources, the IaaS layer can prepare in advance for these changes. This paper proposes a research agenda towards the co-design of IaaS and PaaS layers. Section 2 first discusses the state of the art and motivation for this work. Section 3 details the issues created by the lack of IaaS-PaaS co-design. Then, Section 4 presents our early ideas about opportunities for improvement, while Section 5 outlines the design of interfaces to enable cross-layer cooperation in clouds. Finally, Section 6 concludes this work.

STATE OF THE ART

Many recent research works have been proposed to take into account the consumed energy in data centers. These works can be classified into two categories: (1) energy-aware systems and (2) energy efficient systems. The first category aims at informing cloud providers and users about the consumed energy for running their applications. For instance in [START_REF] Singh | A cloud-based consumer-centric architecture for energy data analytics[END_REF], the authors propose a cloudbased customer-centric architecture that allows consumers to own and control access to their energy usage data and have it analyzed using algorithms of their choice. The second category aims at optimizing the consumed energy for running applications. Different techniques are used to achieve energy efficiency. In [START_REF] Phan | Evolutionary multiobjective optimization for green clouds[END_REF], the authors propose to increase the use of renewable energy by dynamically moving services across data centers. The works in [START_REF] Jaiantilal | Modeling cpu energy consumption for energy efficient scheduling[END_REF][START_REF] Takouna | Dynamic configuration of virtual machine for power-proportional resource provisioning[END_REF] propose to adjust CPU voltage and frequency according to the load for saving energy. In [START_REF] Gadafi | Autonomic energy management in clusters[END_REF], the authors dynamically resize the active server set according to varying workload conditions. Workload consolidation on a limited number of servers is addressed in [START_REF] Duy | Performance evaluation of a Green scheduling algorithm for energy savings in cloud computing[END_REF][START_REF] Meisner | Powernap: eliminating server idle power[END_REF], in order to allow idle servers to be switched off and save energy.

Other works have focused on analyzing the impact of various optimization techniques on energy consumption. For instance, in [START_REF] Orgerie | Demystifying energy consumption in grids and clouds[END_REF], the authors investigate the tradeoffs and limitations of existing energy models for large-scale systems.

ISSUES

Energy-awareness in IaaS and PaaS

PaaS frameworks provide application developers with runtime environments where applications can be easily deployed and managed in the cloud [START_REF] Pierre | ConPaaS: a platform for hosting elastic cloud applications[END_REF][START_REF] Dib | Meryn: open, SLAdriven, cloud bursting PaaS[END_REF]. They do not provide cloud resources directly, but rather make use of an underlying IaaS layer in charge of resource management. Consequently, PaaS does not have access to information about the underlying resources, while IaaS does not have access to information about the running applications. However, both IaaS and PaaS systems target multi-objective optimizations including primarily cost, performance and energy consumption. These parameters are intrinsically related, thus requiring complex trade-offs to be made between them.

Traditionally, energy-awareness and energyefficiency have been addressed at the IaaS level with techniques presented in Section 2. To our knowledge, no research work has been targeting energy-efficiency or energy-awareness at the PaaS level. However, studies on real cloud infrastructures show that servers are considerably underutilized [START_REF] Zhang | Characterizing Task Usage Shapes in Google's Compute Clusters[END_REF], thus wasting large amounts of energy (Acc, 2010). This information about resource utilization is not available at the PaaS level, and nor is the energy consumption of resources. So, the PaaS layer has no means to provide energy-related metrics to users (e.g. energy-efficiency of their VMs, energy budget of their applications), and thus no way to increase their energy-awareness.

Conversely, the IaaS layer has no information about the applications and the users. As previously discussed, such knowledge would be useful to predict future resource utilization, thus leading to better energy-related decisions at the IaaS level.

The case for IaaS-PaaS co-design

gets released by the PaaS layer. The decision to shutdown this VM may have been taken several minutes in advance by the PaaS layer. If this information is not communicated to the IaaS layer, we take the risk that IaaS will invest previous resources (for example by migrating the VM) without seeing any benefit from this action (because the VM gets shut down just after).

Conversely, the PaaS layer may help the IaaS layer in performing its VM management actions. For example, it is often easy at the PaaS level to temporarily redirect one VM's workload to another (by redefining load balancing parameters for example). Offloading a VM for just a few tens of seconds may greatly facilitate IaaS-level management tasks such as VM migration.

In order to avoid counterproductive independent optimizations, we believe that the IaaS and the PaaS should share their energy-related information and coordinate their reconfiguration actions. This coordination aims at allowing system-level optimizations and trade-offs.

OPPORTUNITIES

In order to design efficient cloud frameworks, both performance and energy-wise, it is essential to understand how the cloud stack layers can interact and cooperate. We argue that there is a need to extend the traditional boundaries of each class of cloud computing frameworks to ensure an effective management of resources, while delivering seamless application performance and SLA compliance to higher-level users. In this section we identify a set of design opportunities for both PaaS and IaaS providers to help address these goals. The research directions we investigate are twofold: first, we focus on the bidirectional information sharing opportunities between these cloud layers and the coordinated decisions they can trigger. Second, we discuss security concerns that arise from such interactions.

Workload-aware IaaS frameworks

The emergence of a wide array of open-source IaaS solutions led to a gradual adoption of open standards and interfaces. This was a key incentive for enterprises to deploy private clouds in their data centers, allowing them to take advantage of the flexibility and performance of virtualized infrastructures, without the vendor lock-in shortcomings of proprietary solutions.

Most open IaaS cloud frameworks deliver similar services, typically focusing on resource manage-ment, VM lifecycle and on providing framework usage statistics. However, little has been done to achieve energy-efficient resource usage in such infrastructure clouds, despite the fact that lowering energy consumption is becoming a growing concern for data center management. This context provides a valuable opportunity for IaaS cloud providers to consider techniques that take into account the requirements of higher-level services, such as PaaS offerings, which traditionally build on infrastructure clouds. To this end, we believe that PaaS clouds can deliver a set of insights to the infrastructure level to guide its resource allocation and energy optimization decisions.

Passing workload information

VM allocation strategies typically rely on information concerning physical machine capabilities and their usage across the data center. Nevertheless, several factors may impact the performance of VMs, with critical side effects on application performance, user costs and power consumption. For instance, intensive network traffic between a few VMs may lower the available bandwidth across an entire cluster. On the other hand, PaaS virtual resource management services use the IaaS APIs to specify the number and type of VMs required upon application deployment. No information about the application class or its execution parameters is forwarded to the underlying infrastructure.

In this context, the answer to achieving better resource management at the infrastructure level is to enable IaaS schedulers to analyze and exploit a wider spectrum of parameters associated with the application behavior. Thus, as such workload properties are generally opaque to the underlying VM manager, a possible approach is to enable PaaS services to expose them to the infrastructure layer, as detailed below.

Resource usage patterns. The resource needs of VMs have a potentially heavy impact on the VM migration duration. As an example, it was shown that memory-bound VMs lead to very inefficient migration times [START_REF] Liu | Performance and energy modeling for live migration of virtual machines[END_REF]. IaaS systems may avoid such operations by flagging affected VMs at deployment time. A more fine-grained approach may enable the PaaS services to trigger such flags dynamically, as the application enters a memory-intensive stage, and allow migration outside these intervals.

Execution time. Cloud services can accommodate a large set of application types. For instance, web servers represent long-running jobs with specific access patterns. Oppositely, a short-lived VM, such as a MapReduce job, may raise specific scheduling constraints. Thus, as many MapReduce jobs run for short periods of time, migrating the VMs during ap-plication execution may degrade their performance to unacceptable levels (e.g., migrating the VM may take longer than its actual runtime).

Elasticity needs. Unlike the infrastructure manager, the PaaS services may be able to estimate the application requirements in terms of workload peaks. As elasticity is one of the most appealing cloud features, anticipating application needs and passing them on to the IaaS layer may result in more timely workload adaptation mechanisms, which cannot be achieved by the PaaS services alone.

Using workload information to reduce energy consumption

Such workload knowledge may allow IaaS providers not only to optimize resource allocation decisions, but also to improve the overall energy-consumption levels. This goal is typically achieved by adapting the number of used physical nodes to the workload and shutting down underutilized nodes [START_REF] Orgerie | Demystifying energy consumption in grids and clouds[END_REF]. Taking into account additional PaaS-level aspects may lead to significant and accurate energysaving decisions, as exemplified below.

Resource scheduling. Resource allocation in IaaS clouds is generally based on the VM type and straightforward scheduling policies. Beyond the VM type in terms of required CPUs and memory, an essential parameter is represented by the association of VMs with a specific virtual cluster, and thus, a single application. By taking into account this type of hint made available by PaaS-level services, an IaaS provider may change the allocation process of similar VMs. For example, VMs belonging to the same job may be deployed on the same physical machines to benefit from data locality and avoid unnecessary and performance-degrading network traffic.

Node shutdown policies. A widely studied mechanism to minimize energy consumption is shutting down unused nodes. However, resource managers usually lack any knowledge related to the type of workload or about the real-time evolution of resource requirements. Such information provided by higherlevel services such as PaaS providers may allow for more aggressive shutdown policies, without incurring performance degradations. When the infrastructure manager is aware of the short-term evolution of its workload, it can adopt less conservative policies to reduce the powered-on resources and the total energy usage.

Management of workload peaks. Cloud services have developed on the premise that they can provide on-demand resource adaptation to workload bursts, relieving the client from the burden of overprovisioning resources to cope with unexpected work-load fluctuations. Without any means to foresee the resource requirements of its scheduled jobs, an IaaS cloud is constrained to maintain a pool of unused running physical machines to accommodate possible workload peaks. To minimize the impact of such an approach on energy consumption, PaaS services may anticipate workload trends and, furthermore, forward them to the infrastructure level. Thus, IaaS schedulers can redirect small bursts to public clouds when the cost of switching nodes on is more important than the cost of obtaining remote VMs, both performance and energy-wise. When facing longer-duration workload peaks, the IaaS cloud can automatically adapt the pool of available physical machines, by enabling nodes in advance.

Optimizing PaaS resource management

Whereas advice from application-aware PaaS services may improve the efficiency of resource management decisions, PaaS frameworks can also benefit from infrastructure-level insights.

Exposing low-level knowledge to PaaS clouds

Infrastructure details. PaaS management services interact with the underlying infrastructure through the interfaces exposed by the IaaS provider, which typically hide low-level details from the user. As this information is provided in private clouds for administrative purposes, it can be also made available for higherlevel entities, such as PaaS frameworks running on the infrastructure. For instance, knowledge about cluster utilization and node properties may allow the PaaS services to select the most suitable VM types to execute an application, enabling the use of large VMs when the infrastructure comprises high-performance physical machines. Energy awareness. The first step towards energyefficient services at the platform level is user awareness, in particular if lowering the energy consumption is associated with monetary incentives. The user may be presented with several execution plans and the corresponding energy consumption estimations for the same application. He can thus choose the best execution schedule in terms of energy, regardless of a performance loss or a longer execution time. Nevertheless, the PaaS layer has no means to accurately estimate the energy consumption of a given application. It is at the IaaS level where such predictions can be made, based on the load of the infrastructure and the application requirements.

IaaS-PaaS coordinated decisions

Since PaaS services do not have full control over the management of their virtual resources, the key to an optimized resource usage and energy consumption lies in cross-layer coordination. Despite the fact that all VM allocation and placement decisions are taken at the infrastructure level, such a collaboration can lead to both improved scheduling mechanisms and a reduced incidence of conflicting actions.

Adaptive resource management. Equipped with detailed information about the user job, the IaaS layer can yield multiple scheduling options to handle it. According to user set thresholds, such as an execution deadline or a given energy budget, the IaaS resource manager may be extended to delay jobs according to the availability of energy-efficient or renewable energy sources (e.g., wind energy). PaaS may facilitate certain VM management actions as well. To reduce the time spent during VM migration, the PaaS can be instructed to temporarily decrease the load of the VM, which may entail performance and energy gains.

Preventing conflicting optimizations. Cloud cross-layer communication contributes to another unexpected, yet significant optimization angle for VM management, that is to prevent conflicting optimizations. As an example, the PaaS layer may detect an idle VM and redirect a slice of the job workload to it. Simultaneously, the infrastructure manager can decide to migrate the same idle VM, resulting in both an inefficient migration and a reduced performance for the PaaS application. To prevent such conflicting or counter-productive optimizations, it is essential to define a standardized interface and a set of parameters to be passed between the two layers, which can ensure portability and usability across cloud providers.

Privacy and security issues

As detailed before, cross-layer coordination can help increasing both energy awareness and energy efficiency of the entire cloud system. However, this coordination needs to be clearly defined as it can raise major concerns about the security and privacy of cloud providers and users. Indeed, it has been shown that having access to the energy consumption of a cloud server can allow people to guess with high probability what type of application (among various possible ones) is running in its virtual machines [START_REF] Hlavacs | Energy consumption side-channel attack at virtual machines in a cloud[END_REF]. So, passing the information about energy consumption of servers is not possible for privacy's sake. It is therefore necessary to consider a coordination of these platforms that guarantees the privacy of the users and of the applications, while providing useful information to both layers in order to save energy. To do so, it is required to identify sensitive information that cannot be communicated between the two layers: energy consumption of physical servers, location of virtual machines (allowing users to target specific physical machines in order to be co-located with specific virtual machines). Moreover, the IaaS provider and the PaaS provider can be different, and they may not want to share sensitive information about the utilization of their platform.

PROPOSAL

To facilitate the interaction between cloud layers while preserving the separation and the interoperability across the cloud stack, we argue there is a need for an abstraction layer proposing coordination APIs. For instance, let us take the example of the Apache Libcloud (lib, 2014) library. Libcloud consists of a set of Python libraries providing a unified abstraction across IaaS cloud APIs such as Amazon EC2, OpenStack or OpenNebula. It hides the differences between the various cloud providers by offering generic APIs for compute, storage, authentication, load balancers and DNS services. We believe such interfacing libraries can be the key of cross-layer interoperability in the cloud stack. The development of client libraries that expose APIs to promote collaborative decisions between IaaS clouds and higher-level frameworks may come as an incentive for cloud providers to exploit the potential benefits of inter-layer cooperation.

We aim at defining a set of APIs to complement existing Libcloud features with new mechanisms for designing IaaS-PaaS cooperation strategies. Our goal is twofold. First, we plan to enable IaaS clouds to expose extended infrastructure information to upper layers or applications. In this context, two sets of APIs are required to enable PaaS providers to optimize resource management, as detailed in Section 4.2. On the one hand, the IaaS provider can supply data about the infrastructure status, including cluster utilization and node properties, which can be accessed on-demand through the Libcloud API. On the other hand, the Libcloud suite may include a set of APIs dedicated to exposing infrastructure energy consumption at the higher levels, as the first step towards energy-aware PaaS and IaaS applications.

Our second goal is to add abstract triggers at the Libcloud level to enable PaaS managers to push workload-related information to the IaaS clouds. As the applications executed on VMs rented from IaaS clouds are opaque to the infrastructure level, such APIs are essential to make application information available for the IaaS resource allocation and management mechanisms. Consequently, such APIs should allow PaaS managers to ask the IaaS provider to flag the activity of specific VMs or groups of VMs as CPU-, IO-or memory-intensive, to predict workload peaks or to provide hints regarding the execution time of groups of VMs that belong to the same job.

Such an API-driven approach to cross-layer cooperation achieves several benefits. First, each cloud provider may choose to implement only the required plugins according to the amount of information it intends to share. While private clouds may implement the whole set of energy-related and cooperation APIs, public providers may decide to share only non-sensitive data. Furthermore, an abstraction layer such as Libcloud fosters interoperability and portability between clouds, such that PaaS providers can make data available and collect infrastructure information according to predefined formats, while each IaaS cloud can provide its own interface implementation according to its internal structure and APIs.

CONCLUSION

The growing demand for "green" and energy-efficient data centers is a key driver in the design of cloud computing frameworks, which target not only performance, but also cost efficiency. In particular, designing virtualized environments for energy efficiency has led to a wide spectrum of optimization techniques, most of which are intended for IaaS cloud offerings. Whereas higher-level PaaS services hold key information related to workload properties and execution, they do not have direct control over the resourcemanagement mechanisms that could trigger energy gains. In this work, we argued that there is a need for coordinated actions between IaaS and PaaS cloud layers. We believe that the conventional boundaries of the cloud stack layers can be extended to promote energy-awareness at higher levels, while enhancing low-level resource managers with detailed workload profiles to achieve fine-grained energy saving mechanisms.

We investigated a set of opportunities for both IaaS and PaaS providers to address these objectives and highlighted the tradeoffs of cross-layer information passing. As future work, our goal is to design the proposed API suite by extending Libcloud, in order to enable cross-layer cooperation between open-source IaaS and PaaS cloud implementations and to develop prototypes of cooperating cloud frameworks.

We believe that the separation of the cloud stack in two distinct IaaS and PaaS layers, while having great advantages for portability and separation of concerns, is detrimental in terms of energy awareness. If each layer is allowed to take energy-related decisions independently, these uncoordinated actions can lead to significant resource waste and performance degradation, possibly negating the benefits of energy awareness altogether.For instance, the IaaS layer can decide to migrate a virtual machine (VM) in order to perform a better consolidation for energy-efficiency purposes. Yet, this same VM may end a few seconds later because it


Eric Canc Ès
email: cances@cermics.enpc.fr

Rachida Chakir
email: rachida.chakir@ifsttar.fr

Yvon Maday
email: maday@ann.jussieu.fr

Two-grid methods for a class of nonlinear elliptic eigenvalue problems

Keywords: Nonlinear eigenvalue problem, Spectral and pseudo spectral approximation, Finite element approximation, Ground state computation, Numerical analysis, Two-grid method

In this paper, we introduce and analyze some two-grid methods for nonlinear elliptic eigenvalue problems of the form -div(A∇u)+Vu+ f (u 2 )u = λ u, u L 2 = 1. We provide a priori error estimates for the ground state energy, the eigenvalue λ , and the eigenfunction u, in various Sobolev norms. We focus in particular on the Fourier spectral approximation (for periodic boundary conditions), and on the P 1 and P 2 finite element discretizations (for homogeneous Dirichlet boundary conditions), taking numerical integration errors into account. Finally, we provide numerical examples illustrating our analysis.

Introduction

Nonlinear eigenvalue problems are encountered in various applications in sciences and engineering, including the simulation of Bose-Einstein condensates (Gross-Pitaevskii equation, see e.g. [START_REF] Pitaevskii | Bose-Einstein Condensation[END_REF]), electronic structure calculation (Hartree-Fock method, orbital free and Kohn-Sham Density Functional Theory), and the study of the vibration modes of structures in nonlinear elasticity.

The first results on the numerical analysis of nonlinear eigenvalue problems have been published in [START_REF] Zhou | An analysis of finite-dimensional approximations for the ground state solution of Bose-Einstein condensates[END_REF]. These first results were improved by three of us in Cancès et al. (2010), where optimal a priori error bounds for nonlinear elliptic eigenvalue problems were obtained for the first time. The techniques introduced in Cancès et al. (2010), based on estimates in negative Sobolev norms, have then been applied to a variety of nonlinear eigenvalue problems (see Cancès et al. (2012); [START_REF] Chen | Numerical analysis of finite dimensional approximations of Kohn-Sham models[END_REF]), among which the Kohn-Sham problem [START_REF] Kohn | Self-consistent equations including exchange and correlation effects[END_REF]), which is currently one of the most widely used models in computational physics and chemistry.

As in Cancès et al. (2010), we focus on the nonlinear elliptic eigenvalue problems arising in the study of variational problems of the form

I = inf E(v), v ∈ X, Ω v 2 = 1 , (1.1)
where Ω is a regular bounded domain or a rectangular brick of R d and X = H 1 0 (Ω ), or Ω is the unit cell of a periodic lattice R of R d and X = H 1 # (Ω ), with d = 1, 2, or 3, H 1 # (Ω ) denoting the space of the restrictions to Ω of the H 1 loc , R-periodic functions on R d , and where the energy functional E is of the form

E(v) = 1 2 a(v, v) + 1 2 Ω F(v 2 (x)) dx, with a(u, v) = Ω (A∇u) • ∇v + Ω Vuv.
In all what follows, we assume that

• A ∈ (L ∞ (Ω )) d×d ; A(x) is symmetric for almost all x ∈ Ω ;

∃ α > 0 such that ξ T A(x)ξ α|ξ | 2 , ∀ ξ ∈ R d and almost all x ∈ Ω ;

(1.2) • V ∈ L p (Ω ) for some p > max(1, d/2);

(1.3) • F ∈ C 1 ([0, +∞), R) ∩C 2 ((0, ∞), R), F (0) = 0 and F > 0 on (0, +∞);

(1.4) ∃ 0 q < 2, ∃C ∈ R + such that ∀t 0, |F (t)| C(1 + t q );

(1.5) • F (t)t is locally bounded on [0, +∞).

(1.6)

To simplify the notation, we denote by f the derivative of F. Note that there is no loss of generality in assuming in (1.4) that f (0) = F (0) = 0 since the minimizers of (1.1) are not modified if F(t) is replaced with F(t) + ct.

Problem (1.1) has exactly two minimizers u and -u, one of them, say u, being positive on Ω . In all what follows, u will be the positive minimizer of (1.1). The function u is solution to the Euler-Lagrange equation ∀v ∈ X, E (u)λ u, v X ,X = 0, (1.7)

for some λ ∈ R (the Lagrange multiplier associated with the constraint u L 2 = 1) and equation (1.7), complemented with the constraint u L 2 = 1, takes the form of the nonlinear eigenvalue problem

A u u = λ u, u L 2 = 1, (1.8)
where for all v ∈ X,

A v = -div(A∇•) +V + f (v 2 )
is a linear self-adjoint operator on L 2 (Ω ) with form domain X. Note that E (u) = A u u. It can then be inferred from (1.8) that u ∈ X ∩ C 0 (Ω ), u > 0 in Ω , and λ is the ground state eigenvalue ofA u . An important point is that λ is a simple eigenvalue of A u . These results are classical; their proofs are recalled in Cancès et al. (2010).

We now consider a family of finite-dimensional subspaces (X δ ) δ >0 of X such that ∀v ∈ X, lim

δ →0 min v δ ∈X δ v -v δ H 1 = 0 (1.9)
and the variational approximations of (1.1) consisting in solving

I δ = inf E(v δ ), v δ ∈ X δ , Ω v 2 δ = 1 .
(1.10)

Problem (1.10) has at least one minimizer u δ such that (u, u δ ) L 2 0, which satisfies

∀v δ ∈ X δ , A u δ u δ , v δ X ,X = λ δ (u δ , v δ ) L 2 ,
for some λ δ ∈ R. It is easily seen that (see, e.g., Cancès et al. (2010); [START_REF] Zhou | An analysis of finite-dimensional approximations for the ground state solution of Bose-Einstein condensates[END_REF]) (1.11) or, in words, that the approximate ground state eigenfunction converges to the exact ground state eigenfunction in the H 1 -norm, from which we deduce that I δ and λ δ converge to I and λ , respectively, when δ goes to 0. Optimal convergence rates have been obtained in Cancès et al. (2010) (under stronger assumptions on the nonlinearity F) for finite element and spectral Fourier discretizations.

lim δ →0 u δ -u H 1 = 0,
The numerical simulation of problem (1.10) can be too costly if the approximation space X δ is highdimensional. We will denote by X δ f such a space and call it the fine discretization space. In two-grid methods, problem (1.10) is first solved in a lower-dimensional approximation space X δ c ⊂ X δ f , that we will call the coarse discretization space. Then, the so-obtained solution u δ c is improved by solving a linearized problem in the fine discretization space X δ f . A nice feature of this approach is that, for appropriate choices of the linearized problem and of the coarse discretization space X δ c , the solution u δ c δ f obtained with the two-grid method has the same accurary as the solution u δ f obtained by solving the nonlinear problem (1.10) in the fine discretization space X δ f . Two-grid methods thus allow us to obtain the same accuracy at a much lower price. Such methods were first introduce ind [START_REF] Xu | Local and parallel finite element algorithms based on two-grid discretizations[END_REF] in the framework of nonlinear elliptic boundary value problems. A different two-grid approach is presented in [START_REF] Henning | Two-Level discretization techniques for ground state computations of Bose-Einstein condensates[END_REF]. The evaluation of its computational cost for a given accuracy with respect to our approach needs to be further analyzed. This article is organized as follows. In Section 2, we introduce three different two-grid algorithms to solve (1.1). In Section 3, we provide some abstract a priori error analysis for one of these algorithms. We then show how these abstract results can be applied to spectral Fourier and finite element discretizations in Sections 4 and 5 respectively. Numerical integration errors are dealt with in Section 6. Finally, we give several numerical examples to illustrate our theoretical results in Section 7. The class of nonlinear eigenvalue problems considered in this work is very similar to the one considered in the previous work Cancès et al. (2010) (only some assumptions on the nonlinearity F will differ). For this reason, some of the proofs of the results below are simple adaptations of proofs in Cancès et al. (2010), and will therefore not be detailed for the sake of brevity. Let us mention that some of the results contained in this article have been published (in French) in the PhD thesis of the second author [START_REF] Chakir | Contribution à l'analyse numérique de quelques problèmes en chimie quantique et mécanique[END_REF]).

Two-grid algorithms

Let X δ c and X δ f be coarse and fine discretization spaces such that X δ c ⊂ X δ f ⊂ X. As mentioned above, two-grid methods consist, first in computing a solution of (1.10) in a coarse discretization space X δ c and, second in improving it by solving a linearized problem in the fine discretization space X δ f .

Several two-grid algorithms can therefore be proposed, depending on the type of linear problem we choose to solve in the fine discretization space. In the following, we introduce three of them, the first and third steps of these three schemes being the same.

1. Solve (1.10) in the coarse discretization space X δ c . Recall that the solution u δ c of this problem is such that there exists

λ δ c ∈ R such that (λ δ c , u δ c ) is also solution to the nonlinear eigenvalue problem: find (λ δ c , u δ c ) ∈ R × X δ c such that ∀v δ c ∈ X δ c , A u δc u δ c , v δ c X ,X = λ δ c (u δ c , v δ c ) L 2 .
2. Two-grid scheme 1. Solve the following linear eigenvalue problem in the fine space X δ f :

             find (λ δ c δ f , u δ c δ f ) ∈ R × X δ f such that ∀v δ f ∈ X δ f , a(u δ c δ f , v δ f ) + Ω f (u 2 δ c ) u δ c δ f v δ f = λ δ c δ f Ω u δ c δ f v δ f , u δ c δ f L 2 = 1, (u, u δ c δ f ) L 2 0, λ δ c δ f
is the lowest eigenvalue of the above spectral problem.

(2.1) Two-grid scheme 2a. Solve the following linearized right-hand side problem in the fine space

X δ f :    find u δ c δ f ∈ X δ f such that ∀v δ f ∈ X δ f , a(u δ c δ f , v δ f ) + Ω f (u 2 δ c ) u δ c δ f v δ f = λ δ c Ω u δ c v δ f .
Two-grid scheme 2b. Solve the following linearized right-hand side problem in the fine space

X δ f :    find u δ c δ f ∈ X δ f such that ∀v δ f ∈ X δ f , a(u δ c δ f , v δ f ) = - Ω f (u 2 δ c ) u δ c v δ f + λ δ c Ω u δ c v δ f . 3. Compute the Rayleigh quotient λ δ c δ f for u δ c δ f : λ δ c δ f = A u δc δ f u δ c δ f , u δ c δ f X ,X u δ c δ f 2 L 2 .
(2.2)

In the limit δ f = 0 (that corresponds to X δ f = X), the second step of scheme 1 amounts to computing the ground state (λ δ c 0 , u δ c 0 ) of the self-adjoint operator A u δc , while schemes 2a and 2b amount to solving the boundary value problems

A u δc u δ c 0 = λ δ c u δ c and A 0 u δ c 0 = (λ δ c -f (u 2 δ c
))u δ c , respectively. In this paper, we shall focus on the analysis of the first scheme, both from the theoretical and simulation points of view. The analysis of the other two schemes will be the matter of a forthcoming work. Let us just mention here that on preliminary simulations, schemes 2a and 2b provide similar results as scheme 1 (see [START_REF] Chakir | Contribution à l'analyse numérique de quelques problèmes en chimie quantique et mécanique[END_REF]).

Abstract error analysis of scheme 1

We denote by u the unique positive solution of (1.1), by u δ c a minimizer of the discretized nonlinear problem (1.10) such that (u, u δ c ) L 2 0, and by u δ c δ f the approximation of u computed with scheme 1.

The aim of this section is to establish error bounds on u -

u δ c δ f H 1 , u -u δ c δ f L 2 and |λ -λ δ c δ f
|, in the general framework of assumptions (1.2)-(1.6) and (1.9).

Preliminaries

Our analysis relies on the introduction of the solution u δ c 0 of the two-grid scheme in the limiting case when δ f = 0 (that is for X δ f = X). Recall that u δ c 0 is the positive ground state eigenfunction of A u δc . We denote by λ δ c 0 the associated eigenvalue. The minmax principle gives

λ δ c 0 = inf A u δc v, v X ,X , v ∈ X, Ω v 2 = 1 , (3.1)
and the solution (λ

δ c δ f , u δ c δ f
) provided by the two-grid scheme 1 can then be interpreted as the solution of the variational approximation

λ δ c δ f = inf A u δc v δ f , v δ f X ,X , v δ f ∈ X δ f , Ω v 2 δ f = 1 (3.2)
of problem (3.1) in the discretization space X δ f . Problem (3.2) has at least one minimizer u δ c δ f , which satisfies (2.1), for some λ δ c δ f ∈ R. Note that, when

δ f = δ c , u δ c δ c = u δ c is solution to (3.2).
The following numerical analysis relies on the properties of the mapping v → (λ v , z v ), where λ v denotes the lowest eigenvalue of the self-adjoint operator A v and z v > 0 the associated positive normalized eigenfunction:

   A v z v = λ v z v , z v > 0, z v L 2 = 1.
The function z v is also the minimizer of the problem inf A v w, w X ,X , w ∈ X,

Ω w 2 = 1 ,
(which amounts to minimizing the Rayleigh quotient associated with the self-adjoint operator A v ). In the special case when v = u, we have (λ u , z u ) = (λ , u).

The following technical lemmas will be used throughout the article. For the reader's convenience, we first state all the lemmas, and postpone their proofs until the end of the section. LEMMA 3.1 Under assumptions (1.2), there exist

β 0 ∈ R + and M ∈ R + such that ∀w ∈ X, α 2 ∇w 2 L 2 -β 0 w 2 L 2 a(w, w) M w 2 H 1 . (3.3)
We recall that f denotes the derivative of F and u the unique positive minimizer of (1.1) that satisfies A u u = λ u with λ ∈ R.

LEMMA 3.2 Assume that F satisfies assumptions (1.4)-(1.6). Denoting by r = 6 5-2q ( 6 5 r < 6), there exists a constant C ∈ R + such that for all (v, w, z) ∈ X 3 ,

Ω f (v 2 )wz C 1 + v 2q L 6 w L 6 z L r , (3.4) Ω f (v 2 ) -f (u 2 ) uw C 1 + v 2q L 6 w L 6 u -v L r , (3.5) Ω f (v 2 ) -f (u 2 ) v 2 C 1 + v 2q+1 L 6 u -v L r , (3.6) 0 Ω F(v 2 ) -F(u 2 ) -f (u 2 ) v 2 -u 2 C 1 + v 2q+1 L 6 u -v L r . (3.7)
Besides, in the case where X = H 1 # (Ω ), there exists

C ∈ R + such that for all (v, w) ∈ X 2 , Ω f (u 2 ) -f (v 2 ) w 2 C Ω 1 u |v| (u -v)w 2 (X = H 1 # (Ω ) only), (3.8)
while, in the case where X = H 1 0 (Ω ), for all ε > 0, there exists

C ε ∈ R + such that for all (v, w) ∈ X 2 , Ω f (u 2 ) -f (v 2 ) w 2 ε w 2 L 2 +C ε Ω 1 u |v| (u -v)w 2 (X = H 1 0 (Ω ) only). (3.9) LEMMA 3.3 There exist 0 < M 2 M 1 < ∞ such that ∀v ∈ X, 0 (A u -λ )v, v X ,X M 1 v 2 H 1 , (3.10) and ∀v ∈ u ⊥ := {v ∈ X, (u, v) L 2 = 0} , M 2 v 2 H 1 (A u -λ )v, v X ,X . (3.11)
Moreover, there exists γ > 0 such that, for all w ∈ X such that w L 2 = 1 and (u, w) L 2 0,

γ w -u 2 H 1 (A u -λ )(w -u), (w -u) X ,X .
(3.12)

The properties of the ground state eigenpair (λ v , z v ) of A v are collected in the following lemma.

LEMMA 3.4 There exists a constant

C ∈ R + such that ∀v ∈ X, |λ v | + z v 2 H 1 C 1 + v 2q L 6 , (3.13) ∀v ∈ X, z v -u H 1 C 1 + v 2q L 6 u -v L max(r,2) . (3.14)
In addition,

• in the case when X = H 1 # (Ω ), there exists

C ∈ R + such that ∀v ∈ X, |λ v -λ | C 1 + v 2q L 6 v -u L max(r,2) (X = H 1 # (Ω ) only); (3.15)
• in the case when X = H 1 0 (Ω ), there exists, for any ε > 0, a constant

C ε ∈ R + such that ∀v ∈ X, |λ v -λ u | 2ε +C ε 1 + v 2q L 6 u -v L max(r,2) (X = H 1 0 (Ω ) only). (3.16)
Estimates (3.15)-(3.16) are sufficient for our purpose, but are not optimal; refined estimates are actually given in the proof of Lemma 3.4. For all v ∈ X, we denote by λ 2,v the second eigenvalue of A v . Since λ = λ u is a simple eigenvalue of A u , there is a gap, denoted by g = λ 2,uλ u > 0, between the first and second eigenvalues of A u .

LEMMA 3.5 There exists 0

< η 1 such that for all v ∈ X such that v -u H 1 η, we have λ 2,v -λ v g/2.
PROPOSITION 3.1 There exist η > 0 and 0

< c 0 C 0 < ∞ such that for all v ∈ X such that v-u H 1 η and all w ∈ X such that w L 2 = 1 and (z v , w) L 2 0, we have c 0 w -z v 2 H 1 (A v -λ v )(w -z v ), (w -z v ) X ,X C 0 w -z v 2 H 1 . (3.17)
Proof of Lemma 3.1. For brevity, we only explain in detail the arguments for d = 3, in which case p > 3/2. Under assumptions (1.2), there exists a positive constant M such that ∀w ∈ X, a(w, w)

A L ∞ ∇w 2 L 2 + V L p w 2 L 2p M w 2 H 1 ,
where 1 p = (1p -1 ) -1 < 3. Using Hölder's inequality, we have that for any w ∈ X, a(w, w)

= Ω A∇w • ∇w + Ω V w 2 α ∇w 2 L 2 -V L p w 2-3/p L 2 w 3/p L 6 α ∇w 2 L 2 -C 3/p 6 V L p w 2-3/p L 2 w 3/p H 1 , where C 6 is the Sobolev constant such that for all v ∈ X, v L 6 C 6 v H 1 . Using Young's inequality, we have for all ε > 0 and w ∈ X, C 3/p 6 V L p w 2-3/p L 2 w 3/p H 1 = 1 ε C 3/p 6 V L p w 2-3/p L 2 ε w 3/p H 1 2p -3 2p 1 ε C 3/p 6 V L p w 2-3/p L 2 2p/(2p-3) + 3 2p ε w 3/p H 1 2p/3 = 2p -3 2p 1 ε C 3/p 6 V L p 2p/3 3/(2p-3) w 2 L 2 + 3 2p ε 2p/3 w 2 H 1 .
Choosing 3 2p ε 2p/3 = α 2 , we get

C 3/p 6 V L p w 2-3/p L 2 w 3/p H 1 2p -3 2p 3 pα C 2 6 V 2p/3 L p 3/(2p-3) w 2 L 2 + α 2 w 2 L 2 + α 2 ∇w 2 L 2 , which implies a(w, w) α 2 ∇w 2 L 2 - 2p -3 2p 3 pα C 2 6 V 2p/3 L p 3/(2p-3) + α 2 w 2 L 2 .
Hence, there exists a positive constant

β 0 such that ∀w ∈ X, a(w, w) α 2 ∇w 2 L 2 -β 0 w 2 L 2 .
This completes the proof.

Proof of Lemma 3.2. In this proof, C denotes a non-negative constant independent on v, w and z, but whose value is allowed to change from one line to another. We recall that r = 6 5-2q . Proof of (3.4). It follows from assumption (1.5) that for all (v, w, z) ∈ X,

Ω f (v 2 )wz C Ω (1 + |v| 2q )|w| |z| C 1 + v 2q L 6 w L 6 z L r .
Proof of (3.5). We first write

Ω ( f (v 2 ) -f (u 2 ))uw = Ω w v,u (u -v)w, with w v,u = - f (u 2 ) -f (v 2 ) u -v u. (3.18) Since u ∈ L ∞ (Ω ) and u 0, it holds | w v,u | 2 f (v 2 ) + f ( u 2 L ∞ ) , when |v| < u/2, 4 sup t∈(0,4 u 2 L ∞ ] F (t)t, when u/2 |v| < 2u, | f (v 2 )| + f ( u 2 L ∞ ), when |v| 2u. (3.19)
The above estimate is easily obtained in the case when |v| < u/2 or |v| 2u. When u/2 |v| < 2u, we observe that

f (u 2 ) -f (v 2 ) = v 2 u 2 F (t) dt = v 2 u 2 tF (t) t dt 2   sup t∈(0,4 u 2 L ∞ ] F (t)t   |ln u -ln |v|| .
It follows that when u/2 |v| < 2u, there exists u/2 ξ < 2u such that

| w v,u | 2   sup t∈(0,4 u 2 L ∞ ] F (t)t   ln u -ln v u -v u = 2   sup t∈(0,4 u 2 L ∞ ] F (t)t   u ξ 4 sup t∈(0,4 u 2 L ∞ ] F (t)t.
Thus, (3.19) is proved. This estimate, together with assumptions (1.5) and (1.6), yields which, combined with (3.18), straightforwardly leads to (3.5).

| w v,u | C 1 + |v| 2q ,
Proof of (3.6). For all v ∈ X, we can write

Ω ( f (v 2 ) -f (u 2 ))v 2 = Ω w v,u (v -u), with w v,u = v 2 f (v 2 ) -f (u 2 ) v -u . (3.20)
As u ∈ L ∞ (Ω ), we have (see the proof of (Cancès et al., 2010, Theorem 1))

|w v,u | 12u sup t∈(0,4 u 2 L ∞ ] F (t)t, when |v| < 2u, 2 | f (v 2 )| + f ( u 2 L ∞ ) |v|, when |v| 2u.
We infer from assumptions (1.5) and (1.6) that

|w v,u | C(1 + |v| 2q+1 ). (3.21)
Putting together (3.20) and (3.21), we obtain

Ω ( f (v 2 ) -f (u 2 ))v 2 w v,u L 6/(1+2q) u -v L r C 1 + v 2q+1 L 6 u -v L r .
Proof of (3.7). The left-hand side inequality in (3.7) follows from the convexity of F (assumption (1.4)).

On the other hand,

∀v ∈ X, Ω F(v 2 ) -F(u 2 ) -f (u 2 ) v 2 -u 2 = Ω f (ξ ) -f (u 2 ) (v 2 -u 2 ) with ξ ∈ [min(u 2 , v 2 ), max(u 2 , v 2 )].
Using assumption (1.5) and the boundedness of u, we get

Ω F(v 2 ) -F(u 2 ) -f (u 2 ) v 2 -u 2 C Ω 1 + |v| 2q+1 |u -v| C 1 + v 2q+1 L 6 u -v L r .
Proof of (3.8). We assume here that X = H 1 # (Ω ). Since u is continuous, everywhere positive, and periodic, there exists a positive constant α 0 such that u α 0 > 0. Denoting by

Ω -= x ∈ Ω | |v(x)| < α 0 2 and Ω + = x ∈ Ω | |v(x)| α 0 2 , we have ∀(v, w) ∈ X 2 , Ω ( f (u 2 ) -f (v 2 ))w 2 = Ω - ( f (u 2 ) -f (v 2 ))w 2 + Ω + ( f (u 2 ) -f (v 2 ))w 2 . (3.22) Since, over Ω -, v is such that |v| < u u L ∞ and f is monotonically increasing, we have | f (v 2 )| | f ( u 2 L ∞ )| over Ω -, so that Ω - ( f (u 2 ) -f (v 2 ))w 2 2 f ( u 2 L ∞ ) Ω - u -v u -v w 2 4 f ( u 2 L ∞ ) α 0 Ω - (u -v)w 2 C Ω - (u -v)w 2 .
(3.23)

Denoting by

Ω 1 + = x ∈ Ω | u(x) > |v(x)| α 0 2 and Ω 2 + = {x ∈ Ω | |v(x)| u(x) α 0 }
, and using the fact that f is monotonically increasing, we obtain

Ω + ( f (u 2 ) -f (v 2 ))w 2 = Ω 1 + ( f (u 2 ) -f (v 2 ))w 2 + Ω 2 + ( f (u 2 ) -f (v 2 ))w 2 Ω 1 + ( f (u 2 ) -f (v 2 ))w 2 = Ω 1 + f (ξ )(u 2 -v 2 )w 2 , where ξ ∈ [v 2 , u 2 ] ⊂ α 2 0 /4, u 2 L ∞ . Since F ∈ C 2 ((0, ∞), R) (assumption (1.4)), we get Ω + ( f (u 2 ) -f (v 2 ))w 2 C Ω 1 + (u 2 -v 2 )w 2 C Ω 1 + 2u(u -v)w 2 C Ω 1 + (u -v)w 2 . (3.24)
Combining (3.22), (3.23), and (3.24), we obtain (3.8).

Proof of (3.9). We now consider the case when X = H 1 0 (Ω ). Since F ∈ C 1 ([0, +∞), R) and f (0) = 0, there exists, for any ε > 0, a constant β ε > 0 such that for all 0 t β 2 ε ,

| f (t)| ε. (3.25)
Since f is monotonically increasing, we have for all

(v, w) ∈ X 2 , Ω ( f (u 2 ) -f (v 2 ))w 2 = Ω 1 u |v| ( f (u 2 ) -f (v 2 ))w 2 + Ω 1 u>|v| ( f (u 2 ) -f (v 2 ))w 2 Ω 1 u>|v| ( f (u 2 ) -f (v 2 ))w 2 . (3.26) Denoting by Ω 1,ε = {x ∈ Ω | 0 |v(x)| < u(x) < β ε }, Ω 2,ε = {x ∈ Ω | 0 |v(x)| < β ε /2, u(x) β ε },
and

Ω 3,ε = {x ∈ Ω | β ε /2 < |v(x)| < u(x)},
we split the right-hand side of (3.26) into three parts. Using (3.25) and the boundedness of u, we get

Ω 1,ε ( f (u 2 ) -f (v 2 ))w 2 2ε w 2 L 2 ,
and

Ω 2,ε ( f (u 2 ) -f (v 2 ))w 2 4 β ε f ( u 2 L ∞ ) Ω 2,ε (u -v)w 2 .
We then note that there exists ξ with v 2 ξ u 2 such that

Ω 3,ε ( f (u 2 ) -f (v 2 ))w 2 = Ω 3,ε f (ξ )(u 2 -v 2 )w 2 2 max t∈[β 2 ε /4, u 2 L ∞ ] F (t) u L ∞ Ω 3,ε (u -v)w 2 .
Thus, (3.9) is proved with

C ε = 4 β ε f ( u 2 L ∞ ) + 2 max t∈[β 2 ε /4, u 2 L ∞ ] F (t) u L ∞ . Proof of Lemma 3.3.
> 0 such that ∀v ∈ X, (A u -λ )v, v X ,X η v 2 L 2 -|(u, v) L 2 | 2 0. (3.27) Since w L 2 = 1 and u L 2 = 1, we have w -u 2 L 2 -|(w -u, u) L 2 | 2 w -u 2 L 2 -(1 -(w, u) L 2 ) = 1 2 w -u 2 L 2 ,
which together with (3.27) implies

(A u -λ )(w -u), (w -u) X ,X η 2 w -u 2 L 2 . (3.28)
In view of inequality ( 22) in Cancès et al. (2010), there exists a constant C ∈ R + such that

(A u -λ )(w -u), (w -u) X ,X α 2 ∇(w -u) 2 L 2 -C w -u 2 L 2 . (3.29)
We obtain (3.12) with γ = αη 2(η+2C) by combining (3.28) and (3.29). Proof of Lemma 3.4.

In this proof, C, C 1 , C 2 denote non-negative constants independent of v, but whose values are allowed to change from one line to another. As λ v is the lowest eigenvalue of A v , we infer from (1.4), (3.3) and the boundedness of u that

λ v = A v z v , z v X ,X A v u, u X ,X = a(u, u) + Ω f (v 2 )u 2 M u 2 H 1 +C 1 + v 2q L 2q C 1 + v 2q L 2q . (3.30)
Using (3.3), the fact that z v L 2 = 1, and the positivity of F (t) (which implies that f (t 2 ) 0 for all t ∈ R), we obtain

λ v = a(z v , z v ) + Ω f (v 2 )z 2 v α 2 z v 2 H 1 - α 2 -β 0 ,
which, together with (3.30), readily leads to (3.13).

We now turn to the proof of (3.15) and (3.16). Let v ∈ X. We shall analyze each case of the alternative λ v > λ = λ u or λ v λ = λ u . In the former case, since λ v is the lowest eigenvalue of A v , we have

λ < λ v = A v z v , z v X ,X A v u, u X ,X = λ + Ω ( f (v 2 ) -f (u 2 ))u 2 .
In the latter case, we use this time the fact that λ u is the lowest eigenvalue of A u to get

λ v λ = A u u, u X ,X A u z v , z v X ,X = λ v + Ω ( f (u 2 ) -f (v 2 ))z 2 v .
Therefore, using either (3.5) with w = u (former case), or (3.8)-(3.9) with w = z v and (3.13) (latter case), we obtain that, for all v ∈ X,

|λ v -λ | C 1 + v 2q L 6 u -v L max(r,2) (X = H 1 # (Ω ) only) and |λ v -λ | 2ε +C ε 1 + v 2q L 6 u -v L max(r,2) (X = H 1 0 (Ω ) only). Since λ v z v -λ u u = A v z v -A u u, we have (λ v z v -λ u u, z v -u) = A v z v -A u u, z v -u X ,X =a(z v -u, z v -u) + Ω f (v 2 )z v (z v -u) - Ω f (u 2 )u(z v -u) =a(z v -u, z v -u) + Ω f (u 2 )(z v -u) 2 + Ω ( f (v 2 ) -f (u 2 ))z v (z v -u) = A u (z v -u), (z v -u) X ,X + Ω ( f (v 2 ) -f (u 2 ))z v (z v -u).
Hence, we have

(A u -λ u )(z v -u), (z v -u) X ,X = (λ v -λ u ) Ω z v (z v -u) + Ω ( f (u 2 ) -f (v 2 ))z v (z v -u),
which together with (3.12) implies that

γ z v -u 2 H 1 (λ v -λ u ) Ω z v (z v -u) + Ω ( f (u 2 ) -f (v 2 ))z v (z v -u) = 1 2 (λ v -λ u ) z v -u 2 L 2 + Ω ( f (u 2 ) -f (v 2 ))u(z v -u) + Ω ( f (u 2 ) -f (v 2 ))(z v -u) 2 .
To conclude the argument, we need again to distinguish the two cases X = H 1 # (Ω ) and X = H 1 0 (Ω ). In the former case, we can use (3.5), (3.8), (3.13) and (3.15) to get

γ z v -u 2 H 1 C 1 + v 2q L 6 z v -u 2 L 2 + z v -u L 6 u -v L max(r,2) +C u -v L 2 z v -u L 3 z v -u L 6 C 1 + v 2q L 6 z v -u H 1 u -v L max(r,2) , (3.31)
where we have used that

z v -u L 2 z v L 2 + u L 2 = 2.
In the latter case, from (3.5), (3.9) and (3.16) we have

γ z v -u 2 H 1 1 2 2ε +C ε 1 + v 2q L 6 u -v L max(r,2) z v -u 2 L 2 +C 1 + v 2q L 6 z v -u L 6 u -v L r +ε z v -u 2 L 2 +C ε u -v L 2 z v -u L 3 z v -u L 6 . (3.32)
We can choose ε = γ/4 and get (3.14). This completes the proof.

Proof of Lemma 3.5. We first notice that if λ 2,v λ 2,u , then

λ 2,v -λ v g + λ -λ v ,
so that λ 2,vλ v g/2 follows from (3.15) and (3.16) provided uv H 1 being small enough.

Let us now deal with the case where λ 2,v < λ 2,u . Since

λ 2,u = A u z 2,u , z 2,u C, we have ∀v ∈ X, λ 2,v < C.
On the other hand, using again (3.3), we get

∀v ∈ X, λ 2,v = A v z 2,v , z 2,v α 2 z 2,v 2
H 1 - α 2 -β 0 .
Hence, there exists a constant

C ∈ R + such that ∀v ∈ X, z 2,v H 1 C. (3.33) We now decompose z 2,v as z 2,v = (u, z 2,v ) L 2 u + α 2,v z ⊥ 2,v with α 2,v 0 and z ⊥ 2,v ∈ u ⊥ such that z ⊥ 2,v L 2 = 1. We have α 2 2,v = 1 -(u, z 2,v ) 2 L 2 = 1 -(u -z v , z 2,v ) 2 L 2 1 -u -z v 2 L 2 . (3.34)
We deduce from (3.14) that there exists 0 < η 0 1 such that

∀v ∈ B u,η 0 , u -z v H 1 1/2,
where B u,η 0 is the ball in H 1 with center u and radius η 0 . It then follows from (3.34) that

∀v ∈ B u,η 0 , 1 α 2 2,v 1 + 2 u -z v 2 L 2 .
As λ 2,u is the smallest eigenvalue of A u in u ⊥ , we obtain, using again (3.34) and the above estimate,

∀v ∈ B u,η 0 , λ 2,u A u z ⊥ 2,v , z ⊥ 2,v = 1 α 2 2,v A u (z 2,v -(u, z 2,v ) L 2 u, z 2,v -(u, z 2,v ) L 2 u = 1 α 2 2,v A u z 2,v , z 2,v -λ (u, z 2,v ) 2 L 2 1 α 2 2,v A u z 2,v , z 2,v + |λ | u -z v 2 L 2 = 1 α 2 2,v A v z 2,v , z 2,v + Ω f (u 2 ) -f (v 2 ) z 2 2,v + |λ | u -z v 2 L 2 1 + 2 u -z v 2 L 2 λ 2,v + |λ | u -z v 2 L 2 + Ω f (u 2 ) -f (v 2 ) z 2 2,v , hence λ 2,v -λ 2,u -(λ v -λ )-2λ 2,v u-z v 2 L 2 -1 + 2 u -z v 2 L 2 |λ | u -z v 2 L 2 + Ω f (u 2 ) -f (v 2 ) z 2 2,v .
Therefore, for any v ∈ B u,η 0 , we have

λ 2,v -λ v g -(λ v -λ ) -2λ 2,v u -z v 2 L 2 -1 + 2 u -z v 2 L 2 |λ | u -z v 2 L 2 + Ω f (u 2 ) -f (v 2 ) z 2 2,v . (3.35)
The existence of some 0 < η η 0 such that λ 2,vλ v g/2 for all v ∈ B u,η easily follows from (3.14), (3.15), (3.16), (3.33) and (3.35) also in the case where λ 2,v < λ 2,u .

Proof of Proposition 3.1. Let 0 < η 1 be as in Lemma 3.5, v ∈ X such that uv H 1 η, and w ∈ X such that w L 2 = 1 and (w, z v ) L 2 0. Note that v H 1 u H 1 + 1. Using (3.3), (3.13) and the fact that f is non-negative on R + , we have

∀z ∈ X, (A v -λ v )z, z X ,X =a(z, z) + Ω f (v 2 )z 2 -λ v z 2 L 2 α 2 ∇z 2 L 2 -(λ v + β 0 ) z 2 L 2 α 2 ∇z 2 L 2 -β z 2 L 2 ,
where the constant β is independent of v and z. In particular,

(A v -λ v )(w -z v ), (w -z v ) X ,X α 2 ∇(w -z v ) 2 L 2 -β w -z v 2 L 2 . (3.36)
From Lemma 3.5, we see that for all z ∈ X such that (z, z v ) 0,

(A v -λ v )z, z X ,X g 2 z -(z v , z) L 2 z v 2 L 2 = g 2 ( z 2 L 2 -|(z v , z) L 2 | 2 ).
Therefore, we have

(A v -λ v )(w -z v ), (w -z v ) X ,X g 2 ( w -z v 2 L 2 -|(z v , w -z v )| 2 ) g 2 w -z v 2 L 2 -(1 -(z v , w)) = g 4 w -z v 2 L 2 .
(3.37)

Combining (3.36) and (3.37) provides the lower bound of (3.17). We get the upper bound from the following estimate

(A v -λ v )(w -z v ), (w -z v ) X ,X =a(w -z v , w -z v ) + Ω f (v 2 )(w -z v ) 2 -λ v Ω (w -z v ) 2 C w -z v 2 H 1 ,
where we have used (3.3), (3.4) and (3.13).

3.2 Basic error analysis of scheme 1 LEMMA 3.6 Let u δ c δ f be a solution of (2.1). Under assumptions (1.2)-(1.6), we have

lim 0<δ f δ c →0 u -u δ c δ f H 1 = 0.
Proof. In this proof, C and C ε are constants independent of δ c and δ f . We first notice that

u -u δ c δ f H 1 u -u δ c 0 H 1 + u δ c 0 -u δ c δ f H 1 . (3.38) We know from (1.11) that u δ c H 1 u H 1 + 1 for all δ c > 0 small enough. Using (3.14) with v = u δ c (so that λ v = λ δ c 0 and z v = u δ c 0 ), we obtain u -u δ c 0 H 1 C u -u δ c H 1 ,
which together with (1.11) with δ = δ c implies lim

δ c →0 u -u δ c 0 H 1 = 0. (3.39)
For each δ f > 0, let Π δ f : X → X δ f be the orthogonal projection on X δ f for the H 1 -scalar product: for any

w ∈ X, w -Π δ f w H 1 = min v δ f ∈X δ f w -v δ f H 1 .
Again from (1.11), for any η > 0, there exists δ 0 c > 0 such that for all 0 < δ c δ 0 c , uu δ c H 1 η. Assuming that (u δ c δ f , u δ c 0 ) L 2 0, we deduce from Proposition 3.1 that for all 0 < δ f δ c δ 0 c ,

u δ c δ f -u δ c 0 2 H 1 c -1 0 (A u δc -λ δ c 0 )(u δ c δ f -u δ c 0 ), (u δ c δ f -u δ c 0 ) X ,X =c -1 0 A u δc u δ c δ f , u δ c δ f X ,X -A u δc u δ c 0 , u δ c 0 X ,X c -1 0   A u δc Π δ f u δ c 0 Π δ f u δ c 0 L 2 , Π δ f u δ c 0 Π δ f u δ c 0 L 2 X ,X -A u δc u δ c 0 , u δ c 0 X ,X   =c -1 0 (A u δc -λ δ c 0 ) Π δ f u δ c 0 Π δ f u δ c 0 L 2 -u δ c 0 , Π δ f u δ c 0 Π δ f u δ c 0 L 2 -u δ c 0 X ,X c -1 0 C 0 Π δ f u δ c 0 Π δ f u δ c 0 L 2 -u δ c 0 2 H 1 c -1 0 C 0 1 + Π δ f u δ c 0 H 1 Π δ f u δ c 0 L 2 2 Π δ f u δ c 0 -u δ c 0 2 H 1 c -1 0 C 0 1 + u δ c 0 H 1 Π δ f u δ c 0 L 2 2 Π δ f u δ c 0 -u δ c 0 2 H 1 .
Since lim

δ f →0 u δ c 0 -Π δ f u δ c 0 L 2 = 0 and u δ c 0 L 2 = 1, there exists δ 0 f > 0 such that ∀0 < δ f δ 0 f , Π δ f u δ c 0 L 2 1 2 . It follows that for 0 < δ f δ c δ 0 f , lim δ f →0 u δ c δ f -u δ c 0 H 1 = 0,
which, together with (3.38) and (3.39), leads to the desired result.

LEMMA 3.7 Let P δ f : u ⊥ → u ⊥ ∩ X δ f be the projection operator defined by

∀w δ f ∈ u ⊥ ∩ X δ f , ∀v ∈ u ⊥ , (v -P δ f v, w δ f ) H 1 = 0. We have v -P δ f v H 1 C min v δ f ∈X δ f v -v δ f H 1 .
Proof. For any v ∈ u ⊥ , we have

Π δ f v - (Π δ f v, u) L 2 Π δ f u (Π δ f u, u) L 2 ∈ u ⊥ ∩ X δ f , so that v -P δ f v H 1 v -Π δ f v + (Π δ f v, u) L 2 Π δ f u (Π δ f u, u) L 2 H 1 v -Π δ f v H 1 + (Π δ f v -v, u) L 2 Π δ f u (Π δ f u, u) L 2 H 1 1 + Π δ f u H 1 (Π δ f u, u) L 2 v -Π δ f v H 1 C min v δ f ∈X δ f v δ f -v H 1 . (3.40)
This completes the proof.

In order to state the main result of this section, we need to introduce the following object: for all v ∈ L 2 (Ω ), we denote by ψ v ∈ u ⊥ the unique solution to the adjoint problem: find

ψ v ∈ u ⊥ such that ∀w ∈ u ⊥ , (A u -λ )ψ v , w X ,X = (v, w) L 2 .
(3.41)

The existence and uniqueness of the solution to (3.41) is a straightforward consequence of (3.11) and Lax-Milgram lemma. It follows from (3.11) that

∀v ∈ L 2 (Ω ), ψ v H 1 M -1 2 v X M -1 2 v L 2 .
THEOREM 3.2 Under assumptions (1.2)-(1.6), there exist δ 1 > 0 and C ∈ R + such that for all 0

< δ f δ c δ 1 , γ 2 u -u δ c δ f 2 H 1 E(u δ c δ f ) -E(u) C u -u δ c δ f 2 H 1 (3.42)
and for all r such that 6 5 r = 6 5-2q < 6 :

u -u δ c δ f H 1 C min v δ f ∈X δ f u -v δ f H 1 + u -u δ c L r + u -u δ c δ f L r , (3.43) |λ -λ δ c δ f | C u -u δ c δ f 2 H 1 + u -u δ c L r + u -u δ c δ f L r , (3.44) u -u δ c δ f 2 L 2 = Ω f (u 2 δ c ) -f (u 2 ) u δ c δ f P δ f ψ + (λ -λ δ c δ f ) Ω (u δ c δ f -u)P δ f ψ + (A u -λ )(ψ -P δ f ψ), (u -u δ c δ f ) X ,X + 1 4 u -u δ c δ f 4 L 2 , (3.45)
where we have set ψ = ψ u-u δc δ f .

Proof. Let us recall that for all w ∈ X such that w L 2 = 1,

E(w) -E(u) = 1 2 (A u -λ )(w -u), (w -u) X ,X + 1 2 Ω F(w 2 ) -F(u 2 ) -f (u 2 )(w 2 -u 2 ) . (3.46)
This equality, referred to as (32) in Cancès et al. (2010), can be derived easily. Using (3.12) and the convexity of F, we obtain, for

w = u δ c δ f , E(u δ c δ f ) -E(u) γ 2 u -u δ c δ f 2 H 1 ,
that is the lower bound in (3.42). We now observe that

F(w 2 ) -F(u 2 ) -f (u 2 )(w 2 -u 2 ) =(w 2 -u 2 ) 1 0 f (u 2 + t(w 2 -u 2 )) -f (u 2 ) dt =(w 2 -u 2 ) 2 1 0 t 0 f (u 2 + s(w 2 -u 2 )) ds dt =(w 2 -u 2 ) 2 1 0 (1 -s) f (u 2 + s(w 2 -u 2 )) ds. (3.47)
We are led to split the domain Ω into four parts

Ω 1 = {x ∈ Ω , u(x) < |w(x)| 2 }, Ω 2 = {x ∈ Ω , |w(x)| 2 u(x) < |w(x)|}, Ω 3 = {x ∈ Ω , |w(x)| u(x) < 2|w(x)|} and Ω 4 = {x ∈ Ω , u(x) 2|w(x)
|}, where we remark that, over Ω 2 ∪ Ω 3 , |w(x)| 2 u L ∞ . Hence, from the assumption (1.6) made on F = f , we deduce that u 2 + s(w 2u 2 ) f (u 2 + s(w 2u 2 )) is bounded over Ω 2 ∪ Ω 3 by a constant (say C 3 ).

We infer from (3.10), (3.46) and (3.47) that

E(w) -E(u) M 1 2 w -u 2 H 1 + C 3 2 Ω 2 ∪Ω 3 w 2 -u 2 2 1 0 1 -s u 2 + s(w 2 -u 2 ) ds + 1 2 Ω 1 ∪Ω 4 w 2 -u 2 1 0 f (u 2 + t(w 2 -u 2 )) -f (u 2 ) dt M 1 2 w -u 2 H 1 + C 3 2 Ω 2 u 2 -w 2 -w 2 ln u 2 w 2 + C 3 2 Ω 3 u 2 -w 2 + w 2 ln w 2 u 2 + 1 2 Ω 1 u 2 -w 2 1 0 f (u 2 ) -f (u 2 + t(w 2 -u 2 )) dt + 1 2 Ω 4 w 2 -u 2 1 0 f (u 2 + t(w 2 -u 2 )) -f (u 2 ) dt .
Using thatln(1a) a + 2a 2 , for any a, 0 a 3/4, we first get that, on Ω 2 ,

(u 2 -w 2 ) -w 2 ln u 2 w 2 = (u 2 -w 2 ) -w 2 ln 1 - w 2 -u 2 w 2 2 (w 2 -u 2 ) 2 w 2 8(|w| -u) 2 8(w -u) 2 .
Then, using that ln(1a) -a andln(1a) 4a , for any a, 0 a 3/4, we get that, on Ω 3 ,

u 2 -w 2 + w 2 ln w 2 u 2 = u 2 -w 2 + u 2 ln 1 - u 2 -w 2 u 2 + (w 2 -u 2 ) ln 1 - u 2 -w 2 u 2 4 (w 2 -u 2 ) 2 u 2 16(|w| -u) 2 16(w -u) 2 .
Finally, the facts that f is positive, increasing, and that 0

|w 2 -u 2 | 3(|w| -u) 2 3(w -u) 2 on Ω 1 ∪ Ω 4 , we get E(w) -E(u) M 1 2 w -u 2 H 1 + 8C 3 Ω 2 ∪Ω 3 (u -w) 2 + 3 Ω 1 f (u 2 ) (w -u) 2 + 3 Ω 4 f (w 2 ) (w -u) 2 M 1 2 w -u 2 H 1 + 8C 3 Ω 2 ∪Ω 3 (u -w) 2 + 3 Ω 1 f (u 2 ) (w -u) 2 +C Ω 4 (1 + |w| 2q ) (w -u) 2 . Taking w = u δ c δ f
, we obtain

E(u δ c δ f ) -E(u) M 2 u -u δ c δ f 2 H1 ,
where the constant M 2 depends on the

H 1 -norm of u δ c δ f
, which is itself uniformly bounded when 0 <

δ f δ c δ 1 .
λ δ c δ f -λ = (A u -λ )(u -u δ c δ f ), (u -u δ c δ f ) X ,X + Ω f (u 2 δ c ) -f (u 2 ) (u δ c δ f ) 2 . (3.48)
Using on the one hand (3.4) with

w = u δ c δ f + u δ c , z = u δ c δ f
u δ c and both v = u, and v = u δ c , and on the other hand (3.6) with v = u δ c , we get

Ω f (u 2 δ c ) -f (u 2 ) (u δ c δ f ) 2 = Ω f (u 2 δ c ) -f (u 2 ) ((u δ c δ f ) 2 -u 2 δ c ) + Ω f (u 2 δ c ) -f (u 2 ) u 2 δ c C u -u δ c L r + u -u δ c δ f L r .
Therefore, we obtain that for all 0 < δ f δ c δ 1 ,

|λ δ c δ f -λ | C u -u δ c δ f 2 H 1 + u -u δ c L r + u -u δ c δ f L r . (3.49) Let us now estimate u -u δ c δ f H 1 . We have u -u δ c δ f H 1 u -u δ c * δ f H 1 + u δ c * δ f -u δ c δ f H 1 , (3.50)
where

u δ c * δ f is defined by u δ c * δ f = u δ c δ f + 1 - Ω uu δ c δ f u. (3.51) It is easy to see that v := u -u δ c * δ f ∈ u ⊥ and u δ c * δ f -u δ c δ f = 1 2 u u -u δ c δ f 2 L 2 . (3.52)
We then only need to estimate uu δ c * δ f H 1 . With the previous notation, we have

v -Π δ f v H 1 v H 1 and v -Π δ f v =u -u δ c δ f -u + u Ω uu δ c δ f -Π δ f u -u δ c δ f -u + u Ω uu δ c δ f = -u δ c δ f + u Ω uu δ c δ f -Π δ f -u δ c δ f + u Ω uu δ c δ f = -u δ c δ f + u Ω uu δ c δ f + u δ c δ f -Π δ f u Ω uu δ c δ f =(u -Π δ f u) Ω uu δ c δ f . Hence, v -Π δ f v H 1 c u -Π δ f u H 1 . (3.53)
Due to (3.11) and (3.52), we have

M 2 v 2 H 1 (A u -λ )v, v X ,X = (A u -λ )(v -Π δ f v), v X ,X + (A u -λ )Π δ f v, u δ c δ f -u δ c * δ f X ,X + (A u -λ )Π δ f v, u -u δ c δ f X ,X = (A u -λ )(v -Π δ f v), v X ,X - 1 2 u -u δ c δ f 2 L 2 (A u -λ )Π δ f v, u X ,X + (A u -λ )Π δ f v, u -u δ c δ f X ,X .
(3.54)

For any v δ f ∈ X δ f , we have

(A u -λ )(u -u δ c δ f ), v δ f X ,X = -(A u -λ )u δ c δ f , v δ f X ,X = Ω f (u 2 δ c ) -f (u 2 ) u δ c δ f v δ f -(A u δc -λ )u δ c δ f , v δ f X ,X = Ω f (u 2 δ c ) -f (u 2 ) u δ c δ f v δ f + (λ -λ δ c δ f ) Ω u δ c δ f v δ f = Ω f (u 2 δ c ) -f (u 2 ) (u δ c δ f -u)v δ f + Ω f (u 2 δ c ) -f (u 2 ) uv δ f + (λ -λ δ c δ f ) Ω u δ c δ f v δ f , (3.55)
which together with (3.4) and (3.5) implies that for any v δ f ∈ X δ f ,

(A u -λ )(u -u δ c δ f ), v δ f X ,X C u -u δ c δ f L r + u -u δ c L r v δ f H 1 + |λ -λ δ c δ f |. (3.56)
From (3.54) and (3.56), we have

v|| H 1 C v -Π δ f v H 1 + 1 2 u -u δ c δ f 2 L 2 + u -u δ c δ f L r + u -u δ c L r + |λ -λ δ c δ f | . (3.57)
Therefore, for all 0 < δ f δ c δ 1 , we get from (3.49), (3.50), (3.52), (3.53) and (3.57) that

u -u δ c δ f H 1 C( u -Π δ f u H 1 + u -u δ c δ f L r + u -u δ c L r ).
u -u δ c δ f 2 L 2 = Ω (u -u δ c δ f )(u -u δ c * δ f ) + Ω (u -u δ c δ f )(u δ c * δ f -u δ c δ f ) = Ω (u -u δ c δ f )(u -u δ c * δ f ) + 1 2 u -u δ c δ f 2 L 2 Ω u(u -u δ c δ f ) = Ω (u -u δ c δ f )(u -u δ c * δ f ) + 1 4 u -u δ c δ f 4 L 2 = (A u -λ )ψ, (u -u δ c * δ f ) X ,X + 1 4 u -u δ c δ f 4 L 2 . (3.58)
We have

(A u -λ )ψ, (u -u δ c * δ f ) X ,X = (A u -λ )ψ, (u -u δ c δ f ) X ,X -1 - Ω uu δ c δ f (A u -λ )ψ, u X ,X = (A u -λ )ψ, (u -u δ c δ f ) X ,X ,
which, together with (3.55), (3.58), and the fact that

P δ f ψ ∈ u ⊥ yields u -u δ c δ f 2 L 2 = (A u -λ )ψ, u -u δ c δ f X ,X + 1 4 u -u δ c δ f 4 L 2 = (A u -λ )(u -u δ c δ f ), P δ f ψ X ,X + (A u -λ )(ψ -P δ f ψ), (u -u δ c δ f ) X ,X + 1 4 u -u δ c δ f 4 L 2 = Ω f (u 2 δ c ) -f (u 2 ) u δ c δ f P δ f ψ + (λ -λ δ c δ f ) Ω (u δ c δ f -u)P δ f ψ + (A u -λ )(ψ -P δ f ψ), (u -u δ c δ f ) X ,X + 1 4 u -u δ c δ f 4 L 2 ,
which proves (3.45).

Spectral Fourier discretization

In this section, we consider Ω = (0, 2π) d with d = 1, 2, 3 and X = H 1 # (Ω ), and we make the following assumptions:

V ∈ H σ # (Ω ) for some σ > d/2, (4.1) the function F satisfies (1.4)-(1.6) and is in C [σ ]+2,σ -[σ ]+ε ((0, +∞), R). (4.2)
The positive solution u to (1.1), which satisfies the elliptic equation

-∆ u +Vu + f (u 2 )u = λ u, then is in H σ +2 #
(Ω ) and is bounded away from 0. A natural discretization of (1.1) consists in using a Fourier basis. Denoting for any k ∈ Z d by e k (x) = (2π) -d/2 e ik•x , we have for all v ∈ L 2 (Ω ),

v(x) = ∑ k∈Z d vk e k (x),
where vk is the kth Fourier coefficient of v:

vk := Ω v(x)e k (x)dx = (2π) -d/2 Ω v(x)e -ik•x dx.
The Fourier spectral approximation of the solution to (1.1) is based on the choice

X M =    ∑ k∈Z d ,|k| * M c k e k , ∀k, c k = c -k    ,
where |k| * denotes either the l 2 -norm or the l ∞ -norm of the wave vector k.

Endowing H ρ # (Ω ) with the norm defined by

v H ρ = ∑ k∈Z d (1 + |k| 2 * ) ρ | vk | 2 1/2
, we obtain that for all τ ∈ R, and all v ∈ H τ # (Ω ), the best approximation of v in H ρ # (Ω ) for any ρ τ is

Π M v = ∑ k∈Z d ,|k| * M vk e k .
For all real numbers ρ and τ with ρ τ, we have

∀v ∈ H τ # , v -Π M v H ρ 1 M τ-ρ v H τ .
(4.3)

In this section, we take δ c = M -1 and δ f = N -1 (M N), and u δ c , u δ c 0 and u δ c δ f are denoted as u M , u M 0 and u M N , respectively. It is easy to see that u M 0 ∈ H σ +2 # (Ω ). Aligning the functions u M , u M 0 and u M N in such a way that (u M 0 , u M κ ) L 2 0 for κ = M, N, and using (4.3), we obtain

u M 0 -Π κ u M 0 H 1 1 κ σ +1 u M 0 H σ +2 .
It therefore follows from Lemma 3.6 that

lim 0<M N→∞ u -u M N H 1 = 0. It is proved in Cancès et al. (2010) that u M converges to u in H σ +2 # (Ω ).
In particular, u/2 u M 2u on Ω for M large enough. Besides, u M N is solution to the elliptic equation

-∆ u M N + Π N (Vu M N + f (u 2 M )u M N ) = λ M N u M N .
Thus u M N is uniformly bounded in H 2 # (Ω ), hence in L ∞ (Ω ), and

∆ (u M N -u) =Π N V (u M N -u) + f (u 2 M )u M N -f (u 2 )u + (Π N -I)(Vu + f (u 2 )u) -λ M N (u M N -u) -(λ M N -λ )u. (4.4) Since u M is bounded in L ∞ (Ω ), F ∈ C 1 ([0, ∞), R)∩C 2 ((0, ∞), R
), lim M→0 u-u M H 1 = 0, and lim 0<M N→∞ uu M N H 1 = 0, the right hand side of the above equality converges to 0 in L 2 # (Ω ), which implies that (u M N ) N∈N actually converges to u in H 2 # (Ω ). With a bootstrapping argument, we also deduce from (4.4) that u M N converges to u in H σ +2 # (Ω ). Besides, the unique solution to (3.41) solves the elliptic equation

-∆ ψ v + (V + f (u 2 ) -λ )ψ v = v -(v, u) L 2 u, from which we infer that ψ u-u M N belongs to H σ +2 # (Ω ) and ψ u-u M N H 2 C u -u M N L 2 . Hence ψ u-u M N -Π N ψ u-u M N H 1 1 N ψ u-u M N H 2 C N u -u M N L 2 .
THEOREM 4.1 Under assumptions (4.1)-(4.2), there exists C ∈ R + such that for all N ∈ N,

u -u M H τ C M σ +2-τ for all -σ τ < σ + 2, (4.5) E(u M N ) -E(u) C 1 M σ +3 + 1 N σ +1 2 , (4.6) u -u M N H 1 C 1 M σ +3 + 1 N σ +1 , (4.7) u -u M N L 2 C 1 M min{σ +4,2(σ +1)} + 1 N 1 M σ +3 + 1 N σ +1 , (4.8) |λ -λ M N | C 1 M 2(σ +1) + 1 N 2(σ +1) .
(4.9)

Proof. The proof of (4.5) is detailed in Cancès et al. (2010) 1 .

Let us first come back to (3.48), which we rewrite as,

λ M N -λ = (A u -λ )(u -u M N ), u -u M N X ,X + Ω w M,N (u M -u), (4.10) with w M,N = f (u 2 M ) -f (u 2 ) u 2 M -u 2 (u M + u)(u M N ) 2 = (u M + u)(u M N ) 2 1 0 f (u 2 + t(u 2 M -u 2 )) dt,
where the argument of f , namely (u 2 + t(u 2 Mu 2 )), belongs to H σ +2 # (Ω ) for any t ∈ (0, 1). As, for M large enough, u/2 u M 2u on Ω for M, we also have

u 2 /4 (u 2 + t(u 2 M -u 2 )) 4u 2 . As f ∈ C [σ ]+1,σ -[σ ]+ε ([ u 2 L ∞ /4, 4 u 2 L ∞ ], R)
, we obtain that w M,N is uniformly bounded in H σ # (Ω ) (at least for N large enough). We therefore infer from (4.10) that, for M large enough,

|λ M N -λ | C( u -u M N 2 H 1 + u -u M H -σ ).
(4.11)

1 Note that, as already observed in Cancès et al. (2010), it follows from the fact that the continuous solution u and the discrete ones u M and u M N are bounded away from zero, the assumption that there exist 1 < r 2 and 0 s 5r such that Cancès et al. (2010) is actually not necessary and is thus not made here.

∀R > 0, ∃C R ∈ R + s.t. ∀0 < t 1 R, ∀t 2 ∈ R, |F (t 2 2 )t 2 -F (t 2 1 )t 2 -2F (t 2 1 )t 2 1 (t 2 -t 1 )| C R (1 + |t 2 | s )|t 2 -t 1 | r made in

Finite-element discretization

In this section, we assume that Ω is a rectangular brick of R d with d = 1, 2, 3 and X = H 1 0 (Ω ). By elliptic regularity, the positive solution u to (1.1), which satisfies the elliptic equation

-∆ u +Vu + f (u 2 )u = λ u, is in H 2 (Ω ) ∩ H 1 0 (Ω ) whenever V is in L 2 (Ω ), and is in H 3 (Ω ) ∩ H 1 0 (Ω ) whenever V is in H 1 (Ω )
(use an extension-by-symmetry argument in order to check that there are no vertex or edge singularities, and the fact that

f (u 2 )u 2 ∇u is in L 2 (Ω ) whenever u is in H 2 (Ω )).
Considering a family of quasi-uniform triangulations (T δ ) δ of Ω , we introduce the coarse (X p H ) H (associated to the triangulations indexed by δ = H) and fine (X h ) h (associated to the triangulations indexed by δ = h) finite element subspaces of H 1 0 (Ω ) such that :

• X k δ = v ∈ H 1 0 (Ω ), ∀ K δ ∈ T δ , v |K δ ∈ P k (K δ ) , • k = p or (p, = 1 or 2) • δ = H or h, with 0 < h H, • T h is a sub-triangulation of T H .
As usual, H (resp. h) denote the maximum of the diameters

H K , K ∈ T H (resp. h K , K ∈ T h ).
We denote by I δ ,k the interpolation operator on X k δ . The following estimates are classical (see e.g [START_REF] Bernardi | Discrétisations variationnelles de problèmes aux limites elliptiques[END_REF]; [START_REF] Ciarlet | Handbook of numerical analysis[END_REF]; [START_REF] Ern | Theory and practice of finite elements[END_REF]).

LEMMA 5.1 For any integer n, 0 n k + 1, and for all r and q, 1 r q < +∞, such that ∀K δ ∈ T δ , W n,r (K δ ) is included in C 0 (K δ ), there exists a positive constant c depending only on n, r and q such that, for any function v of W n,r (Ω ), we have :

v -I δ ,k v L ∞ c δ n-d r |v| W n,r , (5.1) v -I δ ,k v W 1,q c δ n-1-d r + d q |v| W n,r .
(5.2) LEMMA 5.2 There exists a positive constant c independent of δ such that, for any v δ ∈ X k δ we have :

v δ L ∞ c ζ (δ ) v δ H 1 where ζ (δ ) =    c v δ H 1 for d=1, c(1 + | log δ |) v δ H 1 for d =2, c δ -1 2 v δ H 1 for d =3. Let u δ ,k be a solution of the minimization problem inf E(v δ ,k ), v δ ,k ∈ X k δ , Ω v 2 δ ,k = 1 such that (u, u δ ,k ) L 2 0. Let
( f (u 2 H,p ) -f (u 2 ))u H,p h, v h, = Ω w H,p (u H,p -u)u H,p h, v h, C w H,p L 6 u H,p h, L 6 v h, L 6 u H,p -u L 2 C u H,p -u L 2 v h, H 1 .
(5.25)

We then infer from (3.45) and (5.25) that

u H,p h, -u 2 L 2 =(λ -λ H,p h, ) Ω (u H,p h, -u)P h, ψ + Ω f (u 2 H,p ) -f (u 2 ) u H,p h, P h, ψ + (A u -λ )(u -u H,p h, ), ψ -P h, ψ X ,X + 1 4 u -u H,p h, 4 L 2 |λ -λ H,p h, | u -u H,p h, 2 L 2 + u -u H,p L 2 u -u H,p h, L 2 + h u -u H,p h, H 1 u -u H,p h, L 2 + 1 4 u -u H,p h, 4 L 2 , Hence u H,p h, -u L 2 C( u -u H,p L 2 + h u -u H,p h, H 1 ).
(5.26)

Inserting this result in (3.43) gives

u -u H,p h, H 1 C( u -Π h, u H 1 + u -u H,p L 2 ),
which leads to (5.16) and (5.20), and then to (5.15) and (5.19). Next, from (5.26), we further deduce (5.17) and (5.21). Finally, (5.18) and (5.22) are consequences of (5.23) and (5.24) respectively.

The effect of numerical integration

Let us now sketch the effect of a practical implementation of the method, and more precisely to the numerical integration of the nonlinear term. For simplicity, we focus on the case when A = I, with periodic boundary conditions and Ω = [0, L) 3 (L > 0). For N g ∈ N\{0}, we perform the numerical integration on the cartersian grid G N g := L N g Z 3 . We now introduce the subspace

W 1D N g = Span e ily |l ∈ 2π L Z, |l| 2π L ( N g -1 2 ) (N g odd) Span e ily |l ∈ 2π L Z, |l| 2π L ( N g
2 ) ⊕ C(e iπN g y/L+e -iπNgy/L )(N g even),

and In order to evaluate the quality of the error bounds obtained in Theorem 5.2, we have performed numerical tests with Ω = [0, 2π] 2 , f (t) = t, V (x) = x 2 + y 2 , using P 1 and P 2 finite elements. We denote the number of degrees of freedom in the coarse and fine grids by DOF H and DOF h , respectively. Fig. 5-8 show the numerical errors using P 1 finite elements for both the coarse grid and the fine grid. These figures agree with the results of Theorem 5.2, except the right figure of Fig. 5 in which the term in h dominates. Fig. 9-10 show the numerical errors using P 1 on the coarse grid and P 2 on the fine grid. Fig. 11-14 show the numerical errors using P 2 finite elements on both the coarse grid and the fine grid. Similar conclusions as for the plane wave approximation hold here which illustrate the various behaviors stated in Section 5.

W 3D N g = W 1D N g ⊗W 1D N g ⊗W 1D N g . It is then possible to define the interpolation projector I N g from C 0 # (Γ , C) onto W 3D N g by [I N g (φ )](x) = φ (x) for any x ∈ G N g . We now consider the following approximate problem inf{E M N g (v N,N g ), v N,N g ∈ X N , Ω |v N,N g | 2 = 1}, (6.1) where E M N g (v N ) = 1 2 Ω |∇v N | 2 + 1 2 Ω I N g (V )v 2 N + 1 2 Ω f (u 2 M )v 2 N .

  The detailed proof of (3.10) and (3.11) can be found in Cancès et al. (2010). Let us prove (3.12). We know from inequality (20) in Cancès et al. (2010) that there exists η

Thus ( 3

 3 .43) is proved. Let us now consider the L 2 estimate and set ψ = ψ u-u

  us recall the main result in Cancès et al. (2010) concerning the finite element discretization. THEOREM 5.1 Assume that V ∈ L 2 (Ω ), the function F satisfies (1.4), (1.5) for q = 1, and (1.6), (5.3) We refer to the proof of Theorem 3 in Cancès et al. (2010) for details. Next for any v h, ∈ X h , there holds Ω

FIG. 4 .

 4 FIG. 4. Numerical errors |λλ M N | (Fourier approximation), as functions of N (left) and M (right) (in log-log scale).

  FIG. 6. Numerical errors uu H,1h,1 L 2 (P 1 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

  FIG. 7. Numerical errors |λλ H,1h,1 | (P 1 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

  FIG.12. Numerical errors uu H,2 h,2 L 2 (P 2 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

  FIG. 13. Numerical errors |λλ H,2 h,2 | (P 2 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

  FIG. 14. Numerical errors |λλ H,2 h,2 | (P 2 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

  FIG. 9. Numerical errors uu H,1 h,2 H 1 (left) and uu H,1 h,2 L 2 (right), as functions of DOF H (in log-log scale).

		10 -2 10 -1			DOF h =82433 DOF h =328705 The line with -1 slope		10 -2 10 -1			DOF h =82433 DOF h =328705 The line with -1 slope
	||Φ-Φ H,1 h,2 || H 1	10 -3				||Φ-Φ H,1 h,2 || L 2	10 -4 10 -3
		10 -4					10 -5
		10 1 10 -5	10 2	10 3 DOF H	10 4	10 5	10 1 10 -6	10 2	10 3 DOF H	10 4	10 5
										DOF h =20737 DOF h =82433 DOF h =328705 DOF h =DOF H The line with -1 slope
						λH,1 h,1 |

FIG. 8. Numerical errors |λλ H,1

h,1 | (P 1 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

  FIG. 10. Numerical errors |λλ H,1 h,2 | (left) and |λλ H,1 h,2 | (right), as functions of DOF H (in log-log scale). FIG. 11. Numerical errors uu H,2 h,2 H 1 (P 2 finite elements), as functions of DOF h (left) and DOF H (right) (in log-log scale).

		10 -1					The line with -1 slope DOF H =353 DOF H =97		10 -2 10 -1
		10 -2							10 -3
	||u-u H,2 h,2 || H 1	10 -3						||u-u H,2 h,2 || H 1	10 -5 10 -4
		10 1 10 -5 10 -4	10 2	10 3	DOF h	10 4	10 5	10 6	10 1 10 -8 10 -7 10 -6	10 2 DOF h =82433 DOF h =328705 DOF h =DOF H The line with -3/2 slope 10 3 DOF H	10 4	10 5
											DOF h =82433 DOF h =328705 The line with -1 slope
								λH,1 h,2 |

Acknowledgements

Financial support from the ANR grant Manif and from the French state funds managed by CALSIMLAB and the ANR within the Investissements d'Avenir programme under reference ANR-11-IDEX-0004-02, are gratefully acknowledged.

We now make use of (3.45), which reads here as

L 2 , (4.12)

with ψ = ψ u-u M N . Reasoning as above, we obtain that the sequence

is uniformly bounded in H σ # (Ω ) (at least for M large enough). Setting σ * = min{σ , 2}, we have

which, together with (3.40), (4.3) and (4.12), implies

Therefore, we have

Let v := u M * Nu, with u M * N being defined as in (3.51). We deduce from (3.54) and (3.55) that

We also have

From (3.10), (4.14) and (4.15), we obtain

which together with (4.3), (4.11) and (4.13), completes the proof of (4.6)-(4.9). and there exist 1 < r 2 and 0 s + r 3 such that ∀R > 0, ∃C R ∈ R + for which

(5.4)

Then there exist δ 0 > 0 and C ∈ R + such that for all 0

(5.8)

If, in addition, V ∈ H 1 (Ω ), (5.4) is satisfied for r = 2, F ∈ C 3 ((0, +∞), R), (5.9)

F ∈ C 3 ((0, +∞), R) and F (t)t 1/2 and F (t)t 3/2 are locally bounded in [0, +∞), (5.10) then there exist δ 0 > 0 and C ∈ R + such that for all 0 < δ δ 0 ,

(5.14)

In this section, we will take X δ c = X p H , and X δ f = X h . Let u H,p 0 ∈ X be the unique solution of

and u δ c δ f = u H,p h, be the solution of the following lineared eigenvalue problem (Two-grid scheme 1): find

Proof. To establish this result, we first remark that

From (5.1), we have lim

Using (5.2) and Lemma 5.2, we obtain

This completes the proof.

The following theorem states the behavior of the two-grid approach in the finite element context.

THEOREM 5.2 If (5.3) and (5.4) are satisfied, then there exist c ∈ R + and h 0 ∈ R + such that for all 0 < h, H h 0 , we have :

(5.18)

If, in addition, (5.9) and (5.10) are satisfied, then there exist c ∈ R + and h 0 ∈ R + such that for all 0 < h, H h 0 and p, = 1 or 2, we have :

(5.22)

Proof. We follow step by step the same lines as in Theorem 4.1. The analysis will be done gradually under the various regularity assumptions on F. We first start with the analysis of the eigenvalues. Proceeding as in (4.10), we get

We have already derived from this equality the generic estimate (3.44), which for q = 1 and r = 2, gives

(5.23)

If V and F satisfies the additional regularity assumptions (5.9)-(5.10), then w H,h belongs to H 1 (Ω ) and

(5.24)

Let us denote by u M N,N g a solution to (6.1) such that (u M N,N g , u M N ) 0. Then u M N,N g satisfies the following Euler-Lagrange equation:

LEMMA 6.1 There exists a positive constant M 3 such that

where

Proof. It is easy to see that

Note that λ M N is the variational approximation in X N of some eigenvalue of A u M . As (u M ) M∈N converges to u in L ∞ (Ω ), A u M -A u converges to 0 in operator norm. So the n th eigenvalue of A u M converges to the n th eigenvalue of A u when N goes to infinity, the convergence being uniform in n. As the sequence (λ M N ) N∈N converges to λ , the non-degenerate ground state eigenvalue of A u , we obtain that for N large enough, λ M N is the non-degenerate ground state eigenvalue of A u M in X N . We conclude the proof by proceeding as in Lemma 3.1.

Following step by step the same lines as in Cancès et al. (2010), we can prove the following result (we omit the details here for the sake brevity). THEOREM 6.1 Assume that V ∈ H σ # (Ω ) for some σ > d/2 and that the function F satisfies (1.4)-(1.6) and is in C [σ ]+2,σ -[σ ]+ε ((0, +∞), R). Then there exists C > 0 such that for all N ∈ N,

Numerical examples

In order to evaluate the quality of the error bounds obtained in Theorem 4.1, we have performed numerical tests with Ω = (0, 2π) and f (t) = t. The Fourier coefficients of the potential V are given by

, from which we deduce that V ∈ H σ # (0, 2π) for all σ < 3/2. The reference values for u and λ are obtained for N = 500. We first fix M and study the behaviors of the numerical errors u

as functions of N. Let us consider for example the case when M = 10. From the left figures of Fig. 1 and Fig. 2 we can see that uu M N H 1 and uu M N L 2 decay respectively as N -2.5 and N -3.5 up to N = 40, while from N = 40 the errors decay slowly and finally reach plateaus, on which the terms in 1 M τ dominate. Then, we fix N and study the numerical errors u decay respectively as M -4.5 and M -5 before reaching plateaus. These results are in agreement with the assertions of Theorem 4.1.

The same conclusion holds for the error on the eigenvalue, both as a function of N or M. An interesting observation is that the two-grid scheme 1 actually leads to two approximations of the eigenvalue, the first one being λ δ c δ f , the second being the Rayleigh quotient (2.2). Our simulations (and this can easily be confirmed by theoretical arguments) show that the rates of convergence of these two approximations are the same. Note, however, that the accuracy is somehow better for the second approximation and that, in addition, the convergence to zero is more monotonic and smoother. We are unfortunately not able to provide an explanation of this fact.


Ke Sun

Loïc Besnard
email: loic.besnard@irisa.fr

Thierry Gautier

Optimized Distribution of Synchronous Programs via a Polychronous Model

This paper presents a distribution methodology for synchronous programs, applied in particular on programs written in the Quartz language. The given program is first transformed into an intermediate model of guarded actions. After userspecified partitioning, the generated sub-models are transformed into Signal processes. Using the multi-clock calculation model of the Signal language, an optimized data-flow network can be automatically constructed. The optimization includes reducing the communication quantity and the computation load, with no change to the interface behaviors.

Optimized Distribution of Synchronous Programs via a Polychronous Model

Ke Sun * , Loïc Besnard † , and Thierry Gautier *

I. INTRODUCTION Synchronous programming languages such as Esterel [1],

Lustre [2] or Quartz [3] are all based on the synchronous hypothesis [4]. Under this assumption, system behaviors are projected onto a discrete sequence of logical instants. As the sequence is discrete, nothing occurs between two consecutive instants. Such temporal abstraction can greatly facilitate safetycritical reactive system design. It enforces deterministic concurrency of the system: Heisenbugs (i.e., bugs that disappear when one tries to simulate/test the system) are avoided; the system behaviors become predictable. It is also the key to a straightforward translation of synchronous programs to hardware circuits [5]. Furthermore, it guarantees semantic precision by using mathematical models, such as Mealy machine, as supporting foundations. These models enable a series of efficient optimization, compilation and verification techniques to be applied on synchronous programs.

Synchronous programs are modeled on centralized architectures with zero-time communications: any signal emitted by some component is instantaneously received by others. However, in real world, safety-critical reactive systems [6] [7] practically operate over distributed architectures with delayed communications. This mismatch results in a wide gap between centralized design model and distributed implementation. To get over this, desynchronization is introduced with the following requirement: how to distribute a synchronous program while preserving the same observable behaviors, i.e., the same interface behaviors. This issue has drawn a considerable attention both in theoretical challenges and in industrial relevance.

For Esterel and Lustre, a series of distribution methods have been proposed [8]. In [9], starting from a synchronous program, a network of communicating codesign finite state machines is constructed. The distribution of synchronous programs via the modeling as finite deterministic automata [10] is presented in [11]. Based on it, the extended method in [12] focused on the automatic deduction of distributed systems from centralized synchronous circuits. Furthermore, how to distribute synchronous programs to fulfill temporal constraints is introduced in [13]. The primary goal of these distribution methods does not include an optimization of the communication and the computation, which is in the focus of this paper.

A related distribution methodology is presented in [14]. Its procedure consists of three steps: 1) the given synchronous program is compiled into an intermediate AIF model [15], which is a common intermediate format for various synchronous languages [16] [17]; 2) the data dependencies within the model are analyzed and represented in a dependency graph, then the generated graph is subsequently partitioned into subgraphs in which partitions can be made horizontal (for a pipelined execution) or vertical (for a parallel execution); 3) a synchronous elastic system [18] is synthesized by generating a distributed component for each sub-graph and establishing a communication infrastructure with synchronous elastic flow (SELF) protocol [18].

Owing to the analogy to synchronous hardware circuit, one can assume that the given synchronous program holds a master clock driving the whole program computation. We refer to such program as a mono-clocked program [19]. Furthermore, it follows the model of computation (MoC) strict synchrony: during each instant, each input channel always reads data and each output channel always produces data. This feature may cause unnecessary communications or computations. In [20], the communication quantity is reduced, according to an evaluation of communication necessity: a producer transmits a value only when it is really required for the calculation of a receiver. However, it may result in highly overestimated computations for evaluating communication necessity. To avoid this, the evaluation can be refined by constructing additional communications, while this operation may augment the communication. An optimization method on computations is presented in [21]: a condition is computed for every variable to determine whether its value is required for current or future computations; when the condition does not hold, the computation of the corresponding variable can be suppressed.

A new optimized distribution methodology is proposed, which is based on a more flexible MoC. In contrast to synchronous languages, the polychronous language Signal [22] [23] [24] is based on the MoC polychrony1 . As its name suggests, a polychronous program makes use of multiple clocks to drive its execution, the system behaviors are defined on a partially ordered set of instants. One can consider that each component in the program holds its own master clock, and there is no longer a master clock for the whole program.

We refer to such program as a multi-clocked program. If there is no relation between the master clocks, the multi-clocked program does not follow a linear timeline, but a nonlinear one. Nonlinear instants t 1 , t 2 are such that there is no order relation between them. The system behaviors over nonlinear instants occur asynchronously. Hence, the Signal language allows one to conveniently model the asynchronous communications between synchronous modules. This exactly accords with the requirement of distributed systems. In addition, nonlinear timeline implies the possibility to avoid unnecessary synchronization, thereby enabling potential optimization [25]. Due to these advantages, Signal is particularly suited as a coordination layer to finally build a data-flow network over desynchronized processing locations.

Based on the multi-clock calculation model of the Signal language, we propose a distribution methodology for the synchronous programs. First, the given synchronous program is compiled into AIF model. Second, according to the user specifications, the generated model is partitioned into submodels. After the partitioning, the generated sub-models are transformed into equivalent Signal processes. Then, the unnecessary constraints are eliminated from the processes to avoid unnecessary synchronization. Finally, within the Signal framework, the minimal frequencies of communication and computation are computed via multi-clock calculation. This operation can efficiently reduce the communication quantity and the computation load. Along this way, an optimized dataflow network over desynchronized processing locations can be constructed. Note that both the methodology in [14] and ours process AIF model. Both methods are thereby independent of particular synchronous languages.

The rest of this paper is organized as follows. Section II briefly reviews the intermediate AIF model, which serves as the starting point for the distribution procedure, and introduces the Signal language and its associated polychronous semantic model. Section III presents the user-specified model partitioning and the transformation from partitioned models to Signal processes. In Section IV, how to achieve the optimization both on communications and on computations is presented in detail. Meanwhile, our methodology is illustrated through case studies in Section V. Finally, we conclude this paper and look forward the perspectives.

II. FOUNDATIONS A. AIF Model & Synchronous Guarded Actions

To process synchronous programs, it is quite natural to compile them into intermediate models at first. In this way, the whole processing can be modularly divided into several steps and the models can be reused for different purposes [8], such as validation, comparison, model transformation and code generation. Furthermore, the processing on the intermediate model is independent of particular synchronous languages. The distribution of synchronous programs converts to two steps: 1) compilation into intermediate models; 2) distribution of intermediate models. The intermediate model devoted to the Quartz language is the Averest Intermediate Format (AIF) model [15]. The compilation of Quartz programs into AIF models has been implemented in the Quartz/Averest framework2 [26]. The essential part of AIF model is AIF model behaviors that are described by a set of synchronous guarded actions.

Guarded actions are designed in the spirit of traditional guarded commands [27], which are well-established intermediate code for the description of concurrent systems. The guarded actions are in the form of γ ⇒ A, where the Boolean condition γ is called guard and A is called action. The guards represent the complex control structure of the synchronous program. The actions represent the assignments to variables. Due to the category of assignments, each guarded action has either the form γ ⇒ x = τ (called guarded immediate action) or γ ⇒ next(x) = τ (called guarded delayed action). Once the guard is evaluated to true, the corresponding action instantaneously starts. Both kinds of assignments evaluate the expression τ at the current instant. Then, the immediate assignment x = τ instantaneously transfers the value of τ to the variable x, whereas the assignment effect of the delayed one next(x) = τ takes place at the next instant. When the value of a variable x cannot be determined by any action, its value is determined by the absence reaction. It determines the value, according to the storage type of the assigned variable: a non-memorized variable is reset to the default value (denoted as def value, the particular value depends on the data type); a memorized variable keeps its previous value, or takes its default value at the initial instant.

Guarded actions describe model behaviors via two parts: the data-flow part computes locals and outputs; the controlflow part computes labels, which denote the pause locations of control-flow. Guarded actions in the control-flow part are in the form of γ ⇒ next(l) = true, where label l is a nonmemorized Boolean local denoting a pause location. If l holds at the current instant, it means that the control-flow reached the pause location of l at the end of the previous instant, then it resumes from this location at the beginning of the current instant. Note that more than one label can hold at the same instant. This enables the description of the parallelism feature [3] of synchronous programs.

B. The Signal Language

In Section II-B1, we introduce the polychronous model as the formal basis of Signal. Then, an overview of the language is given in Section II-B2.

1) Polychronous Model: We start with the following sets: X is a countable set of variables; B = {ff, tt} is a set of Boolean data values where ff and tt respectively denote false and true; V is a non-empty set of data values, B ⊂ V; and T is a dense set equipped with a partial order relation, denoted by ≤. The elements in T are called tags. We now introduce the notion of time domain.

Definition 1 (time domain). A time domain is a partially ordered set (T , ≤) where T ⊂ T that satisfies: T is countable;

T has a lower bound 0 T for ≤, i.e., ∀t ∈ T , 0 T ≤ t; ≤ over T is well-founded; the width of (T , ≤) is finite.

(T, ≤) provides a continuous time dimension. (T , ≤) defines a discrete time dimension that corresponds to the logical instants [24], at which the presence/absence of data can be observed during the system execution. Thus, the mapping of T on T allows one to move from "abstract" descriptions to "concrete" descriptions [28].

A chain (C, ≤) ⊆ (T , ≤) is a totally ordered set of tags that admits a lower bound 0 C . The notation t + 1 means the immediate successor of a tag t in C, which satisfies ∀t ′ ∈ C, t ≤ t ′ ⇒ t + 1 ≤ t ′ . We denote the set of all chains in T by C T .

Definition 2 (event). An event on a given time domain T is a pair (t, v) ∈ T × V, which associates a tag t with a data value v.

All the events whose tags belong to the same chain, can constitute a data-flow. Formally, Definition 3 (signal). A signal s : C → V is a function from a chain of tags to a non-empty set of data values, where C ∈ C T . The domain of s is denoted by tags(s).

Two signals s 1 and s 2 are identical, denoted by s 1 = s 2 , if and only if tags(s 1 ) = tags(s 2 ) and ∀t ∈ tags(s 1 ), s 1 (t) = s 2 (t). S = ∪ C∈CT (s : C → V) is the set of signals over the time domain (T , ≤).

Definition 4 (behavior). Given a finite subset X of a countable set X of variables, a behavior over X is an injective function b : X → S. The domain of b is denoted by vars(b).

The restriction of b over a set of variables X, denoted by b| X , is the behavior defined by vars(b| X ) = X ∩ vars(b) and ∀x ∈ vars(b| X ), (b| X )(x) = b(x).

Two behaviors b 1 and b 2 are compatible, denoted by b

1 b 2 , if and only if b 1 | vars(b2) = b 2 | vars(b1) . The composition of two compatible behaviors b 1 : X 1 → S, b 2 : X 2 → S is a behavior (b 1 |b 2 ) : (X 1 ∪ X 2 ) → S defined by ∀x ∈ X 1 , (b 1 |b 2 )(x) = b 1 (x) and ∀x ∈ X 2 , (b 1 |b 2 )(x) = b 2 (x).
Definition 5 (process). A process p is a set of behaviors defined over the same domain. The domain is denoted by vars(p), which satisfies ∀b ∈ p, vars(b) = vars(p).

Definition 6 (composition of processes). The composition of two processes p 1 and p 2 , denoted as p 1 |p 2 , is a process consisting of the compositions of any two compatible behaviors b 1 ∈ p 1 and b 2 ∈ p 2 :

p 1 |p 2 = {(b 1 |b 2 )|(b 1 , b 2 ) ∈ p 1 × p 2 , b 1 b 2 }.
The notions presented in this section are sufficient to express the semantics of Signal elementary concepts within this polychronous model [24]. In the remainder of this paper, we also use them to formulate the optimization scheme.

2) An Overview of the Signal Language: Signal is a polychronous language processing unbounded series of typed values, called signals. At any tag t, a signal x (corresponding to b(x) in the polychronous model, where b is a behavior) may be present or absent: when present (i.e., t ∈ tags(b(x))), it holds some value; when absent (i.e., t / ∈ tags(b(x))), it holds no value.

The presence status of x is denoted by its associated clock x : tags(b(x)) → {tt}. Furthermore, the Signal language supports clock calculation. The basic operations contain clock union x 1 +x 2 : tags(b(x 1 )) ∪ tags(b(x 2 )) → {tt}; clock intersection x 1 * x 2 : tags(b(x 1 )) ∩ tags(b(x 2 )) → {tt}; and clock difference x 1 -x 2 : tags(b(x 1 ))\tags(b(x 2 )) → {tt}. In order to enable reasoning on clock calculation, we define 0 for the empty clock (i.e., 0 : ∅ → {tt}) and [x] (resp. There are two kinds of definitions:

• A complete definition (:=) is an equation in which the defined signal is assigned with the associated expression. For instance, y := x is a complete definition of y that is assigned with x when x is present; y is synchronized with x:

∀t ∈ tags(b(x)), b(y)(t) = b(x)(t) and y =x.

• A partial definition (::=) is an equation in which the defined signal is assigned with the associated expression when the expression is present. For instance, y ::= x is a partial definition of y that is assigned with x when x is present; when x is not present, the value of y depends on other partial definitions:

∀t ∈ tags(b(x)), b(y)(t) = b(x)(t).
According to the implied clock relations, the Signal expressions can be classified into two families: synchronous expressions (i. that the clock of x is the clock union of x 1 and x 2 , its value is equal to x 1 when x 1 is present, otherwise equal to x 2 when x 1 is absent but x 2 is present. Given the partial definitions x ::= x 1 , x ::= x 2 , where x 1 #x 2 , x is then identical to x ′ when x 1 and x 2 are absent but x and x ′ are present.

  [¬x]) for the clock at which tags a Boolean signal x holds the value true (resp. false). a) Declarative Constraints A Signal program declares constraints on the involved signals, that must be satisfied by both values and clocks. The constraints on clocks are referred to as clock relations, including synchronization relation x 1 =x 2 , i.e., tags(b(x 1 )) = tags(b(x 2 )); inclusion relation x 1 ≤x 2 , i.e., tags(b(x 1 )) ⊆ tags(b(x 2 )); and mutual exclusion relation x 1 #x 2 , i.e., tags(b(x 1 )) ∩ tags(b(x 2 )) = ∅. Besides the explicit clock relations, constraints can be implicitly declared by the Signal equations. Each equation is a definition associating one defined signal with a Signal expression built on operators over signals. The operands of the operators can be signals or expressions.

  e., all the involved signals have the same clock) and polychronous expressions (i.e., the involved signals may have different clocks). The primitive expressions include • Instantaneous function: x := f (x 1 , . . . , x n ) defines a point-wise n-ary function on sequences of values, in which all the signals x, x 1 , . . . , x n are synchronous. • Delay: x := x ′ $1 init def value defines that x andx ′ are synchronous; the current value of x is equal to the previous value of x ′ and equal to the default value def value at the initial tag.• x =x ′ • b(x)(0 Cx ) = def value • ∀t ∈ tags(b(x)), 0 Cx < t + 1 ⇒ b(x)(t + 1) = b(x ′ )(t) where C x = tags(b(x)). • Downsampling: x := x ′when b defines a downsampling of the signal x ′ that occurs only when both x′ and the downsampling condition b hold true. • x =x ′ * [b], where [b] =when b • ∀t ∈ tags(b(x)), b(x)(t) = b(x ′ )(t) • Deterministic merging: x := x 1 default x 2 defines

  • x =x 1 +x 2 • ∀t ∈ tags(b(x 1 )), b(x)(t) = b(x 1 )(t) • ∀t ∈ tags(b(x 2 ))\tags(b(x 1 )), b(x)(t) = b(x 2 )(t) • Completion: x ::= defaultvalue x ′ completes the definition of x.

From the Greek "poly chronos" to mean multiple clocks.

http://www.averest.org