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

Ludovic Apvrille
email: ludovic.apvrille@telecom-paristech.fr

Test Sequence Generation From Formally Verified SysML Models

Test generation has been acknowledged as a costprone activity reducing productivity and time to market. The expected benefits of Model Based Systems Engineering include automated generation of test sequences from models. The paper proposes verification solutions for the System Modeling Language (SysML). In particular, the paper shows how to link test generation to formal verification. The proposed algorithms are implemented by the free software TTool. Two case studies support discussion on conformance and interoperability testing, respectively.

I. INTRODUCTION

The widespread of Model Based System Engineering approaches [START_REF] Madni | Model-based systems engineering: Motivation, current status, and research opportunities[END_REF] has encouraged several communities to develop their own modelling language. For instance, the Object Management Group (OMG) [START_REF]Object Management Group[END_REF] and the International Council for Systems Engineering (INCOSE) [3] have jointly developed and standardized the Systems Modeling Language (SysML) [START_REF]Systems Modeling Language, version 1.5[END_REF]. The benefits and potential of using SysML have been acknowledged in several application domains, in particular avionics [START_REF] Sergent | Benefits of Model Based System Engineering for Avionics Systems[END_REF], [START_REF] Mehnni | Flight Control System Modeling with SysML to Support Validation, Qualification and Certification[END_REF], [START_REF] Peleska | Model-based Avionic Systems Testing for the Airbus Family[END_REF]. SysML is now supported by proprietary [START_REF]Cameo Systems Modeler[END_REF], [START_REF] Rhapsody | [END_REF] and open-source tools like Papyrus 1 or TTool 2 that help automating an important variety of activities (e.g. [START_REF] Saqui-Sannes | EarlyChecking of SysML Models applied to protocols[END_REF]) throughout the design trajectory of complex systems. One of these activities is test sequence generation, an activity that is cost-prone and time consuming, and therefore worth being automated to reduce time to market of complex systems. Automated test generation from SysML models has already been discussed in [START_REF] Gauthier | A SysML Formal Framework to Combine Discrete and Continuous Simulation for Testing[END_REF], [START_REF] Faria | Practical Results from the Application of Model Checking and Test Generation from UML/SysML Models of On-Board Space Applications[END_REF], [START_REF] Lasalle | Development of an Automated MBT Toolchain from UML/SysML models[END_REF]. The paper proposes another SysML-approach based on earlier work with the Formal Description Technique Estelle [START_REF] Saqui-Sannes | Verification by abstraction as a preamble for interoperability test suite generation[END_REF]. In brief, the reachability graph of the SysML model is computed directly from the block instance diagram that defines the architecture of the system and from the state machine diagrams that define the behaviors of the block instances. The transitions of the graph are labelled using the messages exchanged by pairs of block instances. This transforms the reachability graph of the SysML model into a Labeled Transition System. Theories developed for LTS [START_REF] Konur | Towards an Integrated Approach to Verification and Model-Based Testing in System Engineering[END_REF] therefore apply, in particular minimization that outputs a quotient automaton computed with respect to an equivalence relation such as Milner's observational equivalence [START_REF] Dssouli | Chapter Three -Testing the Control-Flow, Data-Flow, and Time Aspects of Communication Systems: A Survey[END_REF]. Taking a quotient automaton as input, all the paths of maximal size leading to a termination state are computed using the Dijkstra technique [START_REF] Ouerdi | EMV Card: Generation of Test Cases based on SysML Models[END_REF] to finally 1 Papyrus: https://www.eclipse.org/papyrus/ 2 TTool: https://ttool.telecom-paristech.fr/ obtain a set a test sequences. Associated algorithms have been implemented in the free and open-source tool named TTool [START_REF] Gauthier | A SysML Formal Framework to Combine Discrete and Continuous Simulation for Testing[END_REF]. The paper is organized as follows. Section II introduces SysML and more precisely the block instance and state machine diagrams. Section III presents a verification by abstraction approach for SysML models. Section IV uses the output of the verification process to generate test sequences. Section V applies the proposed approach to a client/server protocol and discusses interoperability testing. Section VI uses a UAV in charge of taking pictures to address conformance testing. Section VII surveys related work. Section VIII concludes the paper and outlines future work.

II. SYSML

The SysML standard [3] defines nine type of diagrams that may be used inside one model to cover the requirement capture, analysis and design phases in the trajectory of systems. For the test generation approach discussed in the paper, only the design phase is of interest. During the design phase, one defines the architecture of the system and the behaviors of the blocks the architecture is made up of. The SysML standard defines two architectural design diagrams: the Block Definition Diagram (BDD) and the Internal Block Diagrams (IBD) of the SysML standard. In the paper, we use the version of SysML supported by TTool where the BDD and IBD are merged into one diagram: the Block Instance Diagram (BID). Each block instance has a behavior expressed in the form of a SysML state machine diagram.

A. Block Instance Diagram

A Block Instance Diagram is a tripartite graph with one type of nodes, each one defining a block instance, and two types of relations. First, the "composition" relation of SysML (black diamond) enables to say that one block instance is made up of one or several block instances. Second, the "connect" relation connects two ports used for exchanging signals. For simplicity, the below definitions address the "connect" relation, not the "composition" one.

Definition: block instance. A block instance is a 4-uple (id, al, ml, pl, isl, osl) where:

• id is a String that names the block instance.

• al is an attribute list. The attribute types include Integer, Boolean, Timer, and user-defined Records. An attribute may be defined with an initial value. • ml is a method list.

• pl is a port list.

• isl is an input signal list.

• osl is an output signal list. Thus, it heavily depends on the experience of the user of the simulator, or on the random generator of that simulator. Conversely, formal verification relies on mathematics rather than chance.

A. Working Hypothesis

The verification approach discussed in this section requires exploring the state space of the SysML and more precisely generating (a sub part of) the graph of states that may be reached from the initial state of the system. Full generation of the reachability graph of the SysML model explores all the execution paths the system may go through starting from its initial state. The graph enables checking the system against reachability properties (using the model-checker we have integrated into TTool). Also, when a property to investigate is the reachability of a given state, then the graph generation stops as soon as the property is satisfied (on-the-fly modelchecking). As usually, reachability analysis faces the state space explosion problem. Starting from now, the following hypothesis applies: the reachability graph of the SysML model may be entirely built within an acceptable amount of time and memory.

B. Labeled Reachability Graph

In this section, we assume the SysML model is composed of one block instance diagram and one state machine diagram per block instance. Other diagrams of the SysML model are ignored since they do not play any role in the reachability graph construction process.

Definition: Labeled Reachability Graph. The Labeled Reachability Graph (LRG) of a SysML model is a 5-uple LRG = (S, Σ, Θ, ∆, s 0 ) where:

• S denotes the countable number of global states of the SysML model • Σ is a countable set of observable events included in the list of signals exchanged on ports. Thus we have Σ ⊆ "!"A o (p) ∪ "?"A i (p) ∪ "!"S o (p) "?"S i (p)) where:

-"!" and "?" respectively denote an emission and a reception. -A o (p) is an output signal connected to an asynchronous channel. -A i (p) is an input signal connected to an asynchronous channel. -S o and S 1 are output and input signals connected via a synchronous channel. p denotes a list of attributes of the block to which the signal sending/receiving corresponds.

• Θ is a countable set of internal events, e.g. the assignation of an attribute. Thus we have

Θ ⊆ i(a =< expr >) ∨ i()
where a is an attribute. • ∆ ⊆ Sx(Σ ∪ Θ)xS denotes the set of labeled transitions in the reachability graph. A label contains the name of the signal exchanged by one pair of block instances or the identification of the internal action, i.e. i(...). An empty i() transition is sometimes denoted τ.

• s 0 denotes the initial state of the SysML model, i.e. the global state obtained when all the block instances enter their initial states. Generating the LRG consists in considering all possible transitions from all not-yet-handled states -this not-yethandled set contains s 0 when starting the graph generation-. For each not yet handled transition t fire-able from a given state s 1 ∈ notyethandled, we first create δ ∈ ∆ with δ = (s 1 , σ , s 2 ). We then have to compute whether ∃s ∈ S with s ≡ s2. If this is the case, δ = (s 1 , σ , s 2 ). Otherwise, s 2 is added to S and to the not-yet-handled set. This leads us to define the notion of state and of state equivalence for a SysML model. In terms of tool implementation, the above approach is implemented by the model-checker integrated to TTool. The comparison between states is based on hashing techniques taking into account all elements above listed: state machine pointer, attributes values, local clock, message queues. Finally, contrary to many contributions, we are able to generate LRG directly from SysML model, without the need to use a pivot language. This facilitates the back-tracing to models, and avoids formally proving model transformations. This probably has a performance cost when generating the reachability graphs or more generally studying safety properties.

C. LRG Minimization

The labeled reachability graph of the SysML model may have hundred, thousands and even more states and transitions. Definitely, interpreting such a graph is impossible for Human. Also, the purpose of verification is not to check the entire behavior of the system in one operation. The question of which subset of the system can be defined and checked is then asked. The answer proposed by verification by abstraction is as follows: verification will zoom on a subset of the signals exchanged by pairs of blocks involved in the evolution of the system. To refer to the theories developed for Labeled Transitions Systems, these signals of interest play the role of observable events. Other signals become de facto invisible events. Transitions involving one signal exchange are labeled by the name of the signal. Other transitions are labeled by τ. The objective of minimizing the Labeled Reachability Graph is to get rid of the transitions labeled by τ and to keep the ones labeled by a signal name. Depending on the equivalence relation, the quotient automaton resulting from the minimization process may still contain τ transitions.

TTool implements a 3-step algorithm: 1) Replace each ignored actions with a τ action.

2) Remove all τ transitions i.e. merging states s1 and s2 when there is a τ transition between the two. 3) Minimize the graph. The quotient automaton computation algorithm can be sketched as follows. If the first two stages are quite straightforward from an algorithmic point of view, the minimization itself relies on the identification of coarse blocks, i.e. blocks that are bisimulation-equivalent states. In order to identify these blocks, we rely on the partitioning of the graphs -using splittersiterating on the different symbols in the graph, as explained in [START_REF] Kilada | A C++ Implementation of an Efficient Algorithm for Labeled Transition System Minimization Based on Bisimulation Equivalence[END_REF].

IV. TEST SEQUENCE GENERATION

The paper proposes a test sequence generation where the test sequences are generated from the quotient automaton output by the verification process. Since test sequences are built directly from SysML actions and communication labels, the labelling of test sequences directly refers to elements of the SysML model.

A. Refusal Graphs and Test Sequence Generation

A refusal graph [START_REF] Drira | The refusal Graph: a Tradeoff between Verification and Test, 6th International Workshop on Protocol Test Systems[END_REF] is a deterministic Labeled Transition System that emphasizes on the actions that are accepted from a state (and thus on the ones that are refused).

Definition: Refusal Graph. A Refusal Graph is a 5-uple RG = (G, Σ, ∆, g 0 , Re f ) where Re f : G → P(P(Σ)) is a mapping which defines for each state, the sets of actions that may be refused after the sequence leading to this state. To avoid redundancy, refusal sets must be minimal w.r.t. their inclusion set. Also, to avoid describing imaginary systems, only refused parts of the output set are considered.

The algorithm used to generate a Refusal Graph from a Quotient Automaton can be sketched as follows. It follows all possible paths in the input Quotient Automaton, but stops each time an already met state is encountered (cycle). A transition with a given label "l" is created in the RG each time the outgoing transitions of the QA contains this label at least once from the current state of QA.

Using refusal graphs enables using an operational procedure for implementing the concept of canonical tester [START_REF] Richards | Testing Solutions through SysML / UML, Incose[END_REF]. The latter and the original model have the same traces. Further, its synchronization with the Implementation Under Test must not lead to a deadlock situation.

Finally, generating Test Sequences TS from a Refusal Graph basically consists in identifying all possible paths of maximal size i.e. it consists in finding, from termination states, the longest path from the origin state. To do this, we rely on the Dijkstra algorithm that can compute the path length from one state to another in a LTS. Therefore, test sequences correspond to all longest paths from the initial state to a termination state.

Definition: Test Case. A Test Case is a Labeled Transition System TC = (S, Σ, ∆, s 0 , v) where v : S → {pass, f ail}. The algorithm deriving a refusal graph from one quotient automaton is described below. The following case studies demonstrate the usability of our verification approach on a relevant system. The interest of generating graphs and test sequences from SysML model is shown by the clear relation between graph labels and the SysML model.

V. A CASE STUDY OF INTEROPERABILITY TESTING

A. Testing architecture

Layered design of communication protocol commonly uses a 3-layer pattern where two protocol entities rely on one preexisiting communication service to render in turn a value-

B. Client/Server Protocol

This section abstractly defines a client/server protocol where the server may acknowledge or refuse the requests issued by the client. The use-case diagram in Figure 2 delimits the boundary of the communication system and links it to the client and server applications, as well as to underlying, preexisting network.

Figure 3 depicts a successful completion of an inquiry procedure. The service primitive suffixes abbreviate "request", "indication", "response", and "confirm". Two Protocol Data Units are used: xREQ and xOK. For space reasons, the paper does not show the sequence diagrams developed for the inconclusive termination of the request. In brief, the sequence diagrams of Figures 3 and4 are merely modified to replace "ok" by "nok" and "OK" by "NOK".

Figure 5 depicts the communication architecture that use the client/server protocol.

In terms of protocol machines, Figure 6 depicts the protocol machine associated with the Requester. As far as communication architecture validation is concerned, the events to be preserved by the minimization process are the service primitives exchanged at the boundary between the protocol entities (Requester and Responder on Figure 5) and their respective users (Client and Server on Figure 5). We selected service primitives as observable events to decorate the reachability graph in Figure 8. The minimization process with respect to observational equivalence outputs the quotient automaton depicted by Figure 9. From the quotient automaton in Figure 9, we obtained the test sequences depicted by Figure 10. If the system stops in an intermediate state (e.g., a state different from 4 and 8), then

VI. A CASE STUDY OF CONFORMANCE TESTING: A UAV

A. Testing Architecture

Conformance testing is the process of verifying the correctness of an artifact in the development cycle of a system against its model. For a black box testing approach one may use a testing architecture where the tester accesses the Implementation Under Test (IUT) via one or several Points of Control and Observations (PCOs). Assuming the design is a layered one, Figure 11 depicts such type of conformance testing architecture with two PCOs that respectively test the interfaces of the IUT with its upper and lower layers.

B. Informal Specification of the UAV

The UAV can autonomously take off, fly in a stabilized way, and land at its destination or whenever a critical situation is encountered. It takes pictures at given locations. Only the software related to taking pictures is modeled in this case study: the taking off, flying and landing actions are not Pictures can be taken only when the drone is flying. A remote system located in a ground station can send picture order to the drone. A picture order contains the GPS position of the picture to be taken. To know its current position, a drone has an integrated GPS. When a picture GPS point is reached, with regards to a given threshold, the picture is taken, and then stored on a CompactFlash removable storage system. The system needs 2 seconds to take a picture, and between 4 and 5 seconds to store it in on the memory card. Pictures may be remotely downloaded from the ground station using a download order. Pictures can also be read from the CompactFlash once the drone has come back from its mission. Figure 14 depicts the test sequences generated by TTool: one of them concerns the regular execution of the UAV , the other one evaluates that that the UAV works as expected when the memory for storing pictures is full. These test sequences remain abstract test sequences that need to be transformed into concrete test sequences that may be processed by an actual tester. This transformation is not yet implemented by TTool.

C. UAV Modeling and Test Sequence Generation

VII. RELATED WORK

A. Model-Based Testing for Basic Models

How to derive test sequences from fundamental models such as Extended Finite State Machine or Labeled Transition System has extensively been discussed in the literature, years before the MBSE acronym was coined. In particular, the increasing development of networked systems has stimulated research work on protocol testing In [START_REF] Dssouli | Chapter Three -Testing the Control-Flow, Data-Flow, and Time Aspects of Communication Systems: A Survey[END_REF]. Most of the models listed in [START_REF] Dssouli | Chapter Three -Testing the Control-Flow, Data-Flow, and Time Aspects of Communication Systems: A Survey[END_REF] are discrete event Model of Computation. The advent of Cyber Physical Systems has stimulated research work on hybrid models assembling discrete events and continuous paradigms.

B. Model-Based Testing from SysML models

In [START_REF] Gauthier | A SysML Formal Framework to Combine Discrete and Continuous Simulation for Testing[END_REF]

VIII. CONCLUSIONS

Test generation has been acknowledged as one of the most cost-prone activity in the design trajectory of complex systems. The expected benefits of using a model-based approach includes the possibility to reduce that cost by automating a test generation from a model of the system.

The paper proposes a novel approach for generating test sequences expressed in SysML and more precisely in the dialect of SysML supported by the free and open-source tool named TTool. The block instance diagram describing the architecture of the system and the state machine diagrams describing the behaviors of the block instances serve as input to a process that links test generation to formal verification. Indeed, assuming the reachability graph of the SysML model can be computed, TTool generates the latter as a Labeled Transition System whose transitions are labeled by events appropriately selected by the user of TTool. The labeled reachability graph is minimized using Milner's observational equivalence. The resulting quotient automaton serves as starting point to build up a refusal graph, which serves in turn as starting point for generating test sequences.

Taking a UAV as case study, the paper illustrates the test generation approach in the case of conformance testing. Other types of testing, such as interoperability testing [START_REF] Saqui-Sannes | Verification by abstraction as a preamble for interoperability test suite generation[END_REF] also deserve investigations. Whatever the type of testing, the test sequences must be presented in a standardized form, e.g. using the TTCN notation [START_REF]TTCN standard[END_REF].

Fig. 1 .Fig. 2 .

 12 Fig. 1. Global Testing Architecture (left) and Local Testing Architecture (right)

Figure 7

 7 depicts the protocol machine associated with the Responder.

Figure 8

 8 depicts the reachability graph generated from the SysML model by TTool. Making the reachability graph a Labeled Transition Systems, which enables reuse of verification techniques originally developed for LTS, in particular

Fig. 3 .Fig. 5 .

 35 Fig. 3. Service Scenario (successful request)

Fig. 6 .Fig. 7 .Fig. 9 .Fig. 10 .

 67910 Fig. 6. Protocol Machine for the Requester

Fig. 11 .

 11 Fig. 11. Conformance Testing Architecture

Fig. 12 .

 12 Fig. 12. Architecture of the UAV

Figure 12

 12 Figure 12 depicts the architecture of the UAV in the form of a block instance diagram. The reachability graph of the UAV depicted by Figure 12 has 674 states and 927 transitions. The reachability graph is labeled to preserve a limited set of signal exchanges: requesting for the flight to start, taking pictures, and saving them. The minimization process outputs the quotient

Fig. 13 .

 13 Fig. 13. Quotient Automaton from the UAV model

Fig. 14 .

 14 Fig. 14. Test Sequences for the UAV Model

uple (s start , a f ter, condition, Actions, s end ) where: • s start is the initial state of the transition. • a f ter(tmin,tmax) enables firing the transition after at least between tmin and tmax units of time have elapsed.

  connect is a function Port → Port that connects pairs of ports.

	Definition: Block Instance Diagram. A Block Instance
	Diagram is a 3-uple (Blk, connect, assoc) where:
	• assoc is a function (BlkxSignal) → (BlockxSignal) that
	associates one signal of block B1 to one signal of block
	B2, making it possible to compose the state machines
	belonging to B1 and B2. Here, the term compose denotes
	a composition in the usual sense of finite state machine
	composition, not the composition relation supported de-
	picted by a black diamond in SysML.
	Finally, a Block Instance Diagram depicts the architecture of
	a system as a graph of interconnected blocks.
	B. State Machine Diagram
	Each block instance contains one extended state finite state
	machine that supports states, transitions, attribute settings, in-
	puts and outputs operations on signals, and time manipulation.
	Definition: State Machine. An extended finite state machine
	depicted by a SysML state machine diagram is bi-partite graph
	(s0, S, T ) where
	• S is a set of states (s 0 is the initial state).
	• T is a set of transitions.
	Definition: State Transition. A transition in a state machine
	is a 5-• condition is a Boolean expression that conditions the
	execution of the transition. A condition may use attributes
	of its corresponding block.
	• Actions is a ordered set of action. These actions can be
	executed only once the transition has been enabled i.e. the
	a f ter clause has elapsed and the condition equals true.
	• s end is the final state of the transition.
	A state machine cannot contain parallel states, historic state,
	fork and join states. These behaviors can easily be replaced
	by counterparts. e.g. using sub-blocks and synchronous signals
	for each parallel activity.
	III. FORMAL VERIFICATION
	Adding a state machine to each block instance makes the
	model an executable one. Interactive simulation enables early
	checking of the state machines against design errors. It does
	not necessarily explore the entire state space of the model.

• Blk is a set of block instances.

•

  Definition: State. The state of a SysML model is defined as S b , S a where: • S b is the state of block b. The state of a block b is defined by the value of its attributes al, by one state of its state machine (it can be considered as a pointer to the current state of the state machine of b), and by the value of its clock. Indeed, if we assume a global clock applied to all blocks, we need to use a local clock to remember how much of after clauses has elapsed. Two states s 1 and s 2 are said to be equivalent s 1 ≡ s 2 if and only if all their state values are equal.

• S a represents the state of the signal queue a. There is one signal queue for each asynchronous communication. The state of a signal queue is characterized by its ordered list of n messages m 1 (p 1 1 , p 1 2 , ...), m 2 (p 2 1 , p 2 2 , ...), ...).

Definition: State equivalence.

1 :

 1 QA: Quotient Automaton. RG: Refusal Graph. 2: We go though the graph, starting at state #1. Each time we meet an already handled state, we stop handling this path.

	9:	currentQA = toHandleStates.get(0)
	10:	currentRefusal = Assoc.get(currentQA);
	11:	toHandleStates.remove(0)
	12:	metStates.add(currentQA) {For each out transition of
		currentQA, we create a transition in currentRefusal if
		the transition does not yet exist with the same label}
	13:	Loop on tr = currentQA out transitions
	14:	if not (tr.label exists in out transitions of currentRefusal)
		then
	15:	QA.add(newState)
	16:	Assoc(destination state of tr, newState)
	17:	currentRefusal.addOutTransition(to newState, label
		of tr)
	18:	if not (metStates contains destination state of tr) then
	19:	toHandleStates.add(destination state of tr)
	20:	end if
	21:	end if
	22:	Endloop
	23: end while

3: metStates = {} 4: toHandleStates = { QA.getState(0)} {We have a list of associations Assoc{} between states in QA and RF} 5: State currentState = new State(0); 6: RG.add(currentState) 7: Assoc.add(QA.getState(0), currentState) {Main loop} 8: while toHandleStates is non empty do

  , Gauthier et al use discrete SysML/OCL modeling artifact and add a Modelica support to address continuity. The SysML model is annotated with Modelica and OCL. The approach relies on simulation techniques. In [12], Lasalle et al. present MBT, a test sequence generation tool developed for SysML in the framework of VETESSS, a project connected with automotive industry. MBT generates functional tests from UML or SysML diagrams edited with the open-source and Eclipse-based tool Topcased. The SysML model of the System Under Test is the input of a tool that generates abstract test suites and transfers them to another tool in charge of creating concrete tests in the form of test scripts. In [22], Ouerdi et al. share an experience in vulnerability test generation from SysML models of smart cards. The authors


Hugo Perrin

Stefan Karpitschka

Jacco H Snoeijer

Bruno Andreotti

Peeling a soft viscoelastic adhesive: experiments and scaling laws

The functionality of adhesives relies on their response under the application of a load. Yet, it has remained a challenge to quantitatively relate the macroscopic dynamics of peeling to the dissipative processes inside the adhesive layer. Here we investigate the peeling of a reversible adhesive made of a polymer gel, measuring the relationship between the peeling force, the peeling velocity, and the geometry of the interface at small-scale. Experiments are compared to a theory based on the linear viscoelastic response of the adhesive, augmented with an elastocapillary regularization approach. This theory, fully quantitative in the limit of small surface deformations, demonstrates the emergence of a "wetting" angle at the contact line and exhibits scaling laws for peeling which are in good agreement with the experimental results. Our findings provide a new strategy for design of reversible adhesives, by quantitatively combining wetting, geometry and dissipation.

Introduction

Pressure sensitive adhesives, ubiquitous for domestic and industrial applications, have the characteristic property that they do not undergo chemical reactions during the bonding process and their performance life. Animals with adhesive pads are ubiquitous in nature 1 and have inspired numerous designs of artificial reversibly adhesive materials [2][3][4][5] . The insights combining viscoelasticity, capillarity, and multiscale hierarchical topography 6 are crucial to design innovative adhesives. Namely, an effective adhesive material should stick under physical contact with a substrate and must therefore respond highly compliant, similar to a liquid. Its adhesive performance results from the resistance to peeling it off a substrate. Strong adhesion implies a high energy dissipation, produced in most polymeric materials by fingering instabilities [7][8][9][10] , cavitation 11 and fibrillar deformation.

From a theoretical perspective, pioneering models have considered the opposite limit of weak adhesion, for which the debonding is interfacial, reversible, and the adhesives remain weakly deformed 12 . According to these theories [13][14][15][16][17][18][19][20][21] , the dissipation dur- ing debonding can be related to the linear viscoelastic properties of the adhesive 17,22 . The recent review by Creton and Ciccotti 2 gives a comprehensive overview of the development of the field. Most experiments with peeling adhesive tapes 23,24 or bulk fracture 25 disagree quantitatively with theoretical predictions 21,26,27 . This leaves a gap in first principles understanding, and has led to the conclusion that non-linear viscoelastic dissipation and, most often, damage mechanisms in the polymer network should be taken into account 2,28 .

In this paper we investigate the dynamics of peeling for a reversible viscoelastic adhesive, which can be peeled off without exhibiting irreversible plastic damage (Fig. 1). Rather than considering the classical case of peeling a thin, strong adhesive with a flexible backing off a rigid solid 7,23,24,29 , we take the opposite perspective: we use a thick layer of a weak adhesive on a rigid backing and peel off a thin flexible tape of a much stiffer material (cf. Fig. 1). This way, we disentangle the effects of bending elasticity, viscoelastic dissipation, and adhesion energy. Our model adhesive is made of a multi-scale polymer gel, whose strong dissipation is controlled by linear viscoelasticity. The key finding is that the dissipation in the bulk is determined by the singular deformation in the vicinity of the contact line, this singularity being regularised by surface energy. This is in stark contrast to "classical" peeling where the blunt and frequently irregular crack front cannot localize dissipation sufficiently and the thickness of the adhesive matters 2 . To resolve all our experimental findings in detail, we propose a theory based on which we establish new scaling laws for peeling of such reversible adhesives, with multiplicative factors quantitatively determined in the limit of small deformations.

Set-up

Elastomers are reticulated polymers obtained by cross-linking long polymer chains. The longer the chains between cross-links (or entanglement points), the larger their relaxation timescale and their effective viscosity. Here we use an addition-cure adhesive that consists of a Silicone gel (Dow Corning CY52-276) whose reticulated polymer network is formed by polymerizing small multifunctional prepolymers. At the gel point, a fractal polymer network forms, which is composed of branches of all lengths between the prepolymer and the size of the sample. We mix the two components at the stoichiometric ratio of 1:1 (A:B), providing extra cross-links with respect to the gel point and leading to a finite shear elastic modulus G 1.2 kPa. The linear rheology is very accurately fitted by the simple relation

µ(ω) = G (ω) + iG (ω) = G [1 + (iτω) n ] , (1)
with an exponent n 0.55. This bulk rheology obey Kramers-Kronig relation: both the storage and the loss moduli originate from the same relaxation function

G 1 + Γ(1 -n) -1 (τ/t) n ,
where Γ is the gamma function. The power law dependence of the loss modulus G ∼ ω n reflects the architecture of the polymer network, with a continuum distribution of relaxation times 30 . The cross-over timescale, measured τ = 0.13 s for our system, is determined by the length of largest branches of polymers in the network. Once stirred and degassed, the mixture of the two prepolymers is poured into petri dishes to make gel layers of varying thickness (from 2 -10 mm). In order to directly test the hypothesis of a dissipation governed by linear viscoelasticity, a ten times stiffer gel has also been used. It was prepared by adding 5% of Sylgard 184 to the aforementioned gel, leading to G 10 kPa, n 0.44, and τ 8.6 ms.

A thin, virtually inextensible film of bending modulus B is placed on the gel, and is subsequently peeled off at controlled forcing. The experiment is performed by placing the system upside-down, inclined at a variable angle φ with respect to the horizontal 31 , and attaching a mass at the end of the tape (Fig. 1). The forcing is varied by four orders of magnitude by using different masses and also by peeling due to the weight of the plastic sheet alone. The sheet's bending modulus B was varied by one order of magnitude (respectively B 9.7×10 -5 J and 6.9×10 -6 J) by using two types of tape: a 88 µm thick biaxially oriented polypropylene film (BOPP manufactured by Innovia) and a 34 µm thick Mylar film (Polyethylene terephthalate) coated with Aluminium (PET, manufactured by Toray). These correspond to Young's moduli of the order of 1 GPa, which is well separated from those of the gels. The tapes are smooth at optical scales -submicrometer roughness could not be measured. Despite no further annealing was applied, the samples did not present any plastic damage nor any spontaneous curvature. The position of the contact line, where the sheet joins the gel, and the geometry of the ridge formed below this contact line are recorded using a video camera (1024 × 1024 pixels with a resolution of 20 µm/pix).

The key control parameter of the experiment is the peeling force f per unit width, defined as the energy released by gravity when the contact line moves by a unit distance. To quantify this force, we focus on the case where the length R of the freely hanging tape is sufficiently large to be quasi vertical close to its free end. Then, peeling the tape by a length dR, the end of the tape moves downward by dR(1sin φ ) [c.f. Fig. 1(b)]. Hence, we obtain the peeling force per unit width:

f = λ g(1 -sin φ ), ( 2
)
where λ is the mass at the end of the tape per unit width. In this expression we neglected the weight of the sheet, which can been included leading to a general expression (see Appendix 7.1).

Geometric characteristics of the peeling front

The present setup provides an experimental access to the geometric features in the vicinity of the contact line 32 . This is of key importance to understand the viscoelastic dissipation during the peeling, which occurs as the highly deformed zone travels along the gel. For the studied parameters, the peeling front remains straight and stable: no oscillation nor stick-slip motion like those reported by Cortet & al. 33 were observed. The peeling dynamics reaches a steady travelling state after a transient that lasts less than a second for the largest velocities. The geometric aspects can be inferred from Fig. 1 (b). Interestingly, the typical scales of deformation, e in the normal direction and parallel to the gel layer, depend non-trivially on f . In Fig. 2 we plot the normal displacement e as a function of forcing f , at different inclination angles φ . The response is far from linear, with a scaling law e ∼ f 1/2 , even though typical strains are within the linear range of the gel.

This nonlinear response of the normal displacement can be understood in two steps. First, we analyze the part of the elastic gel that is in contact with the tape. The gel layer is soft and thick, ) and the symbols to the type of experiment, as shown in Fig. 6. The marker's size codes for gel thickness. The data collapses according to (6) shown as the solid line (prefactor 0.26).The circled cross symbol corresponds to a series of data obtained with a 7.5 mm thick layer of the most rigid gel, and with the flexible backing in the horizontal case φ 0 • . Inset: the shape of the sheet is parametrised by θ , the shape of the gel by h.

and characterized by its (static) shear modulus G. A surface deflection of amplitude e and horizontal scale induces a normal stress σ ∼ Ge/ . By contrast, the elasticity of the thin, stiff tape is characterized by its bending modulus B, and gives a typical normal stress σ ∼ Be/ 4 . The balance of stress thus does not select the normal deflection e, but provides access to a lateral length scale given by:

= B G 1/3 . (3)
This elasto-bending length is the wavelength of wrinkles that appear when compressing a soft foundation that is covered by a hard, thin skin 34,35 . In the present case of a peeling experiment, where the elastic film is under tension, represents the decay length over which the angle of the tape aligns to the gel layer (Fig. 1b) ( 1.8 mm for the flexible tape and 4.3 mm for the stiff one). To check this prediction, we have directly measured the length exp over which the tape relaxes to its asymptotic angle. The measurements, reported in Fig. 3 shows that it presents subdominant variations when the force is varied over more than three orders of magnitude. Moreover, it is as expected on the order of the elasto-bending length .

As the length exp is only rigorously defined in the limit of small deformations for which a linear relaxation is expected, we have also measured the inclination of the tape at the contact line θ with respect to the horizontal, for different peeling forces (see Fig. 1b for the definition of θ ). A very simple view of the deformation would be to approximate the tape as a triangle, with normal extension e and lateral extension on the side of the gel given by . This would give a geometric relation tan θtan φ = (e/ )(1 + tan θ tan φ ). However, this approximation has the drawback that it does not capture the large deformation asymptotic: when e → +∞, θ → π/2 since the tape is aligned with gravity already at the location of the contact line. Therefore we propose the fit:

tan θ -tan φ = e α , (4)
as it correctly gives the divergence of the normal deformation e when θ → π/2 and is consistent with the triangular model for The inset shows the case where the gel is horizontal (φ = 0). The data are very well described by the fit (4) for all φ , shown as a solid lines. The prefactor α was found of order unity for all cases (see caption), and shows that as defined by (3) indeed sets the lateral scale -with a weak dependence on the inclination angle φ . The result obtained for φ = 0 is shown in the inset of Fig. 4 and is accurately described by tan θ ∼ e/ , for different f and B, confirming that (3) correctly predicts the lateral scale .

φ = 0.
To explain the nonlinear behaviour of the normal displacement e, we now turn to the freely hanging part of the sheet that is not in contact with the gel. The shape is that of a classical elastica, forced by the mass at the end of the tape. Here we parametrise the shape of the tape by θ (S) relating its local angle to the curvilinear coordinate S (Fig. 2, inset). Introducing the tangential unit vector t = (cos θ , sin θ ), the elastica equation can be integrated to

1 2 Bθ (S) 2 + λ g • t(S) -t(R) = 0, (5)
where we used that the end of the tape at S = R is free from torque. At the contact line we can estimate the bending term from the characteristic horizontal and vertical scales:

Bθ 2 ∼ Be 2 / 4 .
The forcing term in ( 5) can be written as λ g(1sin φ ) = f , which becomes exact when the hanging part of the sheet is long. Combined with (3), this gives the nonlinear scaling for the normal

This relation is successfully tested in Fig. 2, where we find a collapse over 4 decades, for different inclinations φ , for two different bending moduli B, for two different gel elastic moduli G, for different gel thicknesses and for different adhesion conditions.

Reversibility, work of adhesion and contact angle

The reversibility of the peeling process is illustrated by a series of experiments without an additional mass at the end of the tape. When the hanging part of the sheet is sufficiently long, peeling can in fact be induced by the weight of the sheet. However, there is a threshold length below which the surface energy due to adhesion is stronger than gravity, so that the sheet spontaneously reattaches to the gel when its end is released. The resulting peeling velocities v are presented in Fig. 5a, where positive v correspond to peeling and negative v to reattachment. At the threshold point where v = 0, the system is at equilibrium: there is an exact balance between the forcing f and the (conservative) work of adhesion Γ. The value for the softer gel and the metallized sheet is Γ = 52 ± 3 mN.m -1 and for the softer gel and the Mylar sheet is Γ = 19 ± 3 mN.m -1 . For the ten times more rigid gel and the metallized sheet Γ = 10 ± 5 mN.m -1 . Hence, this experiment allows us to accurately determine the non-dissipative contribution to the adhesion process. Note that reversibility is, here, not intended in the thermodynamic sense (there is a visco-elastic dissipation) but in the structural sense: the adhesive does not present irreversible damages after peeling. The work of adhesion can be directly related to the shape of the gel on the side that is detached from the sheet. Figure 5b shows that the free surface of the gel is highly curved before contacting the sheet. Just like in a crack à la Griffith 36 or in the JKR adhesion problem 37 , the normal deflection follows a square root shape h ∼ x 1/2 , where x is the horizontal distance measured from the contact line 37 (see Appendix 7.2). This square root singularity is the solution of the mixed problem, where a contacting region changes to a stress-free interface. It must not be confused with the logarithmic deformation of a surface submitted to a localized force, which appears here (or in the JKR problem) as the far field solution far away from the contact line 34 . Remarkably, however, the gel is clearly seen to make a well-defined contact angle θ Y when touching the sheet (Fig. 5b). This strongly resembles the wetting of liquids, for which the work of adhesion can be expressed as

Γ = γ(1 + cos θ Y ), ( 7
)
where γ is the surface energy of the liquid-vapor interface and θ Y is Young's contact angle. The gel could indeed obey similar wetting laws [38][39][40] below the elastocapillary length

L = γ G , (8)
where in this case γ is the surface energy of the gel. This length, approximately 30 µm for the softer gel and ten times lower for the more rigid, indicates the scale where surface energy dominates over bulk elasticity, and below which we expect "wetting" behavior.

Despite subtleties of capillarity of elastic interfaces 41? -43 , this wetting interpretation is indeed consistent with our direct measurement of the gel's contact angle. For the softer gel, we measured θ Y = 65 ± 5 • with the metallized sheet, and a higher contact 6 Relation between the peeling force f and the peeling velocity v, scaled according to (10). The colors correspond to different inclination angles φ , from red to blue (-45 • , 0 • , 45 • ) and the symbols to the type of experiment, as shown in the table in inset. The symbol size codes for gel thickness. The inset represents the data without rescaling. The circled cross symbol corresponds to the 7.5 mm thick layer of the most rigid gel and the flexible backing, for φ 0 • , with a mass.

angle θ Y = 125 ± 5 • with the less adhesive Mylar sheet. Using the measured Γ and the previously reported value for the surface tension of the gel-vapor interface, γ = 39 mN.m -1 44 , the aforementioned relation predicts contact angles of respectively θ Y = 70 • and θ Y = 120 • , consistent with direct optical estimates.

Energy release and dissipation

We now turn to the most important characteristic of the adhesive, namely the relation between the peeling force and the peeling velocity. The unscaled experimental data is shown in the inset of Fig. 6, where we present the mass λ versus the velocity v. For each dataset with given tape and inclination φ , we find a power law with an exponent 0.53 ± 0.04 for the softer gel and an exponent 0.37 ± 0.15 for the ten times more rigid gel. Below we demonstrate that this directly reflects the exponent n of the rheology, as given by (1). Another observation is that the peeling velocity is independent of the gel thickness (represented by the size of the symbols), from which we deduce that dissipation is localized in the vicinity of the peeling front.

To interpret the relation between forcing and peeling velocity, we make use of an energy balance. Since inertial effects are negligible, the viscoelastic dissipation in the bulk of the deformed gel must equal the changes in surface energy (due to the work of adhesion) and the forcing by the weight. We thus need to determine the dissipation inside the bulk of the gel layer due to the rheology (1). A closed form expression for the dissipation can be obtained in the limit of small deformations, assuming that the gel profile (cf. inset Fig. 2) is a traveling wave h(xvt). The resulting balance reads (see Appendix 7.3):

f -Γ = dq 2π qG (qv) k(q) | ĥ(q)| 2 . ( 9
)
The integral on the right hand side represents the dissipation, where ĥ(q) is the Fourier transform of the gel deformation and k(q) is the spatial Green's function relating deformation to the normal stress. The finite thickness H of the gel layer enters in k(q) for qH 1. In our case deformation is limited to H and we may use the half-space approximation k(q) ∼ (2|q|) -1 ; in agreement with the experimentally observed independence of thickness. Naturally, ω = qv sets the characteristic frequency for the dissipation, as can be seen from the argument of G in (9). Furthermore, the independent calibration (1) allows us to write G ∼ (qvτ) n G.

The final step is to insert the deformation profile ĥ(q) in ( 9). Importantly, when using the crack-shape h ∼ ( x) 1/2 , for which ĥ ∼ 1/2 /|q| 3/2 , the dissipation integral in (9) diverges at large q and gives infinite dissipation at small scales. However, the appearance of a wetting condition at scales below ∼ L provides a cutoff: dissipation becomes integrable when the interface exhibits a finite angle, with details depending on θ Y . Using q ∼ L -1 as the regularisation scale, ( 9) becomes:

f -Γ = β τv L n G . (10)
Figure 6 confirms that all data are well-described by (10) with β a multiplicative factor encoding the boundary condition effect.

The precise value of β primarily depends on the geometry of the gel at the peeling front: it is found to be ∼ 10 times larger for θ Y 65 • than for θ Y 125 • . This illustrates the importance of the contact angle, also for dynamics, since a larger peeling ridge leads to stronger dissipation. β exhibits a subdominant non-monotonic dependence on φ that can be traced back to a weak variation of the lateral relaxation scale (see the value of α in Fig. 4).

β φ 0 • φ -45 • φ 45 • θ Y = 65 • ± 5 • 3.3 1.2-1.4 2.2 θ Y = 125 • ± 5 • 0.43 0.07 θ Y = 160 • ± 15 • 0.25

Discussion

While first theories of peeling of pressure sensitive adhesives have described purely static equilibrium situations 12,45 , it was soon realised that the peel force was velocity dependent and associated to the rheology of the adhesive [46][START_REF] Gent | Proc. R. Soc[END_REF][START_REF] Kendall | [END_REF][49] . The problem of growing cracks in bulk viscoelastic media is closely related, and Schapery 26 predicted a power law dependence of the growth rate of the crack on the applied gross strain, where the exponent is given by the exponent of the creep compliance function. Also closely related is the geometry of a rigid cylinder rolling on a viscoaleastic material, for which a strain energy release rate ∼ v 0.55 has been found 50,51 . Maugis and Barquins 52 showed that peeling of urethane strips off a glass surface yields a power law relation between peeling speed and applied force as well. They attributed this behavior to the rheology of the adhesive since the empiric dissipation relation obeyed the same time-temperature superposition behavior as the loss modulus. They also noticed that, if viscoelastic losses are localized at the crack tip, dissipation was supposed to be geometry independent.

However, to the best of our knowledge, a quantitative theory capable of predicting the peel force for a highly localized dissipation was still missing in literature. Historically, this is most likely due to the fact that typical properties of engineered adhesive tapes lead to characteristic length scales which do not admit such strong localisation: Such tapes typically comprise a thin layer of a soft and highly adherent material, so that the stress localisation ∼ Γ/G ∼ O(1 mm) is spread much wider than the typical layer thickness ∼ O(10 µm) 2,28 . Instead, a cohesive zone forms, composed of fibrils and cavitation bubbles, and the strain energy release rate becomes geometry dependent and typically increases with layer thickness 53 . Those cases clearly require a different type of model, taking the nonlinear rheological properties of the adhesive into account 2,7,10 .

The growing interest in reversible adhesives that can be peeled without bulk cavitation or plastic deformations, as inspired e.g. by biology 1,[START_REF] Ghatak | Proceedings of the Royal Society of London A: Mathematical[END_REF] , suggests revisiting the limit of weakly adherent materials and peeling geometries with localised stress and dissipation. However, the crack singularity of adhesion leads in this case not just to a concentration of the stress to the peeling front, but also to a diverging dissipation in the continuum description. This prevents quantitative predictions unless a physical regularisation mechanism is identified. While the power laws that we find in our experiments are similar to those reported previously 51,52 , we could disentangle the various contributions to the peel force by peeling an inextensible flexible tape with different surface energies off a thick elastic layer [START_REF] Ghatak | [END_REF] . Most importantly we have identified solid surface tension and the corresponding "wetting type" boundary condition as the leading order regularisation mechanism for dissipation in reversible peeling.

In conclusion, we have shown that reversible adhesives can indeed obey simple scaling laws for deformation and dissipation, whose origin can be traced back to linear viscoelasticity and, importantly, solid capillarity. Analysing the near-crack-tip geometry, we have shown that the regularisation of elastic singularities by the wetting condition allows one to get quantitative estimates of the viscoelastic dissipation. The theoretical framework proposed here bridges the gap between adhesive peeling and moving contact lines of fluids, opening the possibility of fully quantitative theories. It opens the promising perspective of designing adhesives by coherently tuning their visco-elastic properties, their surface functionalisation, and their meso-and macroscopic architecture.

Appendix: derivation of the model equations

In this appendix section, we derive the model used along the article as a framework of interpretation of experimental results. We provide at each step the general, rigorous equations, and their expansion at the lowest order in deformation and contact line velocity.

Peeling force and boundary condition at the contact line

Here we will derive the expression of the peeling force f , for the experiments with or without and additional mass attached at the end of the tape. The length R of the tape is arbitrary. The driving force of the peeling motion is then due to the weight of the tape, the additional mass, and the bending moment of the tape.

To derive the governing equations for the shape of the free part of the tape, we use curvilinear coordinate S and the corresponding tangent vector t(S) = cos θ (S) e x + sin θ (S) e y along the contour of the tape; θ (S) denotes the angle between the tangent vector and the horizontal. We define that the contact line position corresponds to the curvilinear coordinate S = 0, and the end of the tape is located at S = R. The torque balance on the piece of tape going from S to R reads:

Bθ (S) e z + ρ S g × R S (R -S ) t(S )dS + λ g × R S t(S )dS = 0 (11)
The first term, which involves the bending modulus B of the tape, is the elastic torque. The second term is the gravity torque due to the weight of the tape itself. ρ s is the mass density of the tape per unit area and g = g(sin φ e x + cos φ e y ) is the gravitational vector. The last term is the gravity torque due to the additional mass at the end of the tape, λ , per unit width of the tape.

We define the peeling force f as the energy released by gravity and bending; here we ignore -for now -the surface tensions of tape and gel. We therefore write the free energy F d of the free part of the tape and calculate its variations with respect to shape and contact line motion. The bending energy and gravitational Fig. 8 Predictions of the linear visco-elastic model (equations ( 14)-( 16)) at vanishing velocity v. Typical equilibrium profiles of the slope h on the free side of the contact line as a function of the coordinate x rescaled by for three different contact angles. At small x L , h saturates to the value selected by the balance of surface tensions. For L x , an intermediate asymptote h ∼ x -1/2 emerges. Therefore, h ∼ x 1/2 , reminiscent of the crack opening problem. Finally, at large x , one recovers an outer asymptote h ∼ x -1 or equivalently h ∼ log x, which is the elastic response to a localised line load. Inset: complete solution h(x) showing the relaxation over a scale on the tape side and the intermediate asymptotics featuring a square root behaviour on the free side.

energy of the detached part of the sheet are expressed as:

F d = R 0 dS 1 2 Bθ (S) 2 + 1 2 ρ S R 2 g sin φ - R 0 dSρ S g • S 0 t(S )dS + λ Rg sin φ -λ g • R 0 dS t(S) (12)
The first line of the right hand side is the bending energy, the second line describes the gravitational energy due to the weight of the tape itself and the third line is the gravitational energy of the additional mass at the end of the tape. Note that both gravitational terms result from two contributions each: not only the total hanging contour length R of the tape changes (second term), but also the contact line position moves simultaneously (first term). Here we ignore that the gel deformation also slightly depends on R, which causes an additional motion of the contact line relative to the laboratory frame. Now, the total variation of the free energy with respect to a change of length R gives the peel force f :

f = - ∂ F d ∂ R (13) = -g(ρ S R + λ ) sin φ + λ g • t(R) + ρ S g • R 0 dS t(S)
To obtain f from the experiment, all terms on the right are evaluated directly from the images. We have checked experimentally that the variation of the angle θ at the contact line with respect to R leads to a negligible contribution to f .

Equations coupling the tape to the gel mechanical response

In order to derive the equations governing the shape of the gel surface, we proceed to a double expansion: we consider the rheology as being linear, which is well satisfied by the gel used here, even under stretching; furthermore, we linearize the equations at small deformations of the gel surface, in order to be able to use the Green's function formalism. Although this linearity requirement is badly obeyed in many of our experiments, we argue here that the laws found using this approximation may still be valid, within a multiplicative constant that cannot be calculated quantitatively. The scaling behaviour of our experimental results show that this is indeed the case. We furthermore neglect gravity in the description of the gel.

We first introduce the elastic stress σ at the surface of the gel. On the free surface side of the gel, the shear stress vanishes, and the normal stress is balanced by the solid Laplace pressure ϒ S h . The solid surface tension ϒ S is considered constant, for simplicity (ϒ S = γ = L G). On the tape side, the strain ε along the surface vanishes on the tape side. This implies that, at linear order, the shear stress vanishes, but not the normal stress. The latter is set by the bending of the tape, Bh . Thus,

σ = ϒ S h for x < 0 -Bh for x > 0 (14)
There are four boundary conditions at the contact line: the selection of Young's angle (Eq. 7), the selection of the tape angle θ (0), and the second and third derivatives inherited from the free part of the tape.

Next we relate the normal traction σ to the normal displacement and to the viscoelastic rheology µ(ω) = G (ω) + iG (ω). When taking the Fourier transform in both space (x → q) and time (t → ω), the kernel relating the deformation to the normal stress can be written as the ratio µ(ω)/k(q) 44,56 , where k(q) is the spatial Green's function corresponding to the linear response. Here we implicitly made use that for an incompressible thick layer, the normal stress is decoupled from the tangential displacement [START_REF] Johnson | Contact mechanics[END_REF] . We consider travelling wave solutions which depend on the single variable xvt. The Fourier transforms of the elastic stress σ and the gel profile ĥ are then related by:

σ (q) = µ(qv) k(q) ĥ(q) (15)
The effect of velocity and surface tension are therefore only to affect the space Green's function. At small velocity v, the rheology can be expanded at the lowest order, leading to the static shape relation. In this case and for a thick layer ( gel thickness), the profile slope and the stress are related to each other by a Hilbert transform:

σ (x) = - 2G π ∞ -∞ h (x ) x -x dx (16)
Considering previous work for liquid droplets on the same gel 44 , the range of velocities for which this is valid extends up to 2 • 10 -3 m/s. The above set of equations is complemented by the asymptotic condition h ∼ -x -1 for |x| , which corresponds to load that is localized on |x| .

Figure 8 shows the typical surface profiles obtained by numerical resolution of the visco-elastic model in the limit of small velocities. Far from the contact line, the shape becomes logarithmic, as expected from the response to a localized force 34 . At intermediate distances to the contact line, the normal deflection follows the square root shape h ∼ x 1/2 , as is expected in such a mixed problem where a compliant contact zone connects to a stress-free surface 36,37 . However, due to solid capillarity, the square root singularity is regularized close to the contact line, leading to a well-defined Young's contact angle. The square root shape is retained as an intermediate asymptote.

The dissipation integral

Here we compute the rate energy dissipation P inside the gel during peeling, which is used to derive equation ( 9) in the main manuscript. We define the displacement field inside the gel u i and the stress tensor σ i j . With this, we compute the energy dissipation per unit time (per unit width of the gel) from the usual integral 56

P = d 2 x σ i j ∂ ui ∂ x j = d 2 x ∂ ∂ x j σ i j ui -u j ∂ σ i j ∂ x j = ds σ i j n j ui . (17)
These contain standard manipulations, where we used mechanical equilibrium ∂ σ i j /∂ x j = 0 and brought the area integral to the boundary with normal vector n j . The boundary integral represents the work done by normal and tangential tractions. Since the bottom of the gel is fixed to a rigid support ( ui = 0) the only contribution comes from the free surface. In addition, only the normal traction performs work: the inextensibility of the sheet imposes vanishing tangential displacement, while the tangential stress vanishes on the side that is peeled from the sheet. In the limit of small deformation, the normal displacement can be identified with the gel's profile h(x,t), so that (17) reduces to

P = ∞ -∞ dx σ (x,t) ḣ(x,t), (18)
where we denote the normal stress σ = σ i j n i n j .

Inserting in (18), the dissipation becomes:

P = dx dq 2π
µ(qv) k(q) h c (q)e iq(x-vt) × dq 2π (iq v)h c (q )e iq (x-vt) (

We can introduce the variable x = xvt and perform the x integral using d x e i(q+q ) x = 2πδ (q + q ). (

) 20
This gives a dissipation P = dq 2π µ(qv) k(q) h c (q) dq 2π (iq v)h c (q ) 2πδ (q + q ) (21)

= -v dq 2π µ(qv) k(q) h c (q)(iq)h c (-q) (22)

= -v dq 2π µ(qv) k(q) (iq)|h c (q)| 2 . ( 23
)
This expression is real, owing to the symmetry properties k(q) = k(-q) and µ(-ω) = µ(ω) * , which implies G (ω) = G (-ω) and G (ω) = -G (-ω). Hence, as expected, only G contributes to the dissipation. Finally, we can define the dissipative force P = f d v. This gives the expression used in equation ( 9)

f d = ∞ -∞ dq 2π qG (qv) k(q) |h c (q)| 2 . ( 24

a

  Laboratoire de Physique Statistique (LPS), UMR 8550 CNRS, ENS, Univ. Paris Diderot, Sorbonne Université, 24 rue Lhomond, 75005, Paris, France. b Laboratoire de Physique et Mécanique des Milieux Hétérogènes (PMMH), UMR 7636 CNRS, ESPCI, Univ. Paris Diderot, Sorbonne Université, 10 rue Vauquelin, 75005 Paris, France c Max Planck Institute for Dynamics and Self-Organization, Am Fassberg 17, 37077 Goettingen, Germany d Physics of Fluids Group, Faculty of Science and Technology, University of Twente, P.O. Box 217, 7500 AE Enschede, The Netherlands

Fig. 1

 1 Fig. 1 (a) Image of the peeling experiment. A tape is pulled vertically from the adhesive gel under the influence of a weight (not in the image) hanging at its end. (b) Schematic illustrating the gravitational energy released by a peeling a length dR. The forcing is controlled by varying the attached mass and the inclination φ . The image also defines the normal deflection e, the lateral scale = (B/G) 1/3 on the side where the tape is still attached, and the tape inclination at the contact line θ .

Fig. 2

 2 Fig.2Relation between the gel deformation e, rescaled by the elastobending length = (B/G) 1/3 , and the peeling force f rescaled by and by the gel elastic modulus G. The colors correspond to different inclination angles φ from red to blue (-45 • , 0 • , 45 • ) and the symbols to the type of experiment, as shown in Fig.6. The marker's size codes for gel thickness. The data collapses according to(6) shown as the solid line (prefactor 0.26).The circled cross symbol corresponds to a series of data obtained with a 7.5 mm thick layer of the most rigid gel, and with the flexible backing in the horizontal case φ 0 • . Inset: the shape of the sheet is parametrised by θ , the shape of the gel by h.

Fig. 3

 3 Fig. 3 Measured lateral decay length exp on the side of the tape divided by the elasto-bending length for different forces. The symbols are the same as on the other figures.

Figure 4

 4 presents the angle θ as a function of the deformation e/ for the different sets of experiments -type of tape, configuration (with or without a mass) and global inclination φ .

Fig. 4

 4 Fig. 4 Angle θ at the contact line as a function of the deformation e, rescaled by the elastic length . The data correspond to different types of tape and configurations. The colors correspond to the global inclination angle φ . The solid lines show the best fit by the form tan θ = tan φ +e/(α ). For red points φ = -45 • and α = 0.51, for green points φ = 0 • and α = 0.6, for blue points φ = 45 • and α = 0.18. Inset: slope tan θ of the sheet at the contact line for φ = 0 • . The solid line indicates tan θ ∼ e/ .

Fig. 5

 5 Fig. 5 (a) Relation between the peeling velocity |v| and the peeling force f on semi-logarithmic axes, for the case where no mass is added to the sheet. Negative velocities, for which the tape reattaches, are on the left, for forces below Γ. Data correspond to the metallized film and the softer gel. (b) Side view of the gel in the equilibrium condition f = Γ, showing Young's contact angle θ Y 65 ± 5 • .

  Fig.6 Relation between the peeling force f and the peeling velocity v, scaled according to(10). The colors correspond to different inclination angles φ , from red to blue (-45 • , 0 • , 45 • ) and the symbols to the type of experiment, as shown in the table in inset. The symbol size codes for gel thickness. The inset represents the data without rescaling. The circled cross symbol corresponds to the 7.5 mm thick layer of the most rigid gel and the flexible backing, for φ 0 • , with a mass.

Fig. 7

 7 Fig.7 Theoretical schematics. (a) For the derivation of the peeling force f we only consider the hanging tape (shaded on the scheme) for which the contour coordinate S > 0. (b) For the derivation of gel-tape equilibrium shape we consider the attached part of the tape and the elastocapillary mechanics of the gel. The shaded part of the attached tape represent the system considered for the derivation from S = 0 to S < 0.

Table 1

 1 Table of the parameter β defined by the force balance (10) according to the various experimental configurations: type of tape and inclination φ .

)

  Acknowledgments. B.A. is supported by Institut Universitaire de France. This work was funded by the ANR grants Smart. JHS acknowledges financial support from ERC (the European Research Council) Consolidator Grant No. 616918. SK acknowledges financial support from the Max Planck -University of Twente Center 'Complex Fluid Dynamics -Fluid Dynamics of Complexity'.

1-9


Maud Gits-Muselli

P Lewis White

Carlo Mengoli

Sharon Chen

Brendan Crowley

Gijs Dingemans

Emilie Fréalle

Rebecca Gorton

Malcom Guiver

Ferry Hagen

P Lewis White

Gijs

Dingemans

Halliday

Gemma Johnson

Katrien Lagrou

Martina Lengerova

Willem Jg Melchers

Lily Novak-Frazer

Riina Rautemaa-Richardson

Emeline Scherer

Joerg Steinmann

Mario Cruciani

Rosemary Barnes

J Peter Donnelly

Juergen Loeffler

Dr Alexandre Alanio
email: alexandre.alanio@pasteur.fr

The Fungal PCR Initiative's evaluation of in-house and commercial Pneumocystis jirovecii qPCR assays: towards a standard for a diagnostics assay

Keywords: Pneumocystis jirovecii, Pneumocystosis, PCP, Diagnosis, qPCR, quantification, DNA, Whole 60 nucleic acids, Panel specimens, standardization, efficiency, quantification cycles, Cq, 61 threshold

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

To calculate the qPCR efficiency of each assay, calibration curves were generated using 210 the highest positive BALF extract (BALF9: initial Cq=17 with the mtSSU in-house RT-qPCR 211 assay). After thawing of an aliquot, a serial dilution of the extract (6 times at 1:5 ratio in 212 negative pooled extract) was tested for each PCR assay. values. The ANOVA analysis was implemented as a mixed effects regression model, where Cq 240 was the dependent variable, the molecular target level was used as categorical regressor, and 241 the covariate under evaluation was added as an additional factor. The interaction between 242 the covariate and the relative genomic load was also included in the model. Anonymization 243 was performed prior to all analysis and identifiers were used to build the categorical grouping.

244

The identifier indicated the mixed model categorical grouping variable. As a post-estimation 245 step, the contrasts (predicted mean effects) were calculated for each pair of factor levels, 246 adjusting for each relative genomic load level. They measured the difference between the 247 effects exerted by each member of the factor pair on the signal intensity.

248

For the monocentre study, no specific statistical analysis was applied. The maximum variation in the Cq value (excluding negative values) between assays 276 was 12.4 cycles (approximately 10,000-fold variation), 9.2 and 9.4 cycles (approximately 1000-277 fold variation) for the 1:1, 1:100 and 1:1000 panel specimens, respectively (Supplemental 278 Figure 2). Of note, mtLSU assays using the same primers and probes, 26 but with different 279 reagents and/or protocol, gave a median difference of 3 cycles (i.e. 10-fold variation) 280 (Supplemental Table 5).

F

281

Two major variables impacted significantly on the results, in addition to minor 282 influences associated with variation in reagents and protocols: (i) the type of qPCR, RT-qPCR 283 vs qPCR, and (ii) the gene targeted for amplification (Figure 1). The RT-qPCR assays based on 284 the detection of WNA (detecting both RNA and DNA) were significantly superior to qPCR 285 assays detecting DNA only at all loads (p<0.001, Supplemental Table 6 andSupplemental 286 Figure 3). The Cq values were also significantly associated to the target gene at all loads with 287 the lowest Cq values associated with mtSSU, successively followed by mtLSU, MSG, and then 288 TUB (Figure 1, Supplemental Table 7). There was a consistent and significant correlation The specificity of the 10 assays was high; nine gave a specificity of 100% (upon duplicate 297 testing of the 10 negative BALFs) and only one assay (assay#12) provided one false-positive 298 result (one of the duplicates) giving a specificity of 95%.

299

The sensitivity of the assays for the detection of P. jirovecii in 9 positive BALFs varied 300 between 88.8% (n=5 assays) and 100% (n=2 assays) with 4 assays having a sensitivity of 94.4% 301 (Supplemental Table 8 and Figure 2B). The lowest Cq values were observed with the mtSSU 302 assays amplifying WNA and the highest with the MSG assays amplifying DNA across all positive 303 BALFs (Figure 2B). Thus, when excluding the negative values (Cq=45), the median of the 304 maximum variation of the Cq value between the 10 assays was 8.23 (approximately 500-fold 305 variation) for a given BALF. Even when considering the 5 commercial assays targeting mtLSU 306 as the target, a median of variation as high as 3.1 (about 10-fold variation) was observed 307 (Supplemental Table 8 and Figure 2B). 308 By serially diluting BALF#9 WNA extract five times at 1:5 dilutions, the qPCR efficiency 309 for each assay was calculated. The amplification efficiency calculated from the Cq values were 310 all between 94 and 113% (Supplemental Table 8, Supplemental Figure 4). This showed that 311 this factor can be discarded as a possible explanation for the differences observed above. The present study underlines the necessity of comparing different assays for the 314 quantification of the P. jirovecii loads. Indeed, large variation was evidenced among all the 315 assays compared. Based on data of 9 commercial kits and 10 in house assays tested in the 316 multicentre and monocentre studies (Table 1), this variation ranged from 500-fold (5 qPCR 317 cycles) to 10,000-fold (12 qPCR cycles) between the assays. 318 This resulted mainly from the target gene (mitochondrial multicopy genes vs. single copy 319 nuclear gene) and the qPCR method/template (qPCR/DNA versus RT-qPCR/WNA). Thus, our 320 studies showed that using mtSSU as the target gene and WNA as the starting material yielded 321 the best sensitivity. Predictably, assays targeting multi copy genes provided the best analytical 322 sensitivity, as already shown. 16,19,20 This was observed in our study for the single copy gene 323 (beta-tubulin) compared to the mtSSU gene (about 200-fold variation). However, among the 324 repeated genes, MSG assay gave higher/later Cq values than mtLSU or mtSSU, suggesting a 325 lower number of repeats for the MSG gene. The MSG gene sequence is reported to be highly 326 variable, 27 which can impacts both the analytical sensitivity/limit of detection (risk of false 327 negative) and qPCR efficiency (risk of variable amplification yield). 328 Among the mitochondrial multicopy assays, the mtSSU assay was two times more 329 sensitive on average than the mtLSU assay, possibly as a result of the higher number of repeats 330 of mtSSU compared to mtLSU in the mitochondrial genome. 20 Interestingly, the number of 331 copies of some mitochondrial genes (mtLSU, NAD1 and CYTB) have been shown to vary 332 according to the fungal load (high vs. medium and low fungal loads). 20 In contrast, mtSSU 333 seems to remain stable between these conditions suggesting that a more reliable 334 quantification could be achieved with mtSSU rather than with other mitochondrial targets. 335 Indeed, the quantification of unique/single copy gene (nuclear) will better reflect the quantity "carriage/colonization" and therefore of limited clinical significance. However, if the negative 361 predictive value of the assay to exclude PCP is taken into consideration, it seems preferable 362 to rely on optimal sensitivity to exclude PCP, rather than managing false negative results. 363 Indeed, low P. jirovecii burden can be associated with severe disease in HIV-negative 364 patients. [START_REF] Robin | Molecular demonstration of a Pneumocystis 568 outbreak in stem cell transplant patients: evidence for transmission in the daycare centre[END_REF] Second, prophylaxis could be considered for immunocompromised patients with a 365 low fungal load to avoid subsequent development of full-blown PCP. [START_REF] Mori | A followup study of asymptomatic carriers of Pneumocystis 571 jirovecii during immunosuppressive therapy for rheumatoid arthritis[END_REF] At least, knowing that 366 P. jirovecii is detectable can tilt the balance between continuing co-trimoxazole or stopping 367 the drug because of the fear of side effects. Third, patients with low fungal loads in an hospital 368 environment can be part of the onward transmission chain in clinical outbreaks. 3,[START_REF] Robin | Molecular demonstration of a Pneumocystis 568 outbreak in stem cell transplant patients: evidence for transmission in the daycare centre[END_REF]

369

There are limitations to this study. To set up a quality control with P. jirovecii, one must more that 70% of the cases, [START_REF] Alanio | Diversity of 574 Pneumocystis jirovecii during infection revealed by Ultra-Deep Pyrosequencing[END_REF][START_REF] Gits-Muselli | New short tandem repeat-based molecular 577 typing method for Pneumocystis jirovecii reveals intrahospital transmission between 578 patients from different wards[END_REF][START_REF] Hauser | Genetic diversity of Pneumocystis carinii in HIV-581 positive and -negative patients as revealed by PCR-SSCP typing[END_REF][START_REF] Parobek | Multilocus microsatellite genotyping array for 583 investigation of genetic epidemiology of Pneumocystis jirovecii[END_REF] with known sequence variation in targets genes. [START_REF] Alanio | Diversity of 574 Pneumocystis jirovecii during infection revealed by Ultra-Deep Pyrosequencing[END_REF] 382 Indeed, false negative results can result from the choice of primers in variable regions of the 383 fungus. [START_REF] Gal | A misleading false-negative result of 586 Pneumocystis real-time PCR assay due to a rare punctual mutation: A French multicentre 587 study[END_REF] This is therefore the responsibility of each centre to investigate this when facing a The absolute quantity of Pneumocystis DNA added in the initial pooled of samples was not 388 known. Nevertheless, the main purpose of our study was to compare different assays and not 389 to determine the absolute analytic capacity of each assay. Another possible bias in pooling 390 human samples is to dilute PCR inhibitors present in one of the samples for elaborating the 391 panel. Using an internal control is of outmost importance to analyse qPCR results in a 392 diagnostic lab and is strongly advocated [START_REF] Bustin | The MIQE guidelines: minimum information for 589 publication of quantitative real-time PCR experiments[END_REF][START_REF] Hoorfar | Making Internal Amplification Control Mandatory for 592 Diagnostic PCR[END_REF] and improves the quality of the results, especially 393 to have a strong validation for negative results. [START_REF] White | Evaluation of Aspergillus PCR protocols for testing 595 serum specimens[END_REF] Thus, some PCR assays not controlling for 394 the presence of PCR inhibitors can be falsely over-evaluated in the present study. With limited 395 availability of Pneumocystis culture the current methodology based on surplus human 396 samples is the only one possible to perform such evaluations, despite all the inherent potential 397 bias. However, we do not think the possible limitations listed above had a major impact on 398 our comparative studies.

399

The limited quantity of BALF extract available restricted the number of evaluations 400 possible. Thus, some assays were tested only once in a single centre (DNA MSG, DNA TUB, 401 WNA mtSSU, DNA mtSSU) making any evaluation on the reproducibility for these specific 402 assays impossible. Additionally, only the analytical sensitivity of the amplification was tested, 403 not the nucleic acids extraction step, which has a major impact for the global performance of 404 all the amplification procedure. [START_REF] White | ) BALF5 (Cq1/Cq2) BALF6 (Cq1/Cq2) BALF7 (Cq1/Cq2) BALF8 (Cq1/Cq2) BALF9 (Cq1/Cq2) PCR Efficiency (%) 2[END_REF] We also did not have the opportunity to evaluate the impact

  comparison, we selected 10 P. jirovecii-negative and 9 P. jirovecii -204 positive BALFs. The 50 µl of the BALF extracts left after routine procedure were combined with 205 70 µl of negative extract (see above), to provide sufficient material for the comparison of the 206 kits. Verification was performed using the mtSSU in house RT-qPCR assay which generated 207 initial Cq values ranging from 17 to 30 cycles for the 9 qPCR-positive BALF extracts . All aliquots 208 were stored at -20°C until use.209

  10 different assays based on the in-house assays available locally 215 and the commercial assays available in his country with the same person and the same qPCR 216 instrument as above. For Centre E in-house assays, the same protocols as above were used. 217 Commercial qPCR assays were performed following the manufacturer's instructions. The 218 major surface glycoprotein (MSG) in-house DNA assay was performed following instruction of 219 the detailed process described in Larsen et al. 25 220 Technical information on the 10 assays, consisting of 5 in-house assays and 5 221 commercial kits, are listed in

  study, 20 datasets of the five samples tested by each assay in duplicate 228 were produced. Therefore, the total number of data values included was 200. Each replicate 229 was analysed as an independent result. Negative results were allocated a value of 45 cycles. 230 Descriptive statistics were performed calculating summary indexes of Cq at all target loads, 231 pure (1:1) and 1:100, 1:1000 dilutions. 232 An ANOVA model was used to evaluate the effect of a two main covariables: DNA vs. 233 WNA, and target gene. The effects of these categorical covariates were quantified by adding 234 each of them consecutively to a basic ANOVA model. An ANOVA model was applied as (1) the 235 real/absolute amount of the fungal burdens were unknown and the ANOVA approach 236 overcame this by simply assuming that the panel samples 2, 3, and 5 are different; (2) the 237 linearity of the relationship Cq/logRLM was not guaranteed; (3) the difference between the 238 effects exerted by two levels of an added covariate would unlikely be constant across all Cq 239

  289 associated with Cq values generated by mtSSU assays amplifying WNA through to TUB assays 290 amplifying DNA only (ANOVA Kruskal-Wallis p<0.05) across all nucleic acid loads. The Cq values 291 obtained with the 3 panel specimens and the 20 PCR assays are represented Figure 2A. Indeed, 292 a mean 7 Cqs difference (200-fold variation) was observed between the TUB assay and the 293 assays targeting multicopy genes. 294 295 Monocentre intra-laboratory evaluation of 10 P. jirovecii qPCR assays

  but the lack of sensitivity of such a target gene hampers its use for clinical 337 diagnosis.338Interestingly, even when a given gene is targeted or the same primer/probes utilized, 339 differences can occur. Using the same mtLSU primers and probes and qPCR, 26 the 340 quantification obtained from five centres showed here a median 3-fold variation in the Cq (10-341 fold difference in concentration). In addition, we showed that five commercial mtLSU kits also 342 gave a 3-fold variation in the Cq. Overall, these differences can come from the different primer 343 concentration/amplification kits/enzyme or thermocyclers used. Indeed, technical 344 parameters can affect qPCR efficiency, such as the type of polymerase used, the concentration 345 of oligonucleotides or the method of Cq determination, suggesting that harmonization of the 346 methods requires more than standardization of the target gene or the primer sequences. 347 Although differences were already known among different DNA extraction methods or 348 qPCR assays,[START_REF] Hoarau | Evaluation of quantitative FTD-Pneumocystis jirovecii kit 515 for Pneumocystis infection diagnosis[END_REF][START_REF] Guillaud-Saumur | Comparison of a 518 commercial real-time PCR assay, RealCycler® PJIR kit, progenie molecular, to an in-house[END_REF][16][17][18][19][20][21][22] the superiority of WNA over DNA extraction has never been reported for P.

  349jirovecii, due to lack of prior of investigation. WNA RT-qPCR amplification kits were initially 350 developed for virology to allow cDNA and genomic DNA amplification in the same run to 351 detect both RNA or DNA viruses.[START_REF] Alanio | Diagnosis of Pneumocystis pneumonia : 561 switching from DNA quantification to RNA expression analysis 562[END_REF][START_REF] Alanio | Diagnosis of histoplasmosis : switching 564 from microscopy to RNA detection[END_REF] Currently, most of the automated extractors allow WNA 352 extraction to test for the presence of both RNA or DNA viruses, and DNA from other 353 pathogenic organisms in the same sample extract as part of a syndromic investigation.[START_REF] Alanio | Performance evaluation of multiplex PCR including Aspergillus-not 566 so simple![END_REF] Such 354 strategy has already been tested with success for detecting P. jirovecii DNA in samples 355 dedicated for influenza diagnosis. 24 Thus, since WNA extraction does not interfere with the 356 detection of P. jirovecii and even provides the best sensitivity in the present study, this 357 extraction method could be recommended.358The rationale for supporting the use of the quantitative assays with the highest 359 sensitivity might be contentious as low fungal loads are often considered to represent

  370 consider the difficulty of obtaining cultured material. Realistically, only surplus human 371 samples can be used and these are limited in both number and sample volume. Hence the 372 choice of pooling human samples and the limited number panel size used for comparison. In 373 pooling samples, it is acknowledged that potential biases are introduced. Ideally, it would have 374 been preferable to select a specific patient, to have an idea of the sequence of the P. jirovecii 375 main target genes and to propose individual panel specimens based on well-defined 376 organisms. Unfortunately, this is impossible due to limited sample volume. Since the BALF 377 material is usually limited after performing thorough routine diagnostics testing, such previous 378 analysis, providing 20 aliquots of several individual BALF extracts is practically impossible. The 379 pooling of BALF extracts also assesses the amplification process in the presence of the genetic 380 variability. In addition, it has been shown that PCP episodes are due to multiple genotypes in 381

  of PCP and our quality control did not address this possibility. For 385 example, MGB probes who are very sensitive to polymorphism/mismatch should not be used 386 to avoid the risk of false negative result.[START_REF] Gal | A misleading false-negative result of 586 Pneumocystis real-time PCR assay due to a rare punctual mutation: A French multicentre 587 study[END_REF] 387

  405of qPCR platform (thermocycler and Cq determination), because too many parameters varied 406 among the different centres. There are future plans to implement specific studies on these the impact of extraction, qPCR platform and the use of an internal control on 408 the performance of P. jirovecii qPCR. confirm the large variation between of qPCR assays to quantify P. jirovecii 412 nucleic acids. This large variation, up to 12 cycles (10,000-fold variation) for a given sample, 413 makes reaching a consensus threshold value to distinguish between low and high fungal loads 414 in respiratory samples difficult, but highlights optimal assays in regards to sensitivity. The most 415 sensitive association was shown to be RT-qPCR (amplifying combined RNA and DNA) targeting 416 the mtSSU gene and the less sensitive the qPCR assays (DNA only) targeting Major Surface 417 Glycoprotein and the unique beta-tubulin gene targets using. Interpretative thresholds might 418 be defined in the future, provided centres agree to use the most sensitive method. Other 419 factors not directly related to the PCR process that influence the result, such as the protocol 420 for processing the fluid must, however, be performed. Future harmonization of the results 421 should also take into account the procedure for obtaining broncho-alveolar lavage or other The authors declare no conflict of interest for this specific study, except GJ who is an employee 432 of OLM Diagnostics and GD who is an employee of PathoNostics. AA and SB are inventors in a 433 patent dealing with Pneumocystis RT-qPCR. 434 Page 21 of 68 http://mc.manuscriptcentral.com/tmmy Medical Mycology F o r P e e r R e v i e w O n l y

  Cq variation for each assay regarding the target gene and the type of material amplified (DNA or WNA) in all the three positive panel samples at dilutions 1:1 (A), 1:100 (B), and 1:1000 (C). It shows the gain in the ability to detect the target (lower Cq) in WNA assays compared to DNA and with a gradient from mtSSU to TUB genes (mtSSU>mtLSU>MSG>TUB) 95x272mm (300 x 300 DPI) Cq variation in all assays tested in the multicentre (A) and in the monocentre study (B). For the multicentre study, the panel specimen 1:1 (open circle) 1:100 (open square), 1:1000 (black circle) were tested in duplicates among 20 PCR assays tested (assay#2 was tested by 5 centres and assay#12 by 2 centres, see

Table 1 .

 1 Table 1 and the primers and probes used are available in They consisted in 1 commercial and 2 in-house RT-qPCR assays, and 4

	223	commercial and 3 in-house qPCR assays. The targets were mtLSU (7 assays); mtSSU (2 assays);
	224	and one nucleic gene: Major Surface Glycoprotein (MSG, 1 assays) (Table 1).
	225
	226	Statistical analysis

222

Supplemental

RESULTS 250 Multicentre inter-laboratory evaluation of P. jirovecii qPCR assays 251 The

  effects of shipping on the stability of DNA and WNA in the panel was investigated in Centre

	249 273	five centres (#2 centre B, C, O, #3 #16) (Table 2, Supplemental Figure 2 and Supplemental
	274	Table 3 and 4)
	275
	252	E after panel storage for 15 days at room temperature, after a freezing-thawing cycle. Using
	253	mtLSU DNA and mtSSU WNA detection, the coefficient of variation was 2.5%, 1% and 1% for
	254	mtLSU DNA and 2.2%, 1.4%, 1.9% for mtSSU WNA in samples, respectively (Supplemental
	255	Figure 1). All samples were found to be stable indicating that the distribution process did not F
	256	o r introduce bias in the results, even when shipping samples from Europe to Australia. o r
	257	P There were 16 participating laboratories across eight countries including seven P
	258	e e r European countries and Australia. Two centres tested more than one qPCR assays (two in e e r
		R e v R e v
		i e w i e w
		O n O n
		l l
		y y

Table 1 . Characteristics of the assays used in these studies.

 1 Five centres (B, C, E, J, O, E*) used the same in-house mtLSU assay as described 436 previously 26 and two centres (L, E*) used the same commercial assay.

		Assay	Centre ID	Target PCR Type Multicentre	Monocentre	In House	Commercial	Name of the kit or publication (PMID)	Manufacturer
		ID	(Multicentre			study	study	assays	assays
			study)
		1	A	mtLSU qPCR	X		X	26739439
		2	B, C, E, J, O, E* mtLSU qPCR	X	X	X	20946413
		3	D	mtLSU qPCR	X		X	Unpublished1
		5	E, E*	mtSSU qPCR	X	X	X	27672381
		6	H	TUB	qPCR	X		X	15767014
		7	I	mtLSU qPCR	X		X	Unpublished2
		8	M	mtLSU RT-qPCR X		X	Unpublished3
		9	E, E*	mtLSU RT-qPCR X	X	X	20946413
		10	E, E*	mtSSU RT-qPCR X	X	X	27672381
		11	E*	MSG	qPCR		X	X	11825961
		4	F	mtLSU qPCR	X			X	Pneumocystis jirovecii REAL TIME PCR, PG1905	Ridagene
		12	K, P, E*	mtLSU qPCR	X	X		X	PneumoGenius (PN-600)	Pathonostics
		13	N	mtLSU qPCR	X			X	PneumID (Ref OLM2008)	OLM
		14	K	mtLSU qPCR	X			X	Pneumocystis jirovecii (carinii) Real-TM, P2-50FRT Sacace Biotechnologies
		15	L, E*	mtLSU RT-qPCR X	X		X	FTD Pneumocystis jirovecii	Fast Track Diagnostics
		16	G	MSG	qPCR	X			X	LightMix Modular Pneumocystis jiroveci	Roche
		17	E*	mtLSU qPCR		X		X	RealStar® Pneumocystis jirovecii PCR Kit 1.0	Altona
		18	E*	mtLSU qPCR		X		X	RealCycler PJIR Kit	Progenie molecular
		19	E*	mtLSU qPCR		X		X	Pneumocystis jirovecii qPCR (BE-A994)	Bioevolution
	437	E*: evaluation in Centre E in the monocentre study ; qPCR: quantitative PCR; RT-PCR: reverse-transcribed quantitative PCR. Mitochondrial Large Sub-Unit
	438	(mtLSU); mitochondrial Small Sub-Unit (mtSSU); Major Surface Glycoprotein (MSG); beta-tubulin (TUB).
	439

Table 1 )

 1 . For the monocentre study, 9 human broncho-alveolar lavage fluids (BALF1 to BALF9 highlighted by different symbols) were tested in duplicates with 10 in house and commercial assays. It shows that assay #10 (mtSSU and RT-qPCR) gives the best Cq values.

	Medical Mycology

http://mc.manuscriptcentral.com/tmmy

Table 3 : Detailed Cq value for each assay tested in the multicentric study

 3

	Cq=45 corresponds to a negative qPCR result

http://mc.manuscriptcentral.com/tmmy Medical Mycology

Acknowledgments

http://mc.manuscriptcentral.com/tmmy Medical Mycology

genomic loads are indicated as 0, 1:1000, 1:100, and 1:1. The mean effect on Cq was evaluated after ANOVA as mean difference between two levels of relative genomic load.

Legends to Figures