Eric Razafimahazo
email: eric.razafimahazo@isae-supaero.fr

Pierre De Saqui-Sannes
email: pierre.de-saqui-sannes@isae-supaero.fr

Rob A Vingerhoeds
email: rob.vingerhoeds@isae-supaero.fr

Claude Baron
email: cbaron@laas.fr

Julien Soula
email: julien.soula@cstb.fr

Romain Mège
email: romain.mege@cstb.fr

Mastering Complexity for Indoor Inspection Drone Development

Keywords: Systems Engineering, Method, Drone, Building Indoor Inspection, BIM

Drones have successfully been used for many outdoor uses, notably inspection of buildings. Using drones for indoor inspection of buildings raises new design challenges, in particular with respect to Building Information Modeling and indoor navigation issues since satellite-based localization does not work inside a building. The authors of this paper advocate for using a Systems Engineering approach to address these design challenges. The paper surveys systems engineering methods published in the literature and addresses a method adapted to complex systems which are subject to evolving environment such as drones for indoor inspection of buildings. The chosen method provides several views of the system from a global perspective including mission and operational analyses, to more specific perspectives including requirement, functional, logical and physical analyses.

I. INTRODUCTION

Over the past few years, drones have increasingly been used in transportation, logistics or agriculture to achieve various missions such as firefighting, safety, surveillance, generation of climate data, to name a few [START_REF] Nouacer | Towards a framework of key technologies for drones[END_REF]. In building industry, drones have been used to inspect the facade of high-rise buildings that raise safety risk and danger to Humans [START_REF] Duque | Synthesis of unmanned aerial vehicle applications for infrastructures[END_REF]. In [START_REF] Anwar | Construction monitoring and reporting using drones and unmanned aerial vehicles (UAVs)[END_REF] Anwar, Izhar, Muhammad, and Najam use a drone to monitor a construction site by creating a 3D model at various construction stages. Following decades of using inspection tools based on manual or semi-automated techniques, the building industry can now take many advantages from drones, that are expected to provide time-efficient and cost-effective operations [START_REF] Hasan | Revolutionising building inspection techniques to meet large-scale energy demands: A review of the state-of-the-art[END_REF], relying on sensors that offer precise results.

So far, inspection drones have mostly been used for outdoor inspections of buildings. By contrast, little work has been published on drones for indoor inspections of buildings. Using drones for indoor inspection of buildings raises new design challenges, especially navigation issues since satellite-based localization does not work inside a building [START_REF] Hussein | Key enabling technologies for drones[END_REF], or access to rooms with closed doors. Further, operating within unknown and confined environments such as the interior of an existing building increases the complexity of the mission.

The complexity level reached by indoor inspection drones leads the authors of this paper to propose a Systems Engineering (SE) approach to address the design of these drones. Systems Engineering is an interdisciplinary approach enabling the successful design of engineered systems according to the needs of all the stakeholders involved throughout the system life cycle [START_REF]Systems Engineering Definition[END_REF]. From the definition of customer needs and desires, SE addresses the system from various viewpoints to completely and accurately design the system. The benefits of adopting SE approach are clearly highlighted in the literature [START_REF] Henderson | Value and benefits of modelbased systems engineering (MBSE): Evidence from the literature[END_REF], [START_REF]Systems engineering vision[END_REF]. SE enables an effective collaborative work between the engineering teams [START_REF] Zhang | Application of model-based system engineering method in civil aircraft communication system design[END_REF] as well as managing the increasing complexity of today's systems as technology is rapidly progressing [START_REF] Apvrille | An educational case study of using SysML and TTool for unmanned aerial vehicles design[END_REF]. SE approach addresses the system from the necessary comprehensive perspectives allowing the implication of all involved stakeholders in the decision-making process [START_REF] Arrichiello | Systems engineering and digital twin: a vision for the future of cruise ships design, production and operations[END_REF].

Applying a SE approach requires the adoption of a welldefined method to guide the system engineer through the development of the system. In [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] Brazier, van Langen, Lukosch and Vingerhoeds propose a revised version of the Generic Design Model [START_REF] Frances | Strategic knowledge in design: a compositional approach[END_REF] made up of four design subprocesses: requirements design, artefact design, design process coordination, and simulation, experimentation and implementation. This method encompasses the mission statement, the capture of the needs and the desires of all the involved stakeholders from which requirements will be defined and qualified. The method can be adapted to consider the system from different perspectives, especially operational, functional, logical and physical, to enhance the understanding of the system.

The paper is organized as follows. Section II reviews the current building inspection practices and justifies the need for enhancing indoor inspection systems. Section III surveys related work on Systems Engineering methods. Section IV presents the SE method which is used in this paper. Section V presents its application to the indoor inspection drone. Section VI provides discussion on the results and presents the future work. Section VII concludes the paper.

II. NEED FOR ENHANCING INDOOR INSPECTION SYSTEMS

Inspection activities are carried out during the exploitation life of a building to ensure the safety and security of its users. Inspections concern the building sites themselves as well as monitoring the construction progress. Current techniques are mostly based on hand-held and terrestrial tools [START_REF] Hasan | Revolutionising building inspection techniques to meet large-scale energy demands: A review of the state-of-the-art[END_REF], or even direct visual inspection to detect building defects. Thermal cameras, which detect infrared radiation emitted from object's surface, are used to detect energy-related problems especially heat loss, air leakage and moisture [START_REF] Kirimtat | A review of infrared thermography for the investigation of building envelopes: Advances and prospects[END_REF], as well as cracks and erosions on surfaces [START_REF] Alessia | Advanced damage detection techniques in historical buildings using digital photogrammetry and 3D surface anlysis[END_REF]. The thermographer walks through the building to capture thermal images of the building element or system to be inspected. In [START_REF] Bortolini | Building inspection system for evaluating the technical performance of existing buildings[END_REF] Bortolini and Forcada use mobile techniques combining photos and text or voice annotations to detect and describe the building defects. These techniques rely on smartphone's cameras and mobile applications such as Pick & Go [START_REF] Macarulla | Standardizing housing defects: Classification, validation, and benefits[END_REF], Inspect Anything [START_REF] Snappii | Inspect Anything App[END_REF] and Pocket Inspections [START_REF]Pocket Inspections -Mobile Inspections[END_REF].

These technologies are applicable for both indoor and outdoor inspections, however they are still tedious and time consuming since they require the presence of Humans in field to perform the inspections. Further, they are limited to only items that are accessible to people to minimize safety risks during the inspection.

Building industry can now take many advantages from drones in contrast with traditional inspection tools which are based on manual or semi-automated techniques. Drones are indeed expected to provide time-efficient and cost-effective operations [START_REF] Hasan | Revolutionising building inspection techniques to meet large-scale energy demands: A review of the state-of-the-art[END_REF]. So far, drone-based inspections have mostly been used for outdoor purpose. In [START_REF] Duque | Synthesis of unmanned aerial vehicle applications for infrastructures[END_REF] Duque, Seo and Wacker summarize findings on current visual-based inspection, monitoring, and analysis of infrastructures using drones. Drones are used to detect critical cracks and erosion especially on high-rise buildings facade [START_REF] Chen | Georegistering UAV-captured close-range images to GIS-based spatial model for building fac ¸ade inspections[END_REF], [START_REF] Grosso | Collecting built environment information using UAVs: Time and applicability in building inspection activities[END_REF].

On the other hand, the building industry is processing a digitization of its processes thanks to the implementation of the Building Information Modeling (BIM). According to the European Union BIM Task Group [START_REF]Handbook for the introduction of Building Information Modelling by the European Public Sector[END_REF], 'BIM is a digital form of construction and asset operations. It applies to new build projects; and crucially, BIM supports the renovation, refurbishment and maintenance of the built environmentthe largest share of the sector'. BIM models contain both geometrical and semantic information [START_REF] Temel | Investigation of IFC file format for BIM based automated code compliance checking[END_REF]. Geometric information includes dimensional and spatial data; and semantic information consists in descriptive and procedural characteristics, description of the relationships between the objects of the model, as well as technical properties of the materials (such as mechanical, thermal and acoustic performances). The inspection of buildings can rely on BIM models since they provide information that can be used to manage the maintenance and renovation activities on existing buildings during their exploitation life [START_REF] Gouda | BIM and semantic web-based maintenance information for existing buildings[END_REF].

In this paper, we propose an indoor inspection drone to address the lack of efficient inspection means to evaluate the health of buildings. Drones are chosen since they allow significant time and cost-saving operations, as well as a high flexibility to access hard-to-reach areas which compromise the safety of Humans.

An operator supervises the drone to perform safely and correctly its mission from a ground station which is posted at the foot of the building to be inspected. The drone is intended to collect building information which will be used to inspect both building elements and systems:

• Building elements include structural ones, such as walls, beams, columns, foundations, roofs and floors; as well as non-structural ones which are doors, windows, furniture, stairs, ceiling and floor coverings. • Building systems include electrical, heating, ventilation and air-conditioning (HVAC), plumbing, elevator, and fire protection systems. Disorders in building elements are often visually easy to detect [START_REF] Bortolini | Building inspection system for evaluating the technical performance of existing buildings[END_REF], for instance superficial cracks in structural elements which are usually due to mechanical phenomena such as shrinkage or the differential movement of the bricks and bearing structures composing the element [START_REF] Pereira | Urgency of repair of building elements: Prediction and influencing factors in fac ¸ade renders[END_REF]. In addition, the drone detects surface moisture and vegetation growth which are mostly due to varied phenomena related to damp coming especially from the construction materials used in the building element, or from the rainwater and the condensation.

Regarding the building systems, water leakage on the plumbing system can be identified, which represent the most important percentage of the home insurance claims in France (38% in 2017) [26]. Additional defects related to the safety and the comfort of the users are considered, such as damaged insulation and unsafe sockets on the electrical system, as well as corroded parts and leaks on the HVAC and the fire protection systems. Most of these defects are not visually identifiable since they are concealed in the walls or hidden behind other building elements.

Using drones for indoor inspection raises new design challenges, especially to rely the inspection on BIM, and navigation issues in unknown and confined environment which is continuously evolving with the time. The accessibility in building is by itself a challenge: opening closed doors and windows, move through stairs, avoid furniture, rubbles, etc. Further, such an indoor inspection drone is often embedding multiple, varied payload and therefore massive carrying capacity compromising its energy efficiency and even the flight range. Designing such an inspection drone implies to deal with numerous complex requirements. To master that complexity, one needs to adopt a general approach that encompasses all the aspects of the system throughout all the stages of its life cycle which is the case of the SE approach.

III. SYSTEMS ENGINEERING METHODS

According to the International Council on Systems Engineering (INCOSE), 'Systems Engineering (SE) is a transdisciplinary and integrative approach to enable the successful realization, use and retirement of engineered systems, using systems principles and concepts, and scientific, technological and management methods' [START_REF]Systems Engineering Definition[END_REF]. Applying a SE approach relies on a structured method to guide the system engineer through the development of the system.

In [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] Brazier, van Langen, Lukosch and Vingerhoeds propose a structured systems design and engineering approach made up of four design sub-processes: requirements design; artefact design; design process coordination; and simulation, experimentation, and implementation (see also section IV-B). The method presented in [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] encompasses the mission statement followed by the capture of the stakeholder needs and desires which are translated into requirements. It also includes a logical analysis, and studies the structure of the system, considering its subsystems and their interdependencies. The authors describe their method to be suitable to evolving complex systems such as autonomous vehicles.

In [START_REF] Apvrille | An educational case study of using SysML and TTool for unmanned aerial vehicles design[END_REF] Apvrille, de Saqui-Sannes and Vingerhoeds introduce a method applicable to the development of a broad variety of real-time systems. The method is associated with free SysML [27] tool TTool [28]. It starts with the capture of the modeling assumptions which allow iterative and incremental modeling (starting with a simplified model merging to a more realistic one). Each iteration of the method includes assumptions and requirements capture, analysis, design, simulation and verification. A drone serves as case study to demonstrate the interest of the simulator and the model checker of TTool for early detection of design errors in the life cycle of systems.

In [START_REF] Fei | A model-based system engineering approach for aviation system design by applying SysML modeling[END_REF] Fei, Bin, Rui and Shunhua propose a six-step method covering modeling planing and organization, capture of the stakeholders' needs and definition of the requirements, synthesis of alternative system architectures, integration, and verification and validation (V&V) activities. The method proposed in [START_REF] Fei | A model-based system engineering approach for aviation system design by applying SysML modeling[END_REF] is applied to the Full Authority Digital Engine Control (FADEC) system. The authors claim their method can help designers and engineers to develop complex systems in a consistent manner.

In [START_REF] Li | A functional architecture of prognostics and health management using a systems engineering approach[END_REF] Li, Verhagen and Curran present an architecture design method incorporating the RFLP views: Requirements, Functional, Logical and Physical. These views detail the system requirements, the system functions, the system behaviors and the physical items, respectively. The authors apply the method to build a functional architecture of a Prognostics and Health Management (PHM) system. This method considers system requirements as inputs, and ignores mission definition and stakeholders' needs. It provides several views that enhance the understanding of the system. These views consider the representations of the system from a behavioral and an organizational perspective.

In [START_REF] Hernandez | Model-based systems engineering to design collaborative robotics applications[END_REF] Hernandez and Fernandez-Sanchez present a method for the development of a collaborative robotic system. It has five main steps: identify operational scenarios, specify system capabilities and high level functional requirements, specify quality attributes and system non functional requirements, create system functional and physical architectures. This method facilitates the reuse, the maintenance and the scalability of the robot capabilities and functional interfaces as the robot industry implies constantly new applications presenting new challenges. Compared to previously cited methods, this approach provides an operational view identifying the operational needs and related operational scenarios which can be considered for the design of a drone.

In [START_REF] Sara | Applying model based systems engineering (MBSE) to a standard cubesat[END_REF] Spangelo, Kaslow et al. propose a method for Cube-Sat systems, based on the Object Oriented Systems Engineering Methodology (OOSEM) from the INCOSE. OOSEM is a system-level development method combining object-oriented concepts with systems engineering practices. It starts with the mission definition, and includes requirements analysis, logical and allocated architecture design, guidance for trade studies, verification and validation activities [START_REF]Object-Oriented SE Method[END_REF].

These SE methods highlight complementary concepts needed to describe systems, including mission, needs and requirements, operations, functions, behaviors, structures and architectures. These concepts allow one to consider the system throughout all the stages of its life cycle, which improves the understanding and the quality of the results. Further, the use of models to support the development can be considered, that is applying Model-Based Systems Engineering (MBSE). MBSE is expected to improve system design in terms of completeness, process, integrity, accuracy, maturity, performance, and clarity [START_REF] Henderson | Value and benefits of modelbased systems engineering (MBSE): Evidence from the literature[END_REF]. The next section presents the SE method that is used in this paper to address an indoor inspection drone.

IV. SYSTEMS ENGINEERING METHOD FOR AN INDOOR INSPECTION DRONE

The systems design and engineering approach presented in [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] is suitable for complex systems which are subject to continuously evolving environment such as indoor inspection drones.

A. Qualitative criteria

Relying on the SE methods presented in the previous section, we may enumerate qualitative criteria in order to adapt the method presented in [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] which is used in this paper to manage the complexity of the development of an inspection drone. The method shall:

• C1: start with the mission definition, because the mission makes one drone different from others. C1 includes the capture of the customer's initial needs and the formalisation of the problem expected to be solved by the system [START_REF] Faisandier | Systems opportunities and requirements[END_REF]. • C2: address the needs from all the involved stakeholders at each stage of the system's life cycle. C2 allows to consider the entire system from its genesis to its end [START_REF] Bron | System Requirements Engineering: A SysML Supported Requirements Engineering Method[END_REF]. • C3: consider the requirement analysis. C3 includes the definition of stakeholders' requirements which are derived from their needs and will be translated into system requirements [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF]. • C4: support the traceability of the requirements along the development process. The requirements can be refined as much as we go further into the design of the system (iteration) [START_REF] Apvrille | An educational case study of using SysML and TTool for unmanned aerial vehicles design[END_REF]. • C5: address an operational view of the system. C5 allows to consider the behavior of the system within its opera-tional environment including all the possible operational concepts and modes [START_REF] Bron | System Requirements Engineering: A SysML Supported Requirements Engineering Method[END_REF]. • C6: address a functional view of the system. C6 allows to consider the functions, that are tasks, actions or activities to be performed by the system to achieve its intended mission [START_REF] David D Walden | Systems engineering handbook: A guide for system life cycle processes and activities[END_REF]. • C7: address a logical view of the system. C7 allows to study the behavior of the system to propose a preliminary system architecture [START_REF] David D Walden | Systems engineering handbook: A guide for system life cycle processes and activities[END_REF]. • C8: address a physical view of the system. C8 allows to design a physical architecture of the system containing the physical components and subsystems [START_REF] David D Walden | Systems engineering handbook: A guide for system life cycle processes and activities[END_REF].

On the basis of these criteria, the method presented in [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] can be adapted to encompass comprehensive and complementary viewpoints for a complete and an accurate design of an inspection drone.

B. Adopted Systems Engineering method

The method adopted in this paper is based on the systems design and engineering approach presented in [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF] which is made up of four sub-processes. Requirements design translates the stakeholder needs and desires into requirements. Artefact design describes a design object according to the requirements provided by the requirements design. Design process coordination focuses on the determination of the design strategy, including the progress of the design process and the allocation of resources. Simulation, experimentation and implementation allow one to assess the feasibility of the defined requirements and/or design, and facilitate the exchanges with the stakeholders. These sub-processes are active in parallel for each of the stages of the system's life cycle to cope with the continuous evolution of the system and its environment. Hence existing requirements have to be adapted and new ones can be specified since stakeholder needs and desires evolve during the system's life cycle. This is an iterative, seven-step method which is satisfying the set of criteria presented in previous subsection.

1) Mission Analysis: defines the purpose, the mission and the objectives of the system to be developed to solve the problem which is identified from the initial needs expressed by the customer. The purpose justifies the existence of the system whereas the mission expresses the set of services it provides [START_REF] Faisandier | Systems opportunities and requirements[END_REF]. Initial objectives can be formulated to quantify the mission. These characteristics can be refined as much as needed during all the design process.

2) Requirement Analysis: includes the capture of the needs of stakeholders and their translation into system requirements. Stakeholder needs and desires are captured to define correctly their expectations on the system. The stakeholders can be classified into 4 categories according to the typology proposed by the 'Association Franc ¸aise de l'Ingénierie Système' (AFIS), the French association of systems engineering in [START_REF] Gdtis Afis | Découvrir et comprendre l'ingénierie système[END_REF]:

• The acquirer side who will be concerned by the use of the system,

• The supplier side who will be implied in the realization of the system, • Potentially concerned parties who will be more or less directly impacted by the system, • Regulatory bodies, such as certification and standardization agencies.

These needs and desires of stakeholders are interpreted and translated into stakeholder requirements and then into system requirements afterwards [START_REF] Sotelo | Quality Assurance methodology for system requirement definition[END_REF]. Requirements can be qualified into functional, structural, behavioral and experiential [START_REF] Brazier | Complex systems: design, engineering, governance[END_REF].

• Functional requirements consist in the functions provided by the system. • Behavioral requirements consist in the desired behavior of the system. • Structural requirements concern the components of the system and their relationships. • Experiential requirements concern the effect of the use of the system in the real world.

3) Operational Analysis: aims to study the behavior of the system within its operational environment to deliver its intended services. Operational Analysis addresses the operational concepts describing the characteristics of the system from the operator's viewpoint to identify the system context and interfaces [START_REF] David D Walden | Systems engineering handbook: A guide for system life cycle processes and activities[END_REF]. The different operational modes are studied to consider all the possible operational situations of the system throughout all its operational life. 4) Functional Analysis: identifies the functions and subfunctions of the system to satisfy the previously defined system requirements. The aim is to identify the functions and services that the system has to provide, independently of any technical solutions or technological choices, and their interfaces to build a functional architecture of the system. 5) Logical Architecture Analysis: translates the functional architecture into logical architecture. Logical Analysis focuses on decomposing the system into logical elements to design a preliminary system architecture [START_REF] David D Walden | Systems engineering handbook: A guide for system life cycle processes and activities[END_REF]. The functions are allocated to logical elements which can be some available technical solutions, subsystems or components, but being identified at high level.

6) Physical Architecture Analysis: focuses on designing the physical architecture of the system from the previously defined functional and logical architectures. It is composed of physical components with a given set of properties, and their physical interfaces to which the functions are allocated to describe how will the system be implemented.

7)

Verification & Validation: are transverse activities to every life cycle stage of the system. They can be applied to every engineering element that has contributed to the definition of the system itself (operations, requirements, functions, data flow, system elements, interfaces, and design properties). The purpose of V erif ication is to provide evidence that the 'product is built right', while V alidation is intended to ensure that the 'right product is built'. [START_REF] David D Walden | Systems engineering handbook: A guide for system life cycle processes and activities[END_REF].

V. APPLICATION TO AN INDOOR INSPECTION DRONE

In this section, the SE method presented in the previous section is applied at a high level of abstraction to the indoor inspection system of an existing building.

A. Mission Analysis

From the stated problem that the indoor inspection system is intended to solve as presented in section II, we can define its purpose, its mission, and initial objectives to have a vision on the way forward.

• Problem: There is a lack of efficient inspection systems for the interior of an existing building today which makes the operations tedious and costful. • Purpose: The system provides a cost efficient and time saving solution to inspect the interior of an existing building.

• Mission: The system collects valuable building information from the interior of an existing building in order to detect disorders both on building elements and systems. • Objectives: Initial objectives are as follows:

1) The system detects disorders of the building especially cracks more than 1mm of large; water leakage causing surface moisture, vegetation growth and corrosion. 2) The system reduces the inspection time and the workload by at least 30% compared to the conventional means.

B. Requirement analysis

Involved stakeholders during all the stages of the system life cycle are as follows:

• Acquirer side: Building owner, Building inspection specialist, Civil engineer, Operator. Civile (DGAC), the French National Aviation Authority. Stakeholders can express needs and desires related to the optimization of indoor inspections relying on a cost effective and time saving system. The customer expects the system to provide data which are used to detect the disorders of the building. These defects include cracks more than 1mm of large, water leakage, surface moisture, vegetation growth and corrosion both on building elements and systems. As stated in section II, an operator supervises the system from a ground station while it is performing the inspection to ensure safe and secure operation. Up to two people are required to manage the whole system within its operational environment which can reduce the workload by more than 30% compared to the traditional means.

Requirements are derived from the needs and desires in close interaction with stakeholders. A glossary can be built to

Structural

Stk.S.Rq.01 The system shall include the building information collection system Stk.S.Rq.02 The system shall include the information storage system Stk.S.Rq.03 The system shall include the obstacle avoidance system Stk.S.Rq.04 The system shall include the energy storage system Stk.S.Rq.05 The system shall include the communication system define the words precisely and avoid ambiguity in the requirements statements. Glossary terms are capitalized throughout the requirements to show that they have a specific meaning in the context of the statements. An extract of these high-level stakeholder requirements are listed in table I and are classified into functional, structural, and behavioral requirements. These stakeholder requirements will be translated into system requirements that will serve as inputs for the next steps of the development process. For instance, from the Stk.F.Rq.02 'The system shall detect the Disorders of Building Elements', several system requirements can be derived such as 'The system shall take pictures of the Building Element', 'The system shall analyze the picture of the Building Element', and 'The system shall store pictures of the Building Element'.

C. Operational analysis

Operational analysis considers the system within its operational environment. In addition to the involved stakeholders listed in the previous subsection, external systems interacting with the system during its operational life can be identified.

OM03 Initialization

The system is performing self-tests and calibrations before starting the operation

OM04 Nominal

The system is operated with all its nominal capabilities OM05 Stand-by The system is waiting for mission order or for powering off order OM06 Failure

The system is encountering failure impacting one or more capabilities

OM07 Training

The system is operated with limited capabilities which are sufficient for training purpose

For instance, the system interacts with the building itself which is composed of its elements and systems as stated in section II; the external energy source which provides energy to the system; and the environmental conditions of the operational site which is characterized by the weather, the lighting, the visibility and the wind speed. One can describe operational concepts of the indoor inspection system such as:

• OC01: Move inside the building. The system accesses and moves inside the building to cover the targeted areas where it interacts with the building environment. • OC02: Detect disorders of the building elements and systems. The system performs its intended mission consisting in the inspection of the interior of the building. • OC03: Ensure the security and the safety of the operator and the residents. While performing its mission, the system interacts with the operator and potentially with the residents of the building to be inspected with respect to all the security measures.

From these operational concepts, system use cases can be defined. For instance, from 'Detect disorders of the building elements and systems', use cases are: 'Take picture', 'Analyze picture', 'Identify disorders', and 'Store picture'. Operational modes of the indoor inspection system are presented in table II.

D. Functional analysis

The indoor inspection system can be decomposed into relevant high-level functions to meet the previously defined requirements:

• F01: Provide power. The system provides power to all of its components to achieve the mission. • F02: Manage power. The system manages the power with a twofold objective: to move to and within the building, and to supply power to its internal subsystems. • F03: Control the dynamics. The system controls its dynamic characteristics, that are its speed and its direction, to allow the movement and to perform the intended mission. • F04: Communicate with the operator. The system exchanges information with the operator to receive orders and provide him or her with constant feedback. • F05: Track the localisation of the system. The system is able to track its position within the building to be inspected.

• F06: Collect building information. The system collects a significant amount of information of both building elements and systems. • F07: Analyze information. The system should process the gathered building information to detect disorders on both elements and systems. • F08: Store building information. The system stores the gathered building information which will be provided to the operator. • F09: Adapt to the environmental conditions. The system considers the environmental conditions of the operation site, especially the ability to operate in poorly lit sites and to deal with airstreams. From these high-level functions, the context diagram in figure 1 depicts a high-level functional architecture of the system. The central rectangle represents the indoor inspection system as a set of functions. It is composed of high-level functions represented by the rectangles with a stereotype 'function' for clarification, and linked to the indoor inspection system thanks to composition associations (represented by solid lines with solid diamonds on the composite end). These high-level functions can be refined and decomposed into low-level functions to specify the subsystems composing the system.

E. Logical analysis

The high-level functions presented in the previous section can be associated to logical blocks to build a preliminary logical architecture of the indoor inspection system. It is composed of:

• L01: Power Supply. This block stores energy from the external source and provides energy to the whole system.

• L02: Controller. This block has two main objectives: control the dynamics of the system based on the operator input and sensor environmental data, and to manage the others subsystems work together to perform the intended mission. • L03: Sensor. This block has two main objectives: measure the dynamic characteristics of the system to compare them to the desired states, and to sense the building and its environment to collect the requested building information. • L04: Actuator. This block converts energy into action which is the motion of the system corresponding to controller input. • L05: Structure. This block acts as the skeleton of the system. It holds the other blocks composing the whole system. • L06: Communication. This block is responsible of the exchanges between the system and the operator. • L07: Data Storage. This block is responsible of storing the gathered building information. • L08: Ground Station. This block represents the station from where the operator supervises the mobile part while it is collecting the building information. The preliminary logical architecture of the indoor inspection system is presented in the SysML context diagram in figure 2. The inspection system is represented by the central rectangle to which the logical blocks are linked with composition associations. At this abstraction level, the logical architecture decomposes the system into subsystems represented by each single logical block.

F. Physical analysis

The physical analysis consists in converting the logical architecture into physical architecture which is composed of physical components of the system and their interfaces. To get there, the previous steps of the SE method shall be applied at lower levels of abstraction to detail the subsystems. Further, trade-offs and comparison of different alternative solutions are performed through performance evaluation, sizing, dynamic and structural simulations, and similar.

For instance, imaging sensors such as digital cameras can be selected to get accurate pictures of the building elements to be able to detect the disorders in building elements such as cracks with 1mm of large or surface moisture and corrosion as stated in the requirements. Trade-offs are performed to choose the most convenient reference satisfying the required picture resolution and the weight and power consumption constraints.

Microwave sensors [START_REF] Wavescan | WaveScan -Electromagnetics based sensor solutions for smarter liveable cities[END_REF] can be selected to detect hidden wiring and pipelines concealed behind walls and flooring to identify anomalies such as damaged insulation and leakages on the electric, HVAC, and plumbing systems of the building. For the propulsion system, electric motor with a given characteristics, such as the torque and the supply voltage, can be selected after having sized the necessary force to move the mobile part.

Given the association of sensors, motors and all the components of the system, the power supply system can be specified. For instance Li-Po battery with a given characteristics including the capacity and the voltage are often used for drones depending on the required power.

Verification & Validation activities check the consistency between each step of the method including the traceability of the requirements throughout the design of the functional, logical and physical architectures.

VI. DISCUSSION AND FUTURE WORK

The previous section applies the systems engineering method at a high level of abstraction, which consists in the system level of an indoor inspection drone whose mission is to detect disorders of an existing building. Further, it can be iteratively applied at lower-levels corresponding to each subsystems of the system, until resulting in a system completely and accurately characterized by all of its physical components satisfying the needs of the customer.

This work shows how to master the complexity for the development of such a drone using a SE approach by addressing the system from various comprehensive and complementary perspectives. Besides, the use of SysML models enhances the understanding of the system.

The case study used in this paper is a part of the 'Indoor Multi-Usages Drone Acquisition' (IMUDA) project whose objective is to prototype an indoor multi-usages drone. The mission of the indoor inspection drone can be extended to other usages, such as the evaluation of the technical performances of the building, that can include structural, architectural, energetic, thermal, hygrometric, and acoustic performances. Further, the drone can be used to create BIM models of the building in the frame of the digitization of the building sector, in order to address and support those various usages, especially the inspection one as stated in section II. Studies are currently underway to refine the mission of the multi-usages drone, and to apply the SE method presented in section IV to master the complexity of the multi-usages drone development. Each step of the method must take the BIM constraints and requirements into account.

VII. CONCLUSIONS Introducing an indoor inspection drone is an asset to the building industry since it allows time and cost saving operations compared to the conventional means. Designing such a drone is a highly complex task. To cope with that complexity, this paper advocates for using a Systems Engineering approach. This paper compares several SE methods published in the literature according to a set of qualitative criteria to specify a drone.

An iterative seven-step method is proposed to address an indoor inspection drone. The SE method provides comprehensive and complementary views of the system from a global to very specific perspectives which allow mastering the complexity for the inspection drone development. Further, the use of models to run the proposed method enhances the understanding and the design process of the drone throughout the system life cycle, that is applying a Model-Based Systems Engineering approach. The mission of the inspection drone can be extended to several capabilities, especially the creation of digital models of the building to address multiple usages in the frame of the digitization of the building sector.

Fig. 1 .

 1 Fig. 1. High-level functional architecture of the indoor inspection system

Fig. 2 .

 2 Fig. 2. Preliminary logical architecture of the indoor inspection system

  Rq.01 The system shall operate inside the Building Stk.F.Rq.02 The system shall detect the Disorders of Building Elements Stk.F.Rq.03 The system shall detect the Disorders of Building Systems Stk.F.Rq.04 The system shall evaluate the Technical Performances of the Building Stk.F.Rq.05 The system shall provide feedback to the operator Stk.F.Rq.06 The system shall receive Order from the operator Behavioral Stk.B.Rq.01 The system shall operate within a radius of 100m from the Ground Station Stk.B.Rq.02 The system shall work in Building with Multiple Storeys Stk.B.Rq.03 The system shall ensure the safety of the operator Stk.B.Rq.04 The system shall ensure the safety of the resident Stk.B.Rq.05 The system shall handle abnormal situations Stk.B.Rq.06 The system shall have an energy autonomy more than 45min Stk.B.Rq.07 Up to two people shall handle the system Stk.B.Rq.08 The system shall emit less than 80dB noise Stk.B.Rq.09 The system should be set up within less than 15min Stk.B.Rq.10 The energy storage system shall be accessible within less than one minute

		TABLE I
		STAKEHOLDERS REQUIREMENTS
	ID	Stakeholders requirements
		Functional
	Stk.F.


B M Retailleau

M Benali

G Quéméner

P Delahaye

X Fléchard

E Liénard

Geometry optimisation of a transparent axisymmetric ion trap for the MORA project

Keywords: . Paul trap, Spherical harmonics, Laplace's solver, Boundary element method

Geometry optimisation of a transparent axisymmetric ion trap for the MORA project

Introduction

Precision measurements in nuclear β-decay provide a remarkable tool to improve the accuracy of Standard Model (SM) parameters and to search for New Physics (NP) beyond, at the low energy frontier [START_REF] Gonzáles-Alonzo | New physics searches in nuclear and neutron β decay[END_REF]. In particular, the search for new sources of CP (Charge Parity) violation is one of the requirements to explain the matter-antimatter asymmetry observed in the universe, according to Sakharov's criteria [START_REF] Sakharov | Violation of CP invariance, C asymmetry, and baryon asymmetry of the universe[END_REF]. This search can be achieved in nuclear β-decay by measuring the triple correlation between the parent nucleus spin (J ), electron momentum (p e ), and neutrino momentum (p ν ): D J • (p e × p ν ). This triple correlation is sensitive to T (Time) reversal violation and thus to CP violation thanks to CPT conservation. Such a violation would be quantified by a non-zero D value, experimentally determined from an asymmetry in the βrecoil angular distribution measured in a plane perpendicular to J for two opposite directions of this nucleus orientation. In this context, the new project MORA [START_REF] Delahaye | The MORA project[END_REF] aims to measure the triple-correlation D coefficient in the β-decay of laser-polarised 23 Mg + and 39 Ca + ions confined in a Paul trap applying a quadrupole radio frequency (RF) field. The use of Paul traps is considered as an innovative technique in precision measurements of correlation coeffia e-mail: benali@lpccaen.in2p3.fr b e-mail: quemener@lpccaen.in2p3.fr cients in nuclear β-decay in the SM framework [START_REF] Sternberg | Limit on Tensor Currents from 8 Li β Decay[END_REF][START_REF] Brodeur | V ud determination from light nuclide mirror transitions[END_REF]. One example is the transparent Paul trap lpctrap used in the measurement of the β -ν correlation coefficient, a βν , in the decay of different nuclei [START_REF] Ban | Precision measurements in nuclear β-decay with LPCTrap[END_REF][START_REF] Fabian | Precise measurement of the angular correlation parameter a βν in the β decay of 35 Ar with LPCTrap[END_REF][START_REF] Liénard | Precision measurements with LPCTrap at GANIL[END_REF][START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF].

The central element of the MORA apparatus is a transparent Paul trap (moratrap), which will be used to confine singly charged radioactive ions, coupled to a laser system allowing to polarise the nucleus by optical pumping. As shown in Fig. 1, moratrap, which will be installed in a vacuum chamber, is surrounded by four pairs of electron and recoil ion detectors arranged alternately in an octagonal geometry in the azimuthal plane of the trap and allowing close to 2π azimuthal coverage. Two annular silicon detectors (not visible on the figure) located on the trap axis will monitor the polarisation degree thanks to β asymmetry measurement.

The RF potential generated in the trapping volume or region of interest (ROI) is not perfectly a quadrupole but contains some small amplitude, higher order electric multipole components which disturb the ion's motion. Since the potential in the ROI depends on the electrodes shape and on the applied voltages, an optimisation of the trap geometry is mandatory to reduce the higher order harmonics and to generate an optimised quadrupole potential.

In this work, we will describe an efficient method to optimise the moratrap geometry, using two electrostatic solvers developed at LPC Caen. Our trap is a three dimen- sional transparent Paul trap, inspired from the lpctrap geometry, where the electrodes have been modified from an ideal Paul trap to allow the detection of β-decay products in a large solid angle around the azimuthal plane and obtain an efficient injection of the ion bunches. By optimising the geometry, we aim to reach a quadrupole potential of higher quality, to minimise the ion losses from the trap and to increase the trapping lifetime and the space charge capacity. Such improvements are mandatory to reach the statistics required in the MORA experiment to search for NP.

The potential of an ideal Paul trap is commonly defined in cylindrical (r, z) coordinates as:

V ideal (r, z, t) = V 0 cos(Ωt) r 2 -2 z 2 2 r 2 0 , (1)
where V 0 is the amplitude of the RF voltage, Ω = 2πf its pulsation, and r 0 a distance parameter. 3D Paul traps usually approximate the ideal Paul trap by using electrodes consisting of one ring and two end caps, whose surfaces form truncated hyperboloids (see for example Fig. 1 in [START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF]). In such a configuration, r 0 is the distance of the ring surface to the centre of the trap. In geometries which significantly depart from the ideal trap, one can define an effective trapping radius r eff which contains the region where the quadrupole field still sufficiently dominates, so that ion trajectories are still stable in standard conditions defined for the ideal potential. This region is defined in [START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF] as the region for which the contribution to the potential from harmonics of higher order than the quadrupole one, are still below a few percents. The conditions for stability of an ion of mass m and charge q in an ideal Paul trap, in the absence of DC field, are defined with respect to the Mathieu parameter q z :

q z = 4 q V 0 m r 2 0 Ω 2 .
(

) 2
In the first stability region, |q z |< 0.908 [START_REF] March | Quadrupole Ion Trap Mass Spectrometry, Chemical Analysis, a series of monographs on analytical chemistry and its applications[END_REF]. In the pseudopotential approximation limit, valid for small q z values, one can define a pseudo-potential depth D z (resp. D r ) for the z (resp. r) dimension:

D z = q z V 0 8 = 2 D r . (3)
In this model, the maximal charge density ρ max the Paul trap can hold is

ρ max = 3 ε 0 D z z 2 0 , (4)
where ε 0 is the vacuum permittivity and z 0 = r 0 / √ 2 is the distance of the end caps from the trap centre. From Eqs. ( 2) to ( 4) and considering the fact that experimentally one usually fixes a Mathieu parameter in the middle of the stability diagram q z 0.4 [START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF], the maximal charge capacity Q max a trap can hold directly relates to the product r 0 V 0 via the formula:

Q max = 1 2 πε 0 q z r 0 V 0 (5)
in the case of an ideal trap, and

Q max = 1 2 πε 0 q z r eff V eff (6)
in the case of a Paul trap with a limited quadrupole region of radius r eff , where

V eff = V 0 r 2 eff r 2 0
is the maximum potential at the radius r eff . In order to optimise the charge capacity of the trap for the MORA experiment, the optimisation procedure presented in this article therefore aims at enlarging the trapping region, of radius r eff , given some space constraints on the size of the trap.

Sec. 2 presents the electrostatic solvers necessary to compute the potential from an electrode set and Sec. 3 introduces the harmonics series expansion describing this potential, solution of Laplace's equation. Sec. 4 details the objective function used in the optimisation procedure whose results are presented and commented in Secs. 5 and 6. Finally, the effects of electrode misalignment on the potential quality are investigated in Sec. 7.

Laplace's solver

The shape of the electric potential generated within a trap depends on both the electrode geometry and their corresponding applied voltages. In order to estimate and optimise this potential, it needs to be accurately computed within the ROI (the trapping region). This is achieved by solving Laplace's equation in this ROI. For that purpose, the most common approaches use either the finite element method (FEM) or the finite difference method (FDM) which require to mesh both the electrodes and the free space volume. For our study, we have used a home made C++ Laplace's solver based on the boundary element method (BEM). This solver has been developed by one of us and thoroughly validated with analytic examples and expensive commercial software like simion (FDM) or comsol (FEM). Compare to FEM and FDM, BEM exhibits several advantages: it only requires to mesh the surface of the electrodes and is thus able to achieve higher precision with less computation time and memory requirements. This is especially true for geometries involving electrodes with large aspect ratios (e.g. a very thin thickness for very large length and width) for which the 3D meshing quality would require special care. In addition, BEM can easily deal with open systems as the boundary conditions are inherent to the formalism.

Even if our BEM program also handles dielectric materials, in order to simplify its description, here we shall only concentrate on sets of electrodes with applied voltages. A setup is described by N e electrodes each represented by its surface. The surface of electrode e (e ∈ [1, N e ]) is meshed in n e flat polygonal cells (triangles, quadrangles, ...). The full setup is therefore represented by the set of cells C = {C i : i ∈ [1, N ]} where N = Ne e=1 n e is the total number of cells, each with an a priori unknown associated surface charge density σ i assumed constant over the cell surface. The n e cells of electrode e obviously share the same potential applied on the whole electrode. The potential V i of cell C i centred at position r i is related to the surface charge densities σ j of cells C j (j ∈ [1, N ]) through the superposition principle as:

V i = N j=0 σ j 4πε 0 Cj 1 r i -r d 2 r = N j=0 Q ij σ j , (7)
where the integral over the surface of cell C j only depends on C j shape and on the relative position of the target point r i with respect to C j location, as σ j is assumed constant over the whole surface of cell C j . An analytic formula has been derived for the integral in Eq. ( 7). This formula is too complicated to be discussed here and shall be published in a separate article, but it should be emphasised that special care has been taken to suppress numerical instabilities in its evaluation especially when i = j where numerical divergence may occur. Eq. ( 7) can be written for each cell C i , leading to the following set of N equations:

    V 1 V 2 . . . V N     =     Q 1,1 Q 1,2 • • • Q 1,N Q 2,1 Q 2,2 • • • Q 2,N . . . . . . . . . . . . Q N,1 Q N,2 • • • Q N,N         σ 1 σ 2 . . . σ N     , (8)
where V i are the known potentials applied on electrodes and σ i , i ∈ [1, N ] are the unknown surface charge densities. The matrix elements Q ij are computed for a given electrode assembly/geometry leading to a dense square matrix on the contrary to FEM and FDM which deal with sparse matrices of much larger dimensions. The presence of a dense matrix requires special algorithms to solve Eq. ( 8) for the σ i : a direct solver such as LU decomposition can be used to obtain an exact solution or an approximated solution can also be computed more rapidly with iterative solvers such as CMRH [START_REF] Sadok | CMRH: A new method for solving nonsymmetric linear systems based on the Hessenberg reduction algorithm[END_REF]. Once the charge densities have been determined, the potential and field components can be evaluated at any location r in space without interpolation contrary to the FEM and FDM approaches:

V (r) = N i=0 σ i 4πε 0 Ci 1 r -r d 2 r (9)
E(r) = N i=0 σ i 4πε 0 Ci r -r r -r 3 d 2 r (10)
As for Eq. ( 7), integrals in Eqs. ( 9) and ( 10) are analytically evaluated over the surface of each cell.

For our studies we have used two versions of the solver: a full 3D solver, called electrobem, and a derived version, axielectrobem for axially symmetric problems 1 . In electrobem, a complex 3D setup can be modelled with a thorough set of geometric functions handling predefined shapes (flat polygons, ring, cylinder, cone, torus, ...), rotations, translations and symmetry planes as well as fine tuning of mesh properties. This last point is crucial for the precision of the solution. The software can also import pre-computed meshes from the open source, reliable and user friendly gmsh software [START_REF] Geuzaine | Gmsh: a three-dimensional finite element mesh generator with built-in pre-and postprocessing facilities[END_REF]. Both gmsh and root [START_REF] Brun | root -an object oriented data analysis framework[END_REF] may be used to display the geometry. Electric potentials are applied to electrodes and then, once the electrostatic problem has been solved, a complete set of plotting/mapping functions allows to display and/or export the potential and field components in 1D, 2D and 3D as illustrated for instance in Figs. 2 and3. In axielectrobem, any shape axially symmetric around z-axis can be simulated by predefined shapes or by functions (r = f (z)). The meshing is done along longitudinal z-axis and radial r-axis. In this version, the number of cells necessary to model a given axisymmetric setup is substantially smaller than it would be in the 3D version for the same geometry, resulting in a much smaller computation time. axielectrobem is thus particularly well suited for the optimisation of the axially symmetric moratrap and was used in a first phase of our study, whereas the 3D version served in a second phase to estimate the effects, on the trapping potential quality, of possible mechanical misalignment or machining precision which break the setup axisymmetry. Depending on the version, 2D or 3D, the potential inside the trap is expanded differently in terms of multipole coefficients, as shown in next section and in appendix A.

Laplace's equation: spherical harmonics series expansion

In a source-free region, i.e. in the absence of charges, the electric potential satisfies Laplace's equation ∆V (ρ, θ, ϕ) = 0 with (ρ, θ, ϕ) being the usual spherical coordinates. The solution of this equation [START_REF] Durand | Électrostatique et Magnétostatique[END_REF], may be expressed in terms of a uniformly convergent Laplace series also known as spherical harmonics series expansion in the following way:

V (ρ, θ, ϕ) = ∞ n=0 n m=0 2n + 1 4π (n -m)! (n + m)! ρ R 0 n P m n (cos θ) × [α nm cos(mϕ) -β nm sin(mϕ)] = ∞ n=0 n m=0 ρ R 0 n P m n (cos θ) × [A nm cos(mϕ) -B nm sin(mϕ)] , (11)
where P m n (cos θ) are the associated Legendre functions of the first kind, of degree n and order m, and A nm (resp. B nm ) are the normal (resp. skew) spherical harmonics co-efficients (α nm and β nm being the corresponding unnormalised harmonics). R 0 is the convergence radius of the above series expansion, it must be smaller than the distance between the origin (ρ = 0) and the closest electrode: the region of convergence or the ROI must remain a source free region.

In axisymmetric cases with azimuthal symmetry, the potential does not depend on the ϕ angle and therefore we only keep terms with m = 0 in the above expansion. Dropping the m index and the ϕ dependence, Eq. ( 11) reduces to:

V (ρ, θ) = ∞ n=0 α n 2n + 1 4π ρ R 0 n P n (cos θ) = ∞ n=0 A n ρ R 0 n P n (cos θ) ( 12
)
with P n (cos θ) the Legendre polynomial of degree n and A n the harmonic coefficient of order n (α n is the corresponding unnormalised harmonics). A 0 is a constant term for the potential and therefore it does not appear in the field components which enter the ion equations of motion.

A 1 is called the dipole harmonics, A 2 the quadrupole one and so on for higher orders. The harmonics A 2k+1 with k = 0, 1, 2, . . . are said to belong to the dipole series whereas harmonics A 4k+2 with k = 0, 1, 2, . . . are said to belong to the quadrupole series. Because of symmetry reasons, a potential exhibiting a quadrupole behaviour, is more likely to have harmonics from the quadrupole series (A 2 , A 6 , A 10 , . . . ) dominant over other higher order terms. One should also notice that the general trend of A n is to decrease with increasing n and that the contribution from A n to the potential is weak close to the origin and increases as ρ n . If, in addition to the axial symmetry, a setup exhibits a symmetry plane at z = 0, then only harmonics A n with even order n will be present in the expansion. The harmonic coefficient A n is obtained by integration over the sphere S of radius R 0 by:

A n = 2n + 1 4π 2π 0 dϕ π 0 P n (cos θ) V (R 0 , θ) sin θ dθ.
(13) After integration over ϕ, and assuming a symmetry plane at z = 0, A n becomes:

A n = (2n + 1) π 2 0 P n (cos θ) V (R 0 , θ) sin θ dθ. (14)
In order to determine the A n coefficients for a given electrode configuration, Eq. ( 14) is numerically integrated using a 64-nodes Gauss-Legendre quadrature with the potential V (R 0 , θ) calculated on the circle of radius R 0 at the Gauss-Legendre nodes (θ i ) using the BEM solver described in previous section. The extracted harmonic spectrum {A n } will then be used in the minimisation procedure of an objective function to optimise the trap geometry. A similar approach, combined with FFT (Fast Fourier Transform), is used in 3D to determine the harmonics A nm and B nm , but a sphere of radius R 0 is used instead of a circle to sample the potential V (R 0 , θ, ϕ). The potential is determined on a θ -ϕ grid of 180 × 360 = 64800 points.

Let us now introduce the objective function used in this study.

Objective function

As already mentioned, an ideal Paul trap requires a pure quadrupole potential, i.e. A 2 as large as possible while A n = 0 ∀n > 2. In practice, this is almost impossible to achieve and higher order harmonics are present in the potential. The present study aims at optimising the geometry of the moratrap electrodes by minimising higher order harmonics in order to avoid the loss of trapped ions because of trajectory instabilities induced by these higher order terms. As shown in [START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF], a necessary condition to suppress ions losses is to keep the relative contribution to the potential from higher harmonics to less than a few %, say 2%, in the ROI (trapping region). Neglecting the constant term A 0 which does not contribute to the field and only considering n-even terms because of planar symmetry at z = 0, this condition translates into a relative difference to the quadrupole A 2 as:

nmax n=2 ρ R0 n |A n |-ρ R0 2 |A 2 | ρ R0 2 |A 2 | = nmax n>2 ρ R 0 n-2 A n A 2 0.02, (15)
where the series expansion is truncated to n max whose value will be discussed in next section. Choosing 2% as an upper limit is therefore equivalent to find a root ρ 2% of a polynomial of degree n max -2 in ρ. Optimising the trap performances amounts to maximise both the A 2 term contribution as well as ρ 2% . Both operations can be achieved by minimising the following objective or fitness function:

f (a) = 1 ρ 2% (a) R0 2 A 2 (a) , (16)
where ρ 2% (a) and A 2 (a) explicitly show their dependence on the geometry via the vector a which contains the different free parameters describing the electrodes. To keep the shape of the electrodes as simple as possible and therefore to facilitate their machining at LPC Caen, we have chosen the following parameters for some conical electrodes as illustrated in Fig. 4:

-Z min , the minimal axial distance from the trap r-axis to the electrode. -R min , the innermost radius which is the minimal radial distance from the trap z-axis to the electrode. -the cone angle θ.

the electrode section thickness T h.

Some other parameters can be derived from these quantities like the gap G between two electrodes, their outermost radii as well as the radial and axial angular acceptances respectively Ω R and Ω A . These could help imposing some mechanical constrains in the minimisation procedure. For some exploratory work, we had envisaged and tested more complicated shapes for the electrodes, including e.g. some splines to describe parts of their section, but the whole study showed that simpler shapes could suit our precision needs as we shall see in next section.

Optimisation process and results

Maximising the measured statistics with moratrap requires not only the trapping of as many ions as possible during a period as long as possible, but also a large angular aperture Ω R in the radial direction where the β and recoil ions detectors are located. In addition, the D measurement relies on the knowledge of the ion cloud polarisation which shall be continuously monitored during the experiment [START_REF] Delahaye | The MORA project[END_REF]. This will be performed with two annular silicon detectors located along the trap z-axis before and after the electrodes. The central hole in these Si detectors lets the ion beam entering and exiting the trap. The axial angular acceptance of these detectors directly impacts the precision on the polarisation measurement and is constrained by the trap angular opening Ω A . All the above constrains have been taken care of in the optimisation process.

A simplified version of lpctrap (Fig. 5) was used as a starting point to look for an optimal geometry2 as it already provides a large enough radial acceptance Ω R [START_REF] Fabian | Precise measurement of the angular correlation parameter a βν in the β decay of 35 Ar with LPCTrap[END_REF]. This starting geometry respects axial symmetry and also possesses a planar symmetry at z = 0. However, on the contrary to lpctrap which uses three pairs of tube or cylinder shaped electrodes, in order to satisfy moratrap's requirement of a larger axial acceptance Ω A , conical shapes were chosen for the two innermost pairs of electrodes. A third pair of electrodes (R 5 -R 6 ) was added to screen the axially symmetric trap region from the octagonal set of detectors and associated collimators which might disturb the potential in the ROI (Fig. 1). This outermost set of electrodes is tube shaped, as in the lpctrap geometry. Two Einzel lens triplets, used to focus the incoming beam and to clean the trap after a measurement cycle, were also included in the simulation (Fig. 6). The shape of these lenses and of the outermost tube rings were optimised separately to fit some mechanical constrains and were kept fixed along the whole optimisation process of the two innermost electrode pairs. Because of their distance from the ROI they did not have any significant effect on the potential as will be seen in Sec. 7.

To ease later comparisons with lpctrap performances, it was decided to use, in the simulation, the same applied RF voltage of 60 V on the innermost pair of electrodes even if during moratrap operation we shall use 100 V or more (all other electrodes are grounded during the trapping period).

To summarise, the parameters related to outermost electrodes and Einzel lens triplets were already fixed prior to the whole minimisation process and we have chosen to optimise seven free parameters linked to the geometry of the internal and middle electrodes (Fig. 4): the axial distances Z min,12 and Z min,34 , the radial distances R min,12 and R min,34 , a common cone angle θ and the thicknesses T h 12 and T h 34 . The minimisation of the objec-tive function f (a) (Eq. ( 16)) was performed within rootminuit [START_REF] Brun | root -an object oriented data analysis framework[END_REF][START_REF] James | minuit Function Minimization and Error Analysis: Reference Manual Version 94[END_REF] and to further reduce the parameter space to be explored, some lower positive limits were determined by exploratory simulations and set on the different distances Z min and R min . The optimisation requires the computation of the potential harmonic spectrum. To this purpose, the convergence radius R 0 of the expansion was set to 10 mm which is about 5 times larger than the cloud radius observed with lpctrap. In addition, for the evaluation of ρ 2% from Eq. ( 15), we have chosen to truncate the expansion at n max = 18, as the contribution from the A 18 term, belonging to the quadrupole series, is already very small. Indeed the absolute difference between the potential computed in electrobem and the potential synthesised from harmonics up to order 18 is smaller than 2.5×10 -7 V for the 64800 points used to compute the harmonics while the potential at these points varies from 15 to 37 V. Harmonics coefficients for moratrap up to order 18, extracted from potential sampled on a circle of radius ρ 2% = 7.39 mm using axielectrobem, when applying 60 V on the innermost electrodes and 0 V on all the other electrodes.

Rings Parameter R 1 /R 2 R 3 /R 4 R 5 /R
Positive and negative harmonics have been separated in order to display them with logarithmic scale.

After several iterations, minuit converged to a minimum of the objective function at ρ 2% = 7.39 mm. The optimised geometric parameters of the electrodes corresponding to this ρ 2% are presented in Table 1 and a cross section view of moratrap is shown in Fig. 6. The precision on distances is limited to the machining tolerance (10 µm). Effects of this limited precision are studied in Sec. 7. The radial angular acceptance Ω R , as seen from the trap centre, is 35.9 • . It is limited by the axial position of rings R 3 and R 4 . The axial angular acceptance Ω A is 55.4 • and is limited by the inner radius of ring R 1 or R 2 . This acceptance is large enough to install the two annular detectors (with a hole radius of 6 mm) before and after the Einzel lens triplets at 70 mm from the trap centre [START_REF] Delahaye | The MORA project[END_REF].

The harmonic coefficients A n (n ∈ [0, 18]) of the potential, computed at ρ 2% , are presented in Fig. 7. The ratios A4 A2 , A6 A2 and A8 A2 are respectively around -3.2×10 -7 , -1.9× 10 -2 and -3.3 × 10 -4 , they confirm the dominance of the quadrupole term A 2 . Top panel in Fig. 8 shows the potential (without the constant term n = 0) in the trapping region of moratrap for ρ ρ 2% . Middle and bottom panels demonstrate how small the contribution of these higher multipoles (n > 2) is. Since the octupole term A 4 is of the order of 10 -6 , the largest contribution comes from the dodecapole term A 6 . Bottom panel of Fig. 8 presents in more details the contribution from harmonics with orders greater or equal to 8: they are even lowered by about a factor 10 to 100 depending on the angular position. These results demonstrate that the minimisation procedure led to a quadrupole potential with the required quality, a relative difference to A 2 lower than 2% (see Eq. ( 15)), in a region with a radius as large as 7.39 mm.

In next section we shall compare these performances to those of lpctrap.

Comparison of MORATrap and LPCTrap

It is interesting to compare moratrap performances with those of lpctrap. For this purpose, the potential generated in lpctrap (Fig. 5) was simulated in axielectrobem and the harmonics were determined as seen in Fig. 9, where obviously the octupole term is much larger than in moratrap. One should however notice that A 4 and A 6 have opposite signs, thus partly compensating each others in the low radii region before A 6 takes over at larger radii. The ratios A4 A2 , A6 A2 and A8 A2 are respectively around +5.5 × 10 -3 , -1.4 × 10 -2 and +9.0 × 10 -5 . The corresponding radius ρ 2% is equal to 4.36 mm which is about 40% smaller than for moratrap. This clearly demonstrates that we have succeeded to broaden the trapping region in moratrap and to improve the potential quality compared to lpctrap, given some maximal constraints in dimension for the Paul trap. In addition, due to its tube-shape electrodes, lpctrap has an axial angular aperture of 43.6 • whereas for moratrap Ω A = 55.4 • , i.e. 27% larger. The radial angular acceptances Ω R are similar in both traps, respectively 38.0 • and 35.9 • with a small advantage of 6% for lpctrap. However, it is important to emphasise that in a trap, the amount of particles which can be effectively trapped is directly proportional to the depth of the potential well i.e. the value of A 2 and to the r adius ρ 2% (Eq. ( 6)). One can note that in moratrap, at ρ = 4.36 mm, to reach the same value of the quadrupole term as in lpctrap, would require to apply 138 V instead of 60 V like in the simulations presented throughout this paper. To deepen this comparison it is interesting to make a parallel with an ideal Paul trap, introduced in Sec. 1. For real Paul traps, one can define an equivalent internal radius r 0 as defined for ideal Paul traps in Eq. ( 1), representing the radius at which the RF voltage would be applied to yield the same A 2 value at a given ρ 2% :

r 2 0 = ρ 2 2% V 0 A 2 . ( 17
)
Compared to an ideal Paul Trap, moratrap has an effective internal radius r 0 significantly larger than the one of lpctrap which explains the relative difference in A 2 values at a given radius for both traps. From the A 2 , ρ 2% and V 0 values shown in Figs. 7 and9, one finds an effective internal radius of 19.24 mm for moratrap, and of 12.86 mm for lpctrap. This is consistent with what was found in [START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF] and is summarised in Table 2. It is important however to note that simply scaling lpctrap to the same internal radius as moratrap would not have yielded the same improvements: the 2% radius would have been enlarged to only 6.53 mm, compared to 7.39 mm in the case of moratrap, and the axial angular acceptance Ω A would still have been limited by 27% compared to moratrap. The trap capacity of moratrap, estimated thanks to Eq. ( 6) and considering r eff ρ 2% , has been enlarged by more than a factor of 2 compared to the original lpctrap, and by 40% compared to a scaled version of lpctrap. Table 2. Comparison of moratrap and lpctrap when both traps are assumed to behave like ideal Paul traps exhibiting a pure quadrupole potential. The maximum capacity was estimated for a pseudo potential model with a Mathieu parameter qz = 0.4 and for singly charged positive ions. After rescaling lpctrap effective internal radius to the r0 of moratrap, its radius ρ 2% is still smaller than the one obtained for moratrap and its maximum capacity is about 30% smaller.

To further compare both traps, one can investigate the potential non uniformity as a function of the space coordinates. We can for instance define a radial and an axial non uniformity respectively by:

U Radial (r) = 18 n=4 r R 0 n-2 A n P n (0) A 2 P 2 (0) at z = 0, (18)
U Axial (z) = 18 n=4 z R 0 n-2 A n A 2 at r = 0, (19)
where the Legendre polynomial P n (cos θ) present in Eq. ( 12) is equal to unity along z-axis and to P n (0) along given by Eq. ( 18) and Eq. ( 19), as a function of the radial and axial distance from moratrap centre.

r-axis. These quantities are presented in Fig. 10 where, in agreement with the values of ρ 2% found previously, the region where the radial and axial non uniformities are very weak is wider in the case of moratrap.

Another way to visualise the difference between the two traps is to draw a 2D uniformity defined as:

U (r, z) = 18 n=4 r R 0 n-2 A n A 2 P n (cos θ) . (20)
where the absolute value has been chosen such that a logarithmic scale may be used to help distinguishing more details. Figure 11 shows U (r, z) for ρ 4.36 mm in the case of lpctrap (top panel) and ρ 7.39 mm in the case of moratrap (bottom panel). The dashed black circle on bottom panel represents the limit of the region where the potential non uniformity is smaller than 2% in lpctrap.

All these different definitions of the non uniformity conclude to the improvement of the potential quality and of the trapping region volume in moratrap compared to lpctrap.

In next section, we shall investigate how stable is the potential against misalignment and other mechanical characteristics.

Design sensitivity

It is important to check how the mechanical tolerances in the machining and assembling of the different electrodes of moratrap can affect the trapping potential. For this purpose, the effects on the potential due to mechanical precision, in terms of machining defects or misalignment of one or several electrodes, were investigated using both axielectrobem (2D) and electrobem (3D). axielectrobem was preferred in case of defects not breaking axisymmetry, whereas electrobem was used in other cases. When performing simulations within electrobem, the setup needs to be meshed not only in the axial direction, but also angularly around this axis: a circle is approximated by a polygon and thus any cone, disk or cylinder is segmented in many triangular or quadrangular cells. As mentioned in Sec. 2, the computer RAM usage rapidly increases when working in 3D compared to axially symmetric simulations. On the computer used, the available RAM allowed us to divide the setup in 84 cells in angle while keeping the same axial segmentation as the one used in axielectrobem leading to 84 × 414 = 34776 cells in total. Such a meshing has consequences on the computed harmonics spectrum. Figure 12 shows the ratio of the harmonics from the 2D axisymmetric case to the 3D one. There is a very good agreement between the two: it is better than a few percents except for the octupole term A 4 , which appears to be very sensitive to the angular discretisation and dramatically changes not only in absolute value, but also in sign. Nevertheless, the 3D computation leads to a value of A 4 still about a factor 100 less than A 6 and its contribution to the potential is therefore very limited. When using 3 planes of symmetry in electrobem, one can simulate 1/8 of moratrap with as many cells as a full simulation without symmetry planes. In that case, the 3D harmonics spectrum agrees with the 2D one within a few 10 -6 -10 -5 even for A 4 . However, most of the envisaged defects break at least one of these symmetries. Consequently, the full moratrap has to be simulated in these studies and the results should be compared to those from the perfectly aligned 3D version of moratrap rather than to those from axielectrobem.

The precision on the evaluated potential essentially depends on the characteristics of the electrode mesh and is constant once this mesh has been fixed. A single mesh is used to test different trap configurations with or without misalignment. The results obtained for configurations exhibiting some defects are compared to a reference one where the electrodes are perfectly aligned and machined (in the limit of the mesh discretisation) and therefore one expects that the observed differences are highly significant overall.

The considered defects included machining tolerances of δr and δz 20 µm, axial and radial translations of ∆x and ∆z 200 µm and rotations around r-axis of ∆θ 0.2 • . All these values are larger than what can be mechanically achieved for moratrap. These defects were simulated for either one half of the trap in a whole or for electrodes taken individually while keeping the rest of the trap unchanged.

When dealing with individual electrodes, because of their relatively large distance from the trap ROI, all types of defects related to R 6 , E 4 , E 5 and E 6 have been shown to have very minor influence on the harmonics spectrum and therefore on the potential quality. This was also the case for possible machining defects on R 1 , R 2 , R 3 and R 4 as the computer numerical control machining precision was assumed better than 20 µm. The influence of misalignment along the z direction was studied with axielectrobem by varying electrode(s) axial position. An axial displacement ∆z of half of the trap conserves both the axisymmetry and the planar symmetry: such a displacement shifts the trap centre by ∆z/2 and therefore the potential still exhibits a planar symmetry around the position of this shifted centre. The potential harmonics expansion should therefore be determined around this shifted centre to avoid inducing spurious harmonics with odd orders. In the particular case of ion traps, the ion cloud barycentre will follow the potential centre and, in case of a large shift, this may modify the acceptance and have some consequences both on the detector counting rates and on the measured β decay asymmetry. This shall be addressed in other simulations. Figure 13 presents the harmonics spectrum computed around a point located at z = ∆z/2 = -50 µm for a global translation of all electrodes located at negative z by ∆z = -100 µm. One clearly sees that odd harmonics are not present in this spectrum. In this case, the potential quality is still suitable to trap ions as may be seen in Fig. 14 which shows the evolution of ρ 2% , A 2 and the non uniformity (Eq. ( 15)) versus the translation ∆z. A non uniformity of 2.5% is reached for displacements as large as 0.4 mm and even for such large ∆z, ρ 2% is still about 6.9 mm and A 2 is only lowered by less than 10%.

To describe and further compare the harmonics induced by all the studied translations and rotations of different setup parts, we shall use the spherical harmonics with degree n and order m (Eq. ( 11)) instead of the cylindrical harmonics of order n (Eq. ( 12)) used so far. Indeed, breaking the axisymmetry and/or the planar symmetry may induce some harmonics with different degrees n but also and above all with non zero order m on the contrary to what has been seen previously. As an example, Fig. 15 shows the harmonics corresponding to a setup where half of the trap has been rotated by ∆θ = 0.2 • around x-axis. One should notice that even if the po-Figure 14. Evolution of the radius ρ 2% , the quadrupole term A2 at ρ 2% and the uniformity as a function of ∆z, the translation distance of the left half of the trap (at z < 0) from its nominal position. These results were obtained using axielectrobem. tential centre has probably been shifted away from the reference frame origin, the harmonics expansion is still performed around this origin. The largest harmonics generated by the rotation are however a factor 10 smaller than the main ones and should therefore weakly impact the quality of the trapping potential. To estimate their influence and similarly to what was done in Eqs. ( 15), ( 18) and (19), we shall define a 3D non uniformity. This is more difficult as the associated Legendre functions of the first kind P m n (cos θ) are not bounded to [-1, 1] on the contrary to the Legendre polynomials P n (cos θ): their

U 3D = nmax n=1 n m=0 if n =2 m=1 if n=2 max(P m n (cos θ)) |A nm | + |B nm | |A 2,0 | , (21)
where the numerator and therefore U 3D have been maximised by setting the functions P m n (cos θ), cos(mϕ) and sin(mϕ) to their maximal values and by taking the absolute values of the harmonic coefficients. For the denominator we only keep |A 2,0 | without any spatial dependent function which might become null at some location leading to an infinite value of U 3D . When simulating moratrap in electrobem without any misalignment, the 3D non uniformity U 3D is 2.017 %, which shall serve as our reference. One should notice that this value is weakly above the 2% limit obtained with axielectrobem, this highlights the influence of the 3D meshing. The U 3D values obtained for the tested misalignment are presented in Tables 3 and4 respectively for translations (∆z, ∆x) and rotations (∆θ). These clearly demonstrate that displacements of either rings R 5 and R 6 or Einzel lens triplets have negligible influence on the potential quality within the ROI. The largest non uniformities are obtained for translations ∆x of ring electrode R 4 or for the rotation of half of the trap. Even if these non uniformities were extracted from harmonics series expansions performed around the reference frame origin instead of around the potential centre, they still remain smaller than 2.5%, i.e. at a level sufficient to maintain the trapping efficiency [START_REF] Delahaye | The open LPC Paul trap for precision measurements in β decay[END_REF].

All the reasonable mechanical defects studied have been shown to weakly affect the quality of the potential and the

Conclusion

We have presented a method used to optimise the geometry of the axially symmetric ion trap, moratrap, dedicated to the measurement of the triple correlation parameter D in nuclear β-decay of radioactive ions. Our starting point was the lpctrap geometry, the former used transparent Paul trap, and we succeeded in reducing the contribution to the potential from harmonics of order higher than 2. This is necessary to achieve a longer storage time and a better trapping efficiency. According to the pseudopotential approximation, the moratrap capacity has been enlarged by more than a factor of 2 compared to the initial lpctrap geometry. The optimised geometry exhibits a larger axial angular acceptance than in lpctrap. This larger acceptance is necessary for the online monitoring of the trapped ion cloud polarisation, by detecting β particles along the axis of the trap. Further simulations are currently being performed to deeply investigate the ion cloud dynamics, the trapping efficiency as well as the trapping time. The optimised trap electrodes could easily be machined with a mechanical precision of the order of 10 µm.

The whole setup should be soon tested and then installed at JYFL accelerator facility in University of Jyväskylä (Finland) before its final operation at GANIL -DESIR (Caen, France). The optimisation method described in this study could easily be applied to other axisymmetric ion traps or electrodes setups. Different types of fitness functions could bring many improvements to the method and should be explored in detail.

This work was financially supported by Région Normandie via its Réseaux d'Intérêts Normands. The authors would like to thank their collaborators from the LPC Caen CAD group and workshop for their deep involvement in the design, manufacture and assembly of moratrap.

A Homogeneous polynomials in cylindrical axisymmetric coordinates

For completeness, in this appendix, we present a series expansion for the potential of an axisymmetric system in cylindrical coordinates (r, z). Starting from the coordinates conversion relations ρ 2 = r 2 + z 2 and cos θ = z/ρ and introducing the homogeneous harmonic polynomials H n (r, z) given by an explicit relation rather than by a recurrence one: the potential given by Eq. ( 12) is rewritten as :

V (r, z) = ∞ n=0 A n ρ n 0 H n (r, z). ( 23
)
This series expansion converges for points (r, z) satisfying √ r 2 + z 2 ρ 0 with ρ 0 the convergence radius. Again, a planar symmetry at z = 0 allows to further reduce the multipole expansion to even harmonics only. The first harmonic polynomials are listed up to order 10 in Table 5.

n Hn(r, z)

0 1 1 z 2 1 2 2z 2 -

Figure 1 .

 1 Figure 1. moratrap setup overview showing the location of the ion trap in the centre surrounded by 4 beta detectors (Phoswich) and 4 recoil ion detectors (MCP).

Figure 2 .

 2 Figure 2. gmsh cross section view of moratrap showing the electrodes along with iso-potential lines obtained with electrobem. The trap consists of three ring electrode pairs (R1 -R2, R3 -R4 and R5 -R6) surrounded by two Einzel lens triplets (E1 -E3 and E4 -E6). A more detailed description is given in Sec. 5. The thin black lines on the electrodes delimit the perimeter of the polygonal cells used to solve the electrostatic problem.

Figure 3 .

 3 Figure 3. root view of the potential along moratrap axis (r = 0) computed in axielectrobem for 60 V applied on inner electrodes R1 and R2, all other ones being grounded.

Figure 4 .

 4 Figure 4. Definition of the geometric parameters used to optimise the trap: axial distance Zmin from the radial axis to the closest point on a ring electrode, minimal radial distance Rmin, radial thickness of the electrode T h, gap G between two successive electrodes and cone angle θ from z-axis. ΩA and ΩR are respectively the axial and radial angular openings. The circle of radius R0 used to compute the spherical harmonics is also drawn.

Figure 5 .

 5 Figure 5. Cross section view of lpctrap geometry. There are three pairs of ring/tube electrodes with conical extensions and cylindrical mechanical support. R1 and R2 are the RF electrodes, R3 and R4 are used to inject and extract the ions and the external electrodes R5 and R6 are used to correct the potential generated along the r-axis and to steer the trajectories of the decay products. The circle of radius R0 = 10 mm gives the overall scale.

Figure 6 .

 6 Figure 6. Cross section view of the optimised moratrap geometry. A RF voltage is applied to the innermost ring pair (R1, R2), while the middle and outer electrodes, as well as the Einzel lenses are grounded. z-axis is along the beam injection and extraction. The geometry is symmetric with respect to the x -y plane and is invariant around the z axis. The circle of radius R0 = 10 mm gives the overall scale, but more precise dimensions are given in Table1.

Figure 7 .

 7 Figure 7. Harmonics coefficients for moratrap up to order 18, extracted from potential sampled on a circle of radius ρ 2% = 7.39 mm using axielectrobem, when applying 60 V on the innermost electrodes and 0 V on all the other electrodes. Positive and negative harmonics have been separated in order to display them with logarithmic scale.

Figure 8 .

 8 Figure 8. Contribution to the potential for ρ 7.39 mm in moratrap from different harmonics sets (vertical scale is in volts). Top panel corresponds to harmonics with n 2 and exhibits a strong quadrupole dependence. Middle panel is for n 4 and clearly shows the dodecapole contribution. Bottom picture shows the contribution from all harmonics with n 8.

Figure 9 .

 9 Figure 9. Harmonics coefficients for lpctrap up to order 18, extracted from potential sampled on a circle of radius ρ 2% = 4.36 mm using axielectrobem, when applying a 60 V RF voltage on the innermost electrodes.

Figure 10 .

 10 Figure 10. Radial and axial non uniformities, respectively given by Eq. (18) and Eq. (19), as a function of the radial and axial distance from moratrap centre.

Figure 11 .

 11 Figure 11. Top: Potential non uniformity obtained with Eq. (20), for lpctrap region of interest, i.e. for ρ = √ r 2 + z 2 ρ 2% = 4.36 mm. Bottom: Same figure for moratrap with ρ 2% = 7.39 mm. The dashed black circle corresponds to the lpctrap 2% radius.

Figure 12 .

 12 Figure 12. Comparison between harmonics computed with axielectrobem and electrobem with an angular/polygonal segmentation of 84 cells around z-axis. Except for A4, all harmonics are equal to better than a few percents, especially for low order ones.

Figure 13 .

 13 Figure 13. Harmonics spectrum obtained when translating left half of moratrap by ∆z = -100 µm. The expansion has been computed around a centre located at half of the displacement (z = -50 µm) where the translated potential centre stands, i.e. where the ion cloud barycentre will be.

Figure 15 .

 15 Figure 15. Harmonics spectrum (without the constant term A0,0) obtained at ρ 2% with electrobem for half of the trap rotated by 0.2 • around x-axis (vertical scale is in volts). The series expansion centre is kept at the frame origin which does not correspond to the minimum of the potential after rotation. Harmonics with absolute values below 10 -10 have been set to zero as their contribution to the potential is completely negligible.

  volume of the trapping region in moratrap if the machining process and the alignment of the different electrodes during assembly are kept under control. It was therefore decided to design, machine and assemble the moratrap electrodes at LPC Caen. A picture of this realisation is shown in Fig.16.

Figure 16 .

 16 Figure 16. Picture of moratrap on its mechanical support. It has been machined and assembled at LPC Caen.

H

  n (r, z) = ρ n P n (cos θ) k (n -2k)! (k! ) 2 r 2k z n-2k , (22)

Table 1 .

 1 Optimised moratrap electrodes dimensions.

	6

Table 3 .

 3 3D non uniformity due to ∆z and ∆x misalignment of the half trap and of individual electrodes.

	Translated	∆z	U 3D	∆x	U 3D
	part	[µm]	[%]	[µm]	[%]
	Half trap	100	2.127	200	2.398
	R2	100	2.182	200	2.208
	R4	100	2.142	200	2.461
	R6	100	2.017	200	2.017
	E4	100	2.017	200	2.017
	E5	100	2.017	200	2.017
	E6	100	2.017	200	2.017
	Rotated ∆θ U 3D
	part		[ • ]	[%]
	Half trap	0.2 2.404
	R2		0.2 2.267
	R4		0.2 2.159
	R6		0.2 2.018
	E4		0.2 2.017
	E5		0.2 2.017
	E6		0.2 2.017

Table 4 .

 4 3D non uniformity due to ∆θ misalignment. Here the series expansion centre is kept at the frame origin. Realistic mechanical rotations of at most 0.2 • lead to non uniformities well below 2.5%.

Table 5 .

 5 r 2 31 2 2z 3 -3r 2 z 41 8 3r 4 -24z 2 r 2 + 8z 4 5 1 8 8z 5 -40r 2 z 3 + 15r 4 z 6 1 16 -5r 6 + 90z 2 r 4 -120z 4 r 2 + 16z 6 7 1 16 16z 7 -168r 2 z 5 + 210r 4 z 3 -35r 6 z 8 1 128 35r 8 -1120z 2 r 6 + 3360z 4 r 4 -1792z 6 r 2 + 128z 8 9 1 128 128z 9 -2304r 2 z 7 + 6048r 4 z 5 -3360r 6 z 3 +315r 8 z 10 1 256 -63r 10 + 3150z 2 r 8 -16800z 4 r 6 + 20160z 6 r 4-5760z 8 r 2 + 256z 10 Harmonic polynomials of degree n 10 in axisymmetric cylindrical coordinates (r, z).

In axisymmetric problems, the matrix elements Qij and other surface integrals are computed differently making use of complete elliptic integrals.

In reality some parts of lpctrap electrodes are not axisymmetric at large distance from the ROI, but we simplified their geometry and made them cylindrical in the simulation.