Pierre Huchant

Emmanuelle Saillard

Denis Barthou

Patrick Carribault

Multi-Valued Expression Analysis for Collective Checking

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

Introduction

Collective operations and in particular synchronizations are widely used operations in parallel programs. They are part of languages for distributed parallelism such as MPI or PGAS (collective communications), shared-memory models like OpenMP (barriers) and languages for accelerators such as CUDA (synchronization within thread blocks, cooperative groups and at warp-level). A valid use of collective operations requires at least that their sequence is the same for all threads/processes during a parallel execution. An invalid use (collective error ) leads to deadlocks or undefined memory state that may be difficult to reproduce and debug. Indeed, these languages do not require that all processes reach the same textual collective statement (textual alignment property [START_REF] Zhang | Barrier matching for programs with textually unaligned barriers[END_REF][START_REF] Jakobsson | Replicated synchronization for imperative bsp programs[END_REF]). Finding which collective matches a given collective is needed for collective checking and requires to analyse the different concurrent execution paths of a parallel execution.

Aiken and Gay introduced the concept of structural correctness for synchronizations in SPMD programs, based on the notion of multi-valued and singlevalued variables [START_REF] Aiken | Barrier inference[END_REF]. A variable is said multi-valued if its value is dependent on the process identifier (single-valued otherwise). A program has structurally correct synchronization if all processes have the same sequence of synchronization operations. Thus, if a synchronization is executed conditionally, both branches of the condition have a synchronization or the condition expression is single-valued.

We can extend the notion of structural correctness to collectives. In this paper, we propose a novel method to detect collective errors in parallel programs. It combines an inter-procedural analysis to perform collective matching and a dataflow analysis to detect multi-valued variables. The first pass finds control-flow divergence that may lead to collective deadlocks while the second one filters out the divergences that do not depend on process identifier. We show on several benchmarks and applications that this combination is more accurate than the state-of-the-art analyses and resolves some correctness issues. The analysis has been implemented in the PARallel COntrol flow Anomaly CHecker [START_REF] Saillard | PARCOACH: combining static and dynamic validation of MPI collective communications[END_REF][START_REF] Huchant | PARCOACH Extension for a Full-Interprocedural Collectives Verification[END_REF] (PAR-COACH) framework and tested on benchmarks and real HPC applications, using MPI, OpenMP, UPC and CUDA.

Section 2 describes several deadlock situations in parallel programs. Section 3 presents PARCOACH analysis while sections 4 and 5 describe our multi-valued expression detection and its integration in PARCOACH to find collective errors. Section 6 gives related work on dependence analyses and verification tools. Section 7 shows experimental results and Section 8 concludes our work.

Motivation

This section illustrates four possible deadlock situations due to collectives in MPI, OpenMP, CUDA and UPC as presented in Figure 1 (from a to d).

The MPI code (Fig. 1a) contains two collectives: MPI Barrier and MPI Reduce. The call to MPI Barrier at line 17 is performed by all MPI processes, whereas the call MPI Reduce in g at line 3 may deadlock. Indeed, variable n is multi-valued in the condition expression line 14. Odd-ranked processes evaluate the conditional to true and potentially execute the reduce, while even-ranked ones evaluate it to false, hanging in the MPI Barrier at line 17. On the contrary, variable s is single-valued. Hence, all processes in g execute the reduce or none. The goal of our method is to statically report this situation, identifying the conditional at line 14, and only this one, as a potential cause for mismatched calls.

According to the OpenMP specification, the same explicit and implicit4 barriers (syntactically) should be executed by all threads. In practice, the OpenMP runtimes allow the execution of syntactically different barriers, as long as all threads execute the same number of barriers. The code Figure 1b is written in OpenMP. The #pragma omp parallel directive in function f defines r, containing the thread ID (l.10) and s, as private variables. The first barrier line 12 is encountered by all threads as it is in the global control flow. The barrier line 3 is either executed by all threads of the program or by none of them as s is singlevalued when entering function g at line 14. Because s becomes multi-valued at line 15, the barrier line 17 is conditionally executed. This leads to a deadlock situation if the number of threads is greater than 1 at runtime.

The CUDA code (Fig. 1c) manipulates multidimensional thread IDs through predefined variables such as threadIdx. In CUDA, synchronizations are valid if executed by all threads within the same block. Before the first synchronization, the array tile depends on thread IDs at line 7. As the array is shared among threads, they all share the same version after the synchronization. The synchronization at line 10 is conditionally executed depending on tile[0]. As this value does not depend on thread ID, there is no possible deadlock. Depending on the driver, the third synchronization at line 12 may lead to either a deadlock or an undefined memory configuration. This bug can be difficult to detect for a programmer in a real code as this is a silent synchronization error.

The code Figure 1d is written in Unified Parallel C (UPC), a PGAS language. The predefined variable MYTHREAD specifies thread index. In this code, because of the multi-valued expression at line 3, threads with odd ID will call nine barriers (l.5) while the others will call ten barriers (l.9). Although structurally correct, this code leads to a deadlock. Our analysis reports all collectives in a loop as potentially deadlocking. But note that our analysis would return a false positive if the two loops had the same number of iteration.

A static analysis on the previous codes only detects situations and causes of possible deadlocks. If the codes are executed with only one process, no deadlock can happen. Also, our analysis returns a false positive for some correct but structurally incorrect codes like the example Figure 1e. In addition to our analysis, we use PARCOACH instrumentation of programs explained in [START_REF] Huchant | PARCOACH Extension for a Full-Interprocedural Collectives Verification[END_REF] to handle such situations.

PARCOACH Control-flow Analysis

Our analysis first uses PARCOACH to find all conditionals (flow-divergences) that may lead to the execution of different collective sequences.

PARCOACH static analysis relies on the Parallel Program Control Flow Graph (PPCFG) representation of a program [START_REF] Huchant | PARCOACH Extension for a Full-Interprocedural Collectives Verification[END_REF]. The PPCFG is a program control-flow graph (CFG) where nodes are basic blocks and edges represent the possible flow of control between nodes. To build the PPCFG, PARCOACH reduces each function control-flow graph and replaces each callsite by the callee reduced CFG. Then, with a graph traversal of the PPCFG, PARCOACH computes the possible execution order r of each collective c and the iterated postdominance frontier (PDF + ) for collectives of the same type and order. For a set of collectives C r,c , PDF + (C r,c ) corresponds to the control-flow divergences that may result in the execution or non-execution of a collective in C r,c . Note that to handle communicators in MPI programs, PARCOACH analyses the program for each communicator separately.

Multi-Valued Expression Detection

PARCOACH finds all conditionals potentially leading to different sequences of collectives but reports false positives when conditionals do not depend on a multivalued expression. This section presents our multi-valued expressions detection.

Enhanced SSA Our analysis is based on the Static Single Assignment (SSA) form of the program. In SSA, variables are defined exactly once. Variables assigned in multiple statements are renamed into new instances, one per statement. This makes explicit def/use chains. When multiple control-flow paths join in the CFG, renamed variables are combined with a φ-function into a new variable instance. To capture control-flow dependences we compute an enhanced SSA where φ-functions are augmented with their predicates: φ(y 1 , ..., y k ) is transformed into φ(y 1 , ..., y k , p 1 , ..., p k ) with p i the conditionals responsible for the choice of the y i . These conditionals are determined by computing the PDF + of each argument y i of the φ-function as in [START_REF] Scholz | User-input dependence analysis via graph reachability[END_REF].

For C-like programs, variables that can be referenced with their address (address-taken variables), are only manipulated through pointers with load and store instructions in the SSA form. To compute def/use chains for address-taken variables, we rely on the principles exposed in flow-sensitive pointer analyses such as [START_REF] Sui | SVF: Interprocedural Static Value-flow Analysis in LLVM[END_REF][START_REF] Hardekopf | Flow-sensitive Pointer Analysis for Millions of Lines of Code[END_REF]: First a points-to analysis is computed to handle potential aliases among arrays and pointers. Then, each load q = * p is annotated with a function µ(o) for each variable o that may be pointed-to by p to represent a potential use of o at the load. Likewise, each store * p = q is annotated with o = χ(o) for each variable o that may be pointed-to by p to represent a potential def of o at the store. There is a special case to consider for shared variables. After synchronization (#pragma omp barrier in OpenMP, syncthreads in CUDA), shared variables have the same value for all threads. To create a new SSA instance that no longer depends on the value preceding the barrier, synchronizations are annotated with o = χ() for all shared variables o. Then a context-sensitive Mod-Ref Analysis is performed to capture inter-procedural uses and def as described in [START_REF] Sui | Detecting Memory Leaks Statically with Full-Sparse Value-Flow Analysis[END_REF]. The purpose of this analysis is to capture the variables referenced and/or modified in functions through pointers. Each callsite cs is annotated with µ(o) for each variable o indirectly referenced in the called function. Similarly, each callsite is annotated with o = χ(o) to generate a new instance of o for each variable indirectly modified in the called function. For each address-taken variable referenced or modified in a function, a χ function is inserted at the beginning of the entry node of the CFG and a µ function is inserted at the end of the exit node of the CFG to represent their initial and final values. Finally, all address-taken variables are converted into an SSA form. This results in an augmented SSA with value and control dependences, and additional statements in SSA describing the effects of pointer manipulations. All possible def/use chains are built inside the SSA notation. This simplifies the construction of a dependence graph.

PDCG: Program Data-and Control-flow Dependence Graph

A program data-and control-flow dependence graph (PDCG) is built from the enhanced SSA by connecting the def of each variable with its uses, following the rules presented in Figure 2. The PDCG captures inter-procedural dependences (represented by edges between variables from different functions) but its construction only requires to analyze each function once. This graph is used to find all variables/expressions that are multi-valued. To that end, we identify the source statements that generate processes identifier and spread the dependencies following the edges of the PDCG. The first four rules (from OP to STORE) are based on the work in [START_REF] Sui | SVF: Interprocedural Static Value-flow Analysis in LLVM[END_REF] using similar notations. Our differences are highlighted in grey.

OP and PHI rules correspond to straightforward data-and control-flow dependences. For an operation ∶ z = x op y, the def of x and y at lines ′ and ′′ are connected to the def of z at line . For a φ statement ∶ v 3 = φ(v 1 , v 2 , ..., p i , ..), the defs of the old SSA instances v 1 and v 2 at 1 and 2 are connected to the def of the new SSA instance v 3 at . For each predicate p i , the def of p i at 3 is connected to the def of v 3 at to handle the control-flow dependence.

LOAD and STORE rules take into account alias information for load and store statements. For a load statement ∶ q = * p, the def of each object o at ′′ pointed to by p is connected to the def of q at . We also add a link from the def of p at ′ to the def of q at to denote the dependence of q with the array index. Indeed, this corresponds to the case where * p is in the form of A[e] with e an expression. If e is multi-valued, then q is multi-valued. Similarly, for each store instruction ∶ * p = q annotated with [o 2 = χ(o 1 )], the defs of q and p are connected to o 2 . However, we do not connect o 1 to o 2 since we assume that the old value o 1 is overwritten with o 2 (strong update).

The CALL rule handles inter-procedural dependences. At each callsite cs ∶ r = f (..., p, ...), the def of the effective parameter p is connected to the formal parameter q in f . Furthermore, the def of the return value x in f is connected to the def of r at cs . To handle indirect value-flows for address-taken variables, given a callsite annotated with [µ(o 1 )] [o 2 = χ(o 1 )], the def of o 1 in the calling function is connected to o 3 , the first def of o in f . Similarly, the last def of o in f denoted o 4 , is connected to the def of o 2 at cs .

Rule

Statement (SSA form) Edges in the PDCG Value Flow Dependence OP ∶ z = x@ ′ op y@ ′′ z@ ↩ x@ ′ z@ ↩ y@ ′′ PHI ∶ v3 = φ(v1@ 1, v2@ 2, ..., p@ 3, ...) v3@ ↩ v1@ 1 v3@ ↩ v2@ 2 v3@ ↩ p@ 3

LOAD ∶ q = * p@ ′ [µ(o@ ′′ )] q@ ↩ o@ ′′ q@ ↩ p@ ′ STORE ∶ * p@ 1 = q@ 2 [o2 = χ(o1@ 3)] o2@ ↩ q@ 2 o2@ ↩ o1@ 3 o2@ ↩ p@ 1 CALL cs ∶ r = f (..., p@ 1, ...) [µ(o1@ 2)] q@ 3 ↩ p@ 1 r@ cs ↩ x@ 6 [o2 = χ(o1)] f (..., q@ 3, ...) { o3@ 4 ↩ o1@ 2 [o3@ 4 = χ()] ... [µ(o4@ 5)] return x@ 6 }
o2@ cs ↩ o4@ 5 Optimization

PHI

∶ * p@ 1 = q@ 2 [o2 = χ(o1@ 3)] o4@ ′′ = fuse(o2@ , o3@ ′ , o4@ ′′ ) ′ ∶ * p@ 1 = q@ 2 [o3 = χ(o1@ 3 )] remove(o4@ ′′ ↩ pred@ 4 ) ELIM ′′ ∶ o4 = φ(o2@ , o3@ ′ , pred@ 4 ) RESET cs ∶ reset(buf@ 1 , ...) [µ(o1@ 2)] remove(o2@ cs ↩ o4@ 5)

[o2 = χ(o1)] reset(buf@ 3, ...) { [o3@ 4 = χ()] ... [µ(o4@ 5)] } Collective Checking COND ∶ collective(...) with execution order r collective@ ↩ cond@ ′ For all BB ∈ PDF + (C r,collective@ ) matching: ... br cond@ ′ , label1, label2 Fig. 2: Building Rules for the PDCG.

PHI ELIM and RESET both correspond to edge removal optimizations. After augmenting φ-nodes with their predicates, false control dependences can appear if every operand of a φ-node denotes the same value. This occurs in particular when considering two identical function calls in two branches of an if..then..else construct. Even if these two calls use the same single-valued parameters, the returned value will still depend on the predicate of the conditional (augmented SSA). To tackle this issue, the PHIELIM rule fuses such φ-nodes with their operands and disconnects the predicates. In distributed memory, after a value-sharing collective such as an all-to-all collective, the communicated buffer has the same value for all processes. This implies that this buffer does not depend on processes ranks after such collective, whatever its rank-dependence before the collective. To handle this situation, the RESET rule disconnects the path from the old SSA instance of the buffer to its new SSA instance after a value-sharing collective (o 1 @ 2 ↪ o 3 @ 4 ↪ o 4 @ 5 ↪ o 2 @ cs ). The same rule applies to any value-sharing collective where all processes receive the same result such as MPI Allreduce or MPI Broadcast.

Finally, to detect collectives that may not be executed by all processes/threads, we rely on PARCOACH analysis. Each collective c of execution order r is connected to all conditionals in PDF + (C r,c ) (COND rule).

Collective Errors Detection

We use the PDCG to track values and nodes that depend on processes identifiers, flooding the graph from seed functions returning IDs or variables allowing tasks to identify themselves: MPI Comm rank and MPI Group rank in MPI, omp get thread num for OpenMP. In UPC and CUDA, the seed is a variable: MYTHREAD and threadIdx.*. We use the dependence information from the PDCG to filter out single-valued conditionals from the PDF + of potentially unsafe collectives and thus reduce the number of false positives in PARCOACH. The augmented SSA takes into account value and control dependencies and the pointsto analysis provides the dependencies through aliases. Note that thanks to the PDCG, our analysis can be path sensitive: An expression may be multi-valued or not, depending on the preceding calling context. Algorithm 1 describes our whole collective errors detection.

Step 1 represents PARCOACH control-flow analysis (see sec. 3) while steps 2 and 3 respectively detect multi-valued expressions and build the PDCG (see sec. 4). Finally, step 4 filters out single-valued conditionals and outputs warnings for potential collective errors.

Example Figures 3a and3c Depending on whether its value is odd or even, the execution flows either to the basic block labelled if.then or the basic block labelled if.else. These two control-flow paths join at basic block if.end and a φ-function is inserted to combine the values of reg1 and reg2 into variable n. The predicate cmp1 is added to the φ-function to indicate its value depends on cmp1.

Figure 3b shows the PDCG corresponding to this example. Rectangle nodes represent collectives. Diamond and circle nodes respectively represent definitions of address-taken and top-level variables (variables never referenced by their address). The seed function is MPI Comm rank line 7 and the first multi-valued object is o'2. All library functions have mocked-up CFGs, tagging output val- ues as multi-valued when necessary. The graph highlights the rank-dependent path from o'2 to MPI Reduce in g passing through the conditional cmp2 in f.

In this example, the execution of MPI Reduce depends on the value of cmp in g and the call to g depends on the value of cmp2 in f. Hence, MPI Reduce@ 6 is connected to cmp and cmp2. However, there exists no path from o'2 to cmp as cmp does not depend on processes ranks. The execution of MPI Barrier is not governed by any conditional. MPI Barrier@ 28 is then not connected to any node in the graph and it cannot be reached from a seed statement. Since the only collective highlighted in the graph corresponds to MPI Reduce in g and only one of the two conditionals governing its execution is highlighted, our new analysis only issues a warning for the multi-valued conditional line 22 in f and the call to MPI Reduce in g.

Related Work

This section summarizes works on dependence analyses and gives an overview of existing tools for collective errors detection in parallel programs.

φ-functions with predicates in order to handle control dependences. However, address-taken variables are not transformed into an SSA form.

Collective Error Detection Techniques

Static analyses operate on the source code of an application and have the advantage of not requiring execution. They are usually based on model checking and symbolic program execution, limiting their applicability to small and moderate sized applications (the number of reachable states to consider is combinatorial). TASS [START_REF] Siegel | Automatic formal verification of mpi-based parallel programs[END_REF] and CIVL [START_REF] Siegel | Civl: the concurrency intermediate verification language[END_REF] use this approach. They rely on symbolic execution and require source code annotations to detect collective errors in MPI programs. The OpenMP Analysis Toolkit (OAT) [START_REF] Ma | Symbolic Analysis of Concurrency Errors in OpenMP Programs[END_REF] uses the same method for OpenMP programs by exploring all program paths under symbolic values. SimGridMC [START_REF] Pham | Verifying MPI Applications with Sim-GridMC[END_REF]is a model checker for MPI applications. It uses Dynamic Partial Order Reduction and State Equality techniques to face the state space explosion problem. UPC-SPIN [START_REF] Ali | UPC-SPIN: a framework for the model checking of UPC programs[END_REF] generates finite models of UPC programs in the modeling language of the SPIN model checker. For CUDA programs, we can mention GPUVerify [START_REF] Betts | GPUVerify: A Verifier for GPU Kernels[END_REF] that statically checks that all threads execute the same barriers syntactically. Unlike our analysis, the method does not give a precise feedback in case of a potential error. PARCOACH combines an inter-procedural control-flow analysis with a selective instrumentation to find MPI and OpenMP collective errors. The method is limited to control-flow information and returns many false positives. Our new analysis overcomes this limitation and extends the collective verification to other parallel programming models. The method presented by Zhang and Duesterwald in [START_REF] Zhang | Barrier matching for programs with textually unaligned barriers[END_REF] is the closest to our work. It detects synchronization errors with an inter-procedural barrier matching technique for SPMD programs with textually unaligned barriers. Compared to our analysis, this method is only focused on MPI and OpenMP synchronizations and has no pointer analysis.

Unlike static tools, dynamic tools do not report false positives. However, they are dependent on the input data set and may miss errors (false negatives). PARCOACH instruments non verifiable programs to verify if processes/threads are going to call the same collective at a particular step of execution, preventing deadlocks from happening. This instrumentation is similar to what dynamic tools like MUST [START_REF] Hilbrich | Runtime MPI Collective Checking with Tree-based Overlay Networks[END_REF] or UPC-CHECK [START_REF] Coyle | UPC-CHECK: a scalable tool for detecting run-time errors in Unified Parallel C[END_REF] do. However, the instrumentation starts with the first collectives that may deadlock, avoiding a full instrumentation of programs.

Experimental Results

Our analysis is implemented as a pass in the LLVM framework 3.9 integrated into the open source software PARCOACH 5 .

Figures 4 and5 show the impact of our multi-valued expression analysis on PARCOACH. Figure 4 Filtered out warnings and conditions (higher is better) warnings conditions Fig. 4: Percentage of warnings and conditionals filtered by our multi-valued analysis. 100% means that the analysis proves the program is collective error free.

analysis on 3 HPC applications (MILC 6 , Gadget 7 and MPI-PHYLIP [START_REF] Ropelewski | MPI-PHYLIP: Parallelizing Computationally Intensive Phylogenetic Analysis Routines for the Analysis of Large Protein Families[END_REF]), 3 mini HPC applications (CoMD and miniAMR from the Mantevo project [START_REF] Heroux | Improving Performance via Mini-applications[END_REF] and Hydro 8 ) and 5 widely used benchmarks (HPL 9 , IOR 10 , AMG 11 , NAS IS 12 , and the CUDA benchmarks from Rodinia 13 ). In the figure, warnings are collectives that may lead to deadlocks, and conditions correspond to conditionals governing the execution of unsafe collectives. The initial number of warnings and conditionals found by PARCOACH is given at the top of each bar. 100% for a warning bar means that the application is collective error-free (all warnings are removed by our analysis, the code is proved safe), 0% means that our analysis has no impact. For MILC, 91% of the 498 warnings have been removed. PARCOACH now reports 45 warnings. As shown in the figure, about half conditionals are filtered out by our analysis for most applications and all warnings are removed for Coral AMG OMP, MPI IS, UPC IS and Rodinia. Figure 5 gives the percentage of collectives tagged as potentially deadlocking by our analysis. The total number of collectives is given at the top of each bar. In the figure, seven applications have less than 20% of collectives potentially deadlocking.

To highlight the functionality of our analysis, we created a microbenchmark suite containing programs from multiple sources with correct and incorrect use of MPI collectives 14 . We compare the performance of the method presented in [1] and PARCOACH using our multi-valued analysis (PDCG), SVF and Parfait. Table 1 shows the results. Our analysis always detect collective errors compared to the others. For the remaining false-positive results, a more precise dependence analysis is required. This is left for future work.

Conclusion

This article presents a new static/dynamic method to verify that a parallel program has structurally correct collectives. The analysis resorts to an interprocedural static analysis that can prove in some cases that a program is free of collective error. The method has been applied successfully on different languages and is implemented in PARCOACH. Experiments show that our analysis leads to significant improvement over existing PARCOACH. Furthermore, through a more precise use of alias and control dependences, our static analysis outperforms existing data-flow analyses bringing additional preciseness (removing spurious dependencies) and correctness (adding missing dependencies).

7

 7

Fig. 1 :

 1 Fig. 1: Examples of collective issues and a correct program not verifiable by our analysis.

  show the enhanced SSA for the MPI code Figure 1a. The call to MPI Comm rank is annotated with a χ function to denote the indirect definition of object o'1 pointed-to by r. This generates a new SSA instance denoted o'2. Then the object o'2 pointed-to by r is loaded in reg0.

// STEP 1 .Algorithm 1 :

 11 PARCOACH Control-flow Analysis Input: PPCFG, Output: ⋃ r,c Or,c ▷ Create the set Or,c = P DF + (Cr,c) for each collective name c of execution order r // STEP 2. Multi-Valued Expression Detection Input: SSA, Output: eSSA ▷ Build an enhanced SSA (eSSA) that captures data-and control-flow dependencies // STEP 3. PDCG Construction Input: eSSA, seeds, Output: PDCG ▷ Build a PDCG to find all multi-valued expressions and variables from seed functions and variables // STEP 4. Filter-out single-valued conditionals Input: ⋃ r,c Or,c, P DCG, Output: O for c in collective names of execution order r do for each node n in Or,c do if n is single-valued in PDCG (there is an edge between c and n in the PDCG) then Or,c ← Or,cn ▷ remove the node end O ← O ∪ (c, {Or,c}) end Output nodes in O as warnings return O Collective Error Detection.

  displays the percentage of warnings and conditionals filtered out with our multi-valued analysis compared to the initial PARCOACH

Fig. 5 :

 5 Fig. 5: Percentage of collectives potentially deadlocking.

Table 1 :

 1 Multi-valued detection comparison between the work in[START_REF] Zhang | Barrier matching for programs with textually unaligned barriers[END_REF] and PARCOACH (PAR.) using our PDCG, SVF and Parfait. FP = false positives, FN = false negative.

	Program Name	Origin Description	Deadlock	Zhang et. al [1] PDCG SVF Parfait PARCOACH
	field-sensitive	Hydro Structure with a	no	FP	FP	FP	FP
		multi-valued field
	index-dep	PAR. Use of an array	yes	✓	✓	FN	✓
	phi-cond	PAR. Control-flow dependence	yes	✓	✓	FN	✓
	pointer-instance	PAR. Figure 1b	yes	✓	✓	✓	FP
	pointer-alias	PAR. Use of aliases	yes	FN	✓	✓	FP
	barrierReduce	CIVL Collective mismatch	yes	FN	✓	✓	✓
	barrierScatter	CIVL Collective mismatch	yes	FN	✓	✓	✓
	BcastReduce bad	CIVL Collective mismatch	yes	FN	✓	✓	✓
	mismatch-barrier	PAR. Collective mismatch	yes	✓	✓	✓	✓
	mismatch barrier com PAR. Collective mismatch	yes	✓	✓	✓	✓
	mismatch barrier nb	PAR. Collective mismatch	yes	✓	✓	✓	✓
	MPIexample	PAR. Figure 1a	yes	FN	✓	FN	FN
	noerror barrier	PAR. Correct usage of barrier	no	✓	✓	✓	✓
	not verifiable	PAR. Figure 1e	no	FP	FP	FP	FP
	loop barrier	PAR. Figure 1d	yes	✓	✓	✓	✓

There is an implicit barrier at the end of all worksharing constructs, unless a nowait clause is specified.

PARCOACH is available at https://github.com/parcoach/parcoach

1 define void f () { 2 entry : Fig. 3: Enhanced SSA form of the MPI code Figure 1a and its corresponding PDCG.

Dependence Analyses Techniques

Dependence analyses are the cornerstones of many optimizations/analyses in compilers. For instance, dependences are used for Taint Analysis [START_REF] Arzt | FlowDroid: Precise Context, Flow, Field, Objectsensitive and Lifecycle-aware Taint Analysis for Android Apps[END_REF][START_REF] Tripp | TAJ: Effective Taint Analysis of Web Applications[END_REF][START_REF] Shankar | Detecting Format String Vulnerabilities with Type Qualifiers[END_REF] to determine how program inputs may affect the program execution and exploit security vulnerabilities, Information Flow Tracking [START_REF] Denning | Certification of Programs for Secure Information Flow[END_REF]14,[START_REF] Yin | Panorama: Capturing Systemwide Information Flow for Malware Detection and Analysis[END_REF][START_REF] Sabelfeld | Language-based Information-flow Security[END_REF] to prevent confidential information from leaking, static Bug Detection [START_REF] Laguna | Pinpointing Scale-dependent Integer Overflow Bugs in Large-scale Parallel Applications[END_REF][START_REF] Ye | Accelerating Dynamic Detection of Uses of Undefined Values with Static Value-Flow Analysis[END_REF] or code optimization and parallelization (e.g. the polyhedral model [START_REF] Feautrier | Dataflow Analysis of Array and Scalar References[END_REF]). One of the difficult issues when computing data dependences is to deal with pointers/memory aliases and non scalar variables (e.g. arrays, structures). In SVF [START_REF] Sui | SVF: Interprocedural Static Value-flow Analysis in LLVM[END_REF] the authors annotate load and store instructions with µ and χ functions to transform address-taken variables into an SSA form. However, they do not take into account the possible dependence of the pointer itself (through an array index for instance) when they build the data dependence graph.

Many of the aforementioned analyses only consider data dependences whereas Slowinska et al. [START_REF] Slowinska | Pointless Tainting?: Evaluating the Practicality of Pointer Tainting[END_REF] showed that omitting control dependences can be a huge source of false negative results. In [START_REF] Bao | Strict Control Dependence and Its Effect on Dynamic Information Flow Analyses[END_REF], the authors introduced the concept of Strict Control Dependences to reduce the number of false positives in Taint Analyses and Lineage Tracing. In Parfait [START_REF] Cifuentes | Parfait: Designing a Scalable Bug Checker[END_REF] the authors propose to extend


Torbjørn Cunis

Jean-Philippe Condomines

Laurent Burlion

Local stability analysis for large polynomial spline systems

Keywords: Nonlinear analysis, stability analysis, switching functions, polynomial methods, Lyapunov function

published or not. The documents may come

Introduction

Recently, several works on polynomial fitting have been led and provide a constructive method for determining models based on analytical computation due to their continuous and differentiable nature. Computing the exact region of attraction for this kind of nonlinear dynamics is very hard if not impossible. Therefore, researchers have focused on determining polynomial Lyapunov functions for polynomial systems building upon sum-of-squares [START_REF] Topcu | Local stability analysis using simulations and sum-of-squares programming[END_REF][START_REF] Topcu | Robust Region-of-Attraction Estimation[END_REF][START_REF] Khodadadi | Estimation of region of attraction for polynomial nonlinear systems: A numerical method[END_REF] including extensions to rational and composite Lyapunov-functions [START_REF] Valmorbida | Region of attraction estimation using invariant sets and rational Lyapunov functions[END_REF][START_REF] Tan | Stability region analysis using polynomial and composite polynomial Lyapunov functions and sum-of-squares programming[END_REF][START_REF] Chesi | Rational Lyapunov functions for estimating and controlling the robust domain of attraction[END_REF]. However, when polynomials are unsuitable to represent system dynamics, piecewise defined polynomials such as splines [START_REF] Coen | A new approach to linear regression with multivariate splines[END_REF] provide more tractable models, requiring extended effort when determining local stability: it is well known for example, that stability of the subsystems does not guarantee stability of the entire system [START_REF] Johansson | Computation of piecewise quadratic Lyapunov functions for hybrid systems[END_REF]. Therefore, approximation techniques have been developed for the estimation of the region of attraction of piecewise systems. Early work was limited to a priori given, multiple quadratic Lyapunov functions [START_REF] Amato | Estimation of the domain of attraction for a class of hybrid systems[END_REF] and could only provide very rough estimates of the region of attraction. In [START_REF] Chen | Stability Analysis and Region-of-Attraction Estimation Using Piecewise Polynomial Lyapunov Functions: Polynomial Fuzzy Model Approach[END_REF] though, analysis of polynomial fuzzy models is considered employing again composite Lyapunov functions; taking each the point-wise extremum, this approach This paper was not presented at any IFAC meeting. Corresponding author T. Cunis.

Email addresses: tcunis@umich.edu (Torbjørn Cunis), jean-philippe.condomines@enac.fr (Jean-Philippe Condomines), laurent.burlion@rutgers.edu (Laurent Burlion).

provides directly a continuous function. It is worth noting that this might come at the cost of a large number of decision variables and that there are no relaxations with respect to the respectively active subdomains. [11] proposes a further approach using multiple Lyapunov functions for switching hybrid systems with polyhedral subdomains which share a boundary in the origin. Sum-of-squares complexity was discussed [START_REF] Chesi | On the Complexity of SOS Programming and Applications in Control Systems[END_REF]; where subdomains are considered for stability, the complexity increases with the number of bounding constraints.

The present paper focuses on a new formulation of the region of attraction estimation for large piecewise systems of local polynomial dynamics and polynomial domains, such as switching systems and multivariate splines, within the sum-of-squares framework. We present preliminary results and extend previous work to switching systems (in Sec. 2). The main result, an algorithm for splines, is discussed in Sec. 3, and is applied to an engineering example in Sec. 4. The appendix illustrates an extension to multiple Lyapunov functions. For the implementation of the constraints, in particular the polynomial containment problem (Lemma 1), we rely on the semidefinite programming techniques of [START_REF] Topcu | Local stability analysis using simulations and sum-of-squares programming[END_REF][START_REF] Parillo | Semidefinite programming relaxations for semialgebraic problems[END_REF]. A concise discussion of the V-s-iteration is given in [START_REF] Chakraborty | Nonlinear region of attraction analysis for flight control verification and validation[END_REF].

Preliminaries

Consider the autonomous system ẋ = f (x) given by the k ∈ N ordinary differential equations

ẋ = f i (x) , if x ∈ Φ i (1) for 1 ≤ i ≤ k, where x ∈ R n , f i ∈ R [x]
n , and Φ 1 , . . . , Φ k are intersections of polynomial inequalities Assumptions f (0) = 0, i.e., the origin is a stationary point of f . Lemma 1 Let p, q 1 , . . . , q k ∈ R [x]; we have

Ω ϕ≤x0 = def {x ∈ R n | ϕ(x) ≤ x 0 } with ϕ ∈ R [x], x 0 ∈ R n ,
Ω q1≤0 ∩ • • • ∩ Ω q k ≤0 ⊆ Ω p≤0 if there exist s 1 , . . . , s k ∈ Σ [x] such that k i=1 s i q i -p ∈ Σ [x].
Lemma 1 is an explicit formulation of [1, Lemma 2]. We then write with s = (s 1 , . . . , s k ),

s i Ω qi≤0 ⊆ Σ Ω p≤0 (2)
and say "s proves" the inclusion. We also say that "q i , p solve" the inclusion.

Problem 2 Let p ∈ Σ [x]; solve the optimisation max β ,γ ∈R V∈R[x] β such that Ω p≤β ⊆ Ω V≤γ (3) Ω * V≤γ ⊆ {x | ∇V(x) f (x) < 0 } (4)
with V(•) positive definite 2 and V(x) = 0.

If (V, β , γ ) solve Problem 2 for p, then Ω V≤γ is an invariant subset of the region of attraction [1, Lemma 1] and admits the largest inscribing region Ω p≤β . Simultaneously searching for an optimal function V(•) while proving invariance of

Ω V≤γ involves bilinear terms. If the degree of V ∈ R [x] is restricted and f ∈ R [x]
, the V-s-iteration [START_REF] Topcu | Local stability analysis using simulations and sum-of-squares programming[END_REF][START_REF] Chakraborty | Nonlinear region of attraction analysis for flight control verification and validation[END_REF] solves Problem 2 by alternatingiteratively searching for V, γ , and β maximal such that s 1 , s 2 prove (3) and ( 4).

Remark 3

The sum-of-squares formulation is limited to nonnegativity; we thus make use of that p(x

) < 0 if p(x) ≤ -|x| 2 2 for p ∈ R [x],
x = 0, and > 0.

1 A1, . . . , A k form a set partition of a body K if and only if they are pairwise interior-disjunct and i Ai = K.

2 A continuous function ϕ is said to be positive definite (p.d.) if ϕ(•) > 0 everywhere except the origin and ϕ(0) = 0.

We write

Ω V,f, = def x ∇Vf (x) ≤ -x T x for V ∈ R [x], f ∈ R [x]
m , and > 0.

If there exists a single function V : R m → R p.d. such that ∇Vf i (x) < 0 for all x ∈ A * and all 1 ≤ i ≤ k, then A is also invariant for the piecewise system and V is a Lyapunov function of the switching dynamics. However, this requirement is unnecessary strict when it comes to the subsystems that are not active [START_REF] Johansson | Computation of piecewise quadratic Lyapunov functions for hybrid systems[END_REF]; indeed, the following suffices:

Corollary 4 Let V : R m → R be continuous p.d. with V(0) = 0 and A = Ω V≤α for some α ∈ R; if ∀x ∈ A * . (x ∈ Φ i ⇒ ∇Vf i (x) < 0) (5)
for all

1 ≤ i < j ≤ k, then A = i (A ∩ Φ i ) is invariant. Lemma 1 encodes (5) into a polynomial sum-of-squares problem, recalling Φ i is an intersection of polynomial inequalities. Now, i (Ω V≤γ ∩ Φ i ) is invariant if s 2,i Ω * V≤γi ∩ Φ i ⊆ Σ Ω V,fi, , (6)
s 2,i ⊂ Σ [x], for all 1 ≤ i ≤ k and γ = min{γ 1 , . . . , γ k }.

Remark 5

The idea of Corollary 4, and of the paper, can be extended to piecewise defined functions

V(x) = V i (x) , if x ∈ Φ i , for 1 ≤ i ≤ k.
We further illustrate this in the appendix.

Spline systems

While we have not taken further assumptions on the Φ i , in the present literature the invariant set is commonly assumed to cover all domains. 3 For large spline systems with bounded domains and only local stability, each subdomain taken into account whilst not part of the invariant set adds an inactive boundary to the computational load. We therefore present an adapted algorithm to efficiently compute a region of attraction estimation for spline systems.

h ij = E(i, j) is defined and h ij (x) > 0. 4
We define some further notation: for i ∈ I, let adj[i] = def j h ij = E(i, j) is defined be the set of adjacent domains of Φ i , that is,

Φ i = j∈adj[i] {x | h ij (x) ≤ 0 }; V K (•)
denotes the function-candidate in the K-th iteration and I K ⊂ I will denote the set of domains covered by the invariant set With this notation, we can state Algorithm 1 computing the optimal estimate Ω V≤γ for a spline structure: If, for any iteration K, the invariant set Ω V K ≤γ is contained in the subdomains I K , it suffices to check invariance only of these subdomains and with respect to the boundaries in between, instead of proving Eq. ( 6) for all i ∈ I. Now, in order to examine whether an i ∈ I is covered by the optimal invariant set of V K , we preliminary compute the invariant set for some I ⊂ I -{i } and evaluate the distance of i ; only if i is "closer" than the boundary of the preliminary invariant set, i ∈ I K .

Ω V K ≤γ ; for i ∈ I, dist K [i ] is the distance of Φ i with respect to V K , defined as dist K [i ] = sup{γ | Ω V K ≤γ ∩ Φ i = ∅ } ; (7
The algorithm consists of three cascaded loops; the outer for loop over K of the basic V-s-iteration ("K-iteration"), an inner repeat-until loop determining I K ("I K -loop"), and inner-most for loops over the elements of I K . While the restriction to I K in line 8 reduces the problem size, the I K -loop itself adds to the run time. On the other hand, if V K and adj[i] ∩ I K remain unchanged, Eq. ( * i ) yields the same value γ i ; that is, after each i next added to I K , it suffices to re-execute line 8 for any i ∈ I K ∩ adj[i next ]. We refer to the thus modified algorithm as Algorithm 1b.

Proposition 7

After each repetition of the K-iteration in Algorithm 1(b), the following hold:

(1)

Ω V K ≤γ ⊂ i∈I K Φ i ; (2) Ω V K ≤γ is invariant.
PROOF. After each iteration of the I K -loop, we have that

∀i ∈ I K Ω * V K ≤γpre ∩ Φ i ⊆ {x | ∇Vf i (x) < 0 } ; (8) as Φ i ⊆ j∈I Ω hij ≤0 for any I ⊂ adj[i], adj[i] ∩ I K ⊂ adj[i],
γ pre ≤ γ i , and s 2,i proves ( * i ), 5 for all i ∈ I K ; and

Ω V K ≤γmin ⊂ i∈I K Φ i ; (9) 5 Using A1 ∩ A2 ⊆ A 1 ∩ A 2 for A1,2 ⊆ A 1,2 . Algorithm 1 Estimate invariant set Ω V≤γ ⊂ i∈I Φ i with I = I Kmax and V = V Kmax . 1: for K = 1 to K max 2: if K > 1 then 3: find V K p.d. solving s 1 Ω p≤β ⊆ Σ Ω V K ≤γ ( † †) ∀i ∈ I K . s 2,i Ω V K ≤γi ∩ ( * * ) j∈adj[i]∩I K Ω hij ≤0 ⊆ Σ Ω V K ,fi, 4: end 5:
I K := I K-1 6: repeat 7:
for i ∈ I K 8:

find γ i := max γ≥0 γ s.t. s 2,i ⊂ Σ [x] solves s 2,i Ω V K ≤γ ∩ ( * i ) j∈adj[i]∩I K Ω hij ≤0 ⊆ Σ Ω V K ,fi, 9:
end 10:

γ pre := min{γ i | i ∈ I K } 11: for i ∈ adj[I K ] 12: compute dist K [i ] as min γ≥0,x V K (x)=γ γ s.t. j∈adj[i ] h i j (x) ≤ 0 ( ‡ i ) 13:
end 14:

γ min := min i ∈adj[I K ] dist K [i ] 15:
if γ pre ≥ γ min then 16: for i ∈ I K 22: [START_REF] Coen | A new approach to linear regression with multivariate splines[END_REF]. Since γ pre < γ min for termination of the I K -loop, [START_REF] Amato | Estimation of the domain of attraction for a class of hybrid systems[END_REF] implies

i next := arg min i ∈adj[I K ] dist K [i ]
find β := max β≥0 s.t. s 1 ∈ Σ [x] solves s 1 Ω p≤β ⊆ Σ Ω V K ≤γ ( †) 23: end 24: end as γ min ≤ dist K [i ] for i ∈ I K and ( ‡ i ) implies
Ω V K ≤γ ⊂ i∈I K Φ i with γ = γ pre . Thus, Ω V K ≤γ ∩ Φ i = ∅ for all i ∈ I -I K 6 and therefore, Ω V K ≤γ = i∈I (Ω V K ≤γ ∩ Φ i ) is invariant by Corollary 4 for (8) holds. 2
As all loops are limited by either the number of elements in I or K max , Algorithm 1b terminates and Ω V≤γ is an invariant set of the given spline system. Initially, we have I 0 = {i ∈ I | 0 ∈ Φ i }, and, assuming a singleton I 0 = {i 0 }, V 1 = x T P x as solution to the polynomial Lyapunov equation for f i0 . 7

Asymptotic run time estimation

In order to compare the run time of the proposed approach for splines to the basic approach of the previous section, we count the total number of executions of line 8 and the number of decision variables s 2,i involved each time. 8 Here, we assume that the spline structure SP has k subdomains; every domain has (in average) M adjacent cells; the resulting invariant set covers R subdomains; and the number of iterations is chosen as K max = R. Consider now the following, distinct cases: in the worst case, the initial invariant set Ω V1≤γ covers all R subdomains; whereas in the average case, Ω V K ≤γ grows in each repetition of the K-iteration into one further subdomain i next . In both cases, I 0 is taken as singleton. Clearly, the asymptotic run time of the basic approach in both cases is equivalent to

T basic M (k) = R k(M + 1).
Algorithm 1b, in the worst case, repeats line 8 in the first iteration for R times, afterwards once each iteration: proving ( * i ) requires a single decision variable the first time (adj[i 0 ] ∩ {i 0 } = ∅) and m * + 1 decision variables from there on, where m * is the number of adjacent domains of i in I 1 and each repetition of the I K -loop adds one i next ; in the following R -1 iterations, line 8 is executed R times each with M + 1 decision variables; that is,

T worst M (R) = 1 + R M m * =1 (m * + 1) + (R -1) R (M + 1) .
In the average case, Algorithm 1b executes line 8 once in every iteration for each of the r * domains in I K with m * ≤ M decision variables; additionally, since I K+1 = I K ∪ {i next }, line 8 is executed for every adjacent subdomain of i next in I K (i.e., at most M times) with m * + 1 decision variables; that is,

T avg M (R) = 1 + R r * =2 r * (M + 1) + R M m * =1
(m * + 1) . 7 If I0 is not a singleton, P can be found as solution to

A T i P + P Ai < j∈adj[i]∩I 0 hij(x) ∀i ∈ I0
with Ai = ∂fi/∂x. 8 Line 8, being inside of all three loops, is the major difficulty if line 12 is efficiently computed using MATLAB's fmincon or a similar numerical method.

If M = 3 (M = 4), 9 the worst-case asymptotic run time is less than T basic M for R < 0.987 k (R < 0.981 k) and the average-case is less than 1 2 T basic M for R < 0.944 k (R < 0.933 k).

Application example

The short-period motion of a transport aircraft might be given as autonomous system

ẋ1 = x 2 (10) ẋ2 = f M (x 1 , x 2 ) , (11)
where x 1 is the angle of attack, x 2 is the pitch rate, and f M (•) is a 3rd-order piecewise polynomial model of the aerodynamic pitch coefficient defined as 5-by-5 rectangular spline (boundaries depicted in Fig. 1) with f M (0, 0) = 0. 10 After scaling the system to x = Dx with D ∈ R 2×2 , we compute the invariant set of the origin using Algorithm 1b. Such a problem has been discussed in [START_REF] Chakraborty | Nonlinear region of attraction analysis for flight control verification and validation[END_REF] for a polynomial model and in [START_REF] Cunis | Piecewise Polynomial Modeling for Control and Analysis of Aircraft Dynamics beyond Stall[END_REF] for a once-piecewise polynomial model. Algorithm 1b finds the optimal invariant set shown in Fig. 1, covering 20 of 25 domains, after 124 repetitions of the K-iteration and a run time equivalent of T = 10 098. The basic approach runs the same number of repetitions with a run time equivalent of T basic = 13 020.

In each repetition of the K-iteration, line 8 is executed multiple times subject to both the elements already in I K-1 and those added to I K during the inner I K -loop. The number of linear matrix inequality (LMI) variables to solve the sum-of-squares problems ( * i ) then varies with the number of elements in adj[i] ∩ I K . Table 1 gives details of the computations, including the number of 9 M = 3 and M = 4 are tantamount to planar systems with triangular and rectangular domains, respectively.

Conclusion

Extensions of sum-of-squares techniques such as the V-s-iteration for piecewise polynomial systems quickly grow infeasible for large systems, which spline models embody. In this article, we therefore presented an adapted algorithm for splines, relaxing the problem to the subdomains that are actually covered by the region of attraction estimate. We have proven correctness of our approach and demonstrated a worst case run time superior to the basic approach for regions of attractions spanning more than 90 % of the subdomains of a regular planar spline system. While this ratio will shrink for higher dimensions, so does the number of subdomains increase.

1 Notation

 1 forming a set partition of R m . The interior, boundary, and closure of A ⊆ R m are notated by intA, ∂A, and clA, respectively. A * = def A -{0}. The set of sum-of-squares polynomials is notated by Σ [x].

Definition 6 A

 6 spline system is a triple SP = I, f (•) , E , where I ⊂ N are the domains, f : I×R m are the piecewise nonlinear dynamics, and E : I × I → {R m → R} is the weighted switching relation, where the dynamics switch from f i to f j with i, j ∈ I for x ∈ R m if and only if

  ) at last, we extend adj[•] to 2 I with adj[I] = def i∈I adj [i] -I for I ⊂ I.

17 :I

 17 K := I K ∪ {i next } 18: end 19: until I K = I or γ pre < γ min 20: γ := γ pre 21:

Fig. 1 .

 1 Fig.1. Invariant set after K fin = iterations for aircraft short-period motion (solid: invariant set Ω V K ≤γ ; dotted: inscribing ellipsoid Ω p≤β ; dashed: subdomain boundaries).

Table 1

 1 LMI variables in ( * i ) averaged for each repetition K, comparing Algorithm 1b and the basic approach. 11 Details of the application example: number of elements in IK ⊆ I; run time equivalent; number of executions of line 8; and average number of LMI variables.

			Alg. 1b			basic
	K	#I K	T K	#l8	LMI	#I K	T K	#l8	LMI
	1	12	86	28	208	25	105	25	250
	2	12	44	12	356	25	105	25	378
	3	13	56	15	359	25	105	25	378
	4	14	60	16	360	25	105	25	378
	5	15	65	17	363	25	105	25	378
	6	20	111	29	365	25	105	25	378
	7	20	86	20	374	25	105	25	378
	. . .	. . .	. . .	. . .	. . .	. . .	. . .	. . .	. . .
	K fin	20	86	20	374	25	105	25	378

10

See

[START_REF] Cunis | Piecewise Polynomial Model of the Aerodynamic Coefficients of the Generic Transport Model and its Equations of Motion[END_REF]

for details.

Either by looking for global stability[START_REF] Papachristodoulou | Robust Stability Analysis of Nonlinear Hybrid Systems[END_REF] or choosing boundaries crossing the origin[START_REF] Amato | Estimation of the domain of attraction for a class of hybrid systems[END_REF] 11].

For a spline structure to behave like a spline, we tacitly understand that E(•, •) is irreflexive as well as symmetric with hij = E(i, j) = -hji if defined for i, j ∈ I.

Assuming that intΦi,j are disjunct if i = j.

Acknowledgements

The implementation of Algorithm 1 utilises the SOSOPT / SOSAnalysis toolboxes of the University of Minnesota.

A Appendix

The proposed algorithm can be modified for multiple Lyapunov functions:

Then, Eq. (3) holds if

The continuity condition (A.2) cannot be represented by sum-of-squares immediately. Papachristodoulou and Prajna [START_REF] Papachristodoulou | Robust Stability Analysis of Nonlinear Hybrid Systems[END_REF] suggested to add an equality constraint for each polynomial boundary, however leading to increased conservativeness for large systems.

Proposition 9 Eq. (A.2) holds if and only if

for Φ i , Φ j pairwise disjunct.

PROOF. Follows directly from clΦ

In line 3 of Algorithm 1(b), Eq. (A.4) holds for

Now, when adding a subdomain i next to the current I K (line 17), we search for

Remark 10 Instead of (A.3), we might require

for all 1 ≤ i ≤ k when searching for V 1 , . . . , V k .

Eq. (A.6) imposes a less strict constraint on those V i with β i > β .


Chloé Huetz
email: chloe.huetz@u-psud.fr

Samira Souffi
email: samira.souffi@u-psud.frvictor.adenis@u-psud.fr

Victor Adenis

Neural Code: Another breach in the wall? Response to the target article by Romain Brette: Is coding a relevant metaphor for the brain?

Brette presents arguments that query the existence of the neural code. However he has neglected certain evidence that could be viewed as proof that a neural code operates in the brain. Albeit these proofs show a link between neural activity and cognition, we discuss why they fail to demonstrate the existence of an invariant neural code.

By questioning the existence of the neural code, Romain Brette opens again a strong debate between representational views of the brain (cognitivism and computationalism) and sensorimotor/enaction theories [START_REF] O'regan | A sensorimotor account of vision and visual consciousness[END_REF]Noë, 2001, Varela et al 1991), his preference being the latter. According to his view, all cognitive functions, particularly action and perception are viewed as means to interact with the world, without the need to build internal representations of it. Neural activity during perception should be viewed as the result of the organism's interaction with the world taking into account all possible influences, such as its internal state and its actions resulting in a given percept.

Therefore, as the brain does not manipulate representations, it is senseless to try to decipher any code supposed to encrypt representations in neural activity. The results of three research fields focusing on proving that a particular neural code is at play should be addressed by Brette's review to strengthen his point.

First, in sensory physiology, research on tuning curves has been extended to naturalistic stimuli and is divided into two complementary approaches: encoding and decoding.

Based on models of the stimulus-response function, these approaches rely on the idea that neural activity encodes some features of the external world. Successful reconstructions of complex stimuli based on neural responses (decoding), or successful predictions of responses to new stimuli (encoding) are viewed as proofs that the neural code has been cracked. Interpreting these results in light of Brette's arguments seems necessary. Initially, the stimulus reconstruction method (decoding) was performed either with simple artificial stimuli [START_REF] Bialek | Reading a neural code[END_REF] or in peripheral sensory systems [START_REF] Warland | Decoding visual information from a population of retinal ganglion cells[END_REF][START_REF] Rieke | Naturalistic stimuli increase the rate and efficiency of information transmission by primary auditory afferents[END_REF]. More recently, studies have reconstructed natural stimuli from cortical responses [START_REF] Miyawaki | Visual image reconstruction from human brain activity using a combination of multiscale local image decoders[END_REF][START_REF] Naselaris | Bayesian reconstruction of natural images from human brain activity[END_REF], Akrabi et al 2019), opening the spectacular expectation to read subjects' percepts. In the auditory modality, encoding models were used to investigate neural selectivity to a variety of acoustic properties such as phonetic features [START_REF] Mesgarani | Phonetic feature encoding in human superior temporal gyrus[END_REF], pitch [START_REF] Oxenham | How we hear: The perception and neural coding of sound[END_REF], timbre and rhythm [START_REF] Woolley | Functional groups in the avian auditory system[END_REF]. To achieve good performance, the stimulus/response models used in decoding/encoding approaches rely on features such as trial averaging, statistics of natural stimuli, starting time of the stimulus. Thus, the right interpretation should be that an "ideal observer" with an a priori knowledge of the experimental design can infer the stimulus (in the decoding approach) or the neural response (in the encoding approach). Noteworthy, this field has led to an interesting drift from the idea of a fixed relationship between stimulus and neural responses to a more dynamic model, and is now tackling the mechanisms by which sensory responses are modulated by learning, context and history [START_REF] Williamson | Input-Specific Gain Modulation by Local Sensory Context Shapes Cortical and Thalamic Responses to Complex Sounds[END_REF][START_REF] Fritz | Active listening: task-dependent plasticity of spectrotemporal receptive fields in primary auditory cortex[END_REF][START_REF] Holdgraf | Rapid tuning shifts in human auditory cortex enhance speech intelligibility[END_REF].

Second, the field of neuroprosthetic devices offers demonstrations of causal links between neural code and brain functions. The most successful of these devices, cochlear implant (CI), operates with blunt stimulations of auditory nerve terminals. Despite a large current spread in the tympanic ramp, CI allows implanted subjects to have percepts and recover speech understanding. Even though there are huge differences between normal cochlea and CI, the fact that CIs restore hearing can be viewed as a proof that the neural code at play in the periphery has been deciphered and is successfully implemented in a prosthetic device. However, the CI settings that leads to speech comprehension considerably differ from one subject to another, as do the strategies leading to the largest evoked responses in auditory cortex [START_REF] Adenis | ECAP growth function to increasing pulse amplitude or pulse duration demonstrates large inter-animal variability that is reflected in auditory cortex of the guinea pig[END_REF].

Thus, in contrast with the genetic code that is invariant across cells and species, the neural code (understood as changes in neural activity in adaption to a CI) is probably specific for each individual and/or each type of neuron. In line with sensorimotor theories, success of CI shows that the brain is using a new input in a way it can interact again with the environment, which might be the basis of hearing restoration.

A third important field investigates the effect of disrupting a particular feature of neural activity on a cognitive skill. In the visual system, disruption of physiological activity in the primate area MT during presentation of moving stimuli bias the perceptive judgment of a behaving animal [START_REF] Salzman | Cortical microstimulation influences perceptual judgements of motion direction[END_REF][START_REF] Salzman | Neural mechanisms for forming a perceptual decision[END_REF] thus making the first link between neural code (understood as a pattern of activity of specific neurons) and behavioral performance. More recently, studies performed in the hippocampus have showed that disrupting the replay of spiking patterns occurring across neuronal ensembles during the sharp wave ripples profoundly alters the memory of previously acquired information [START_REF] Girardeau | Selective suppression of hippocampal ripples impairs spatial memory[END_REF][START_REF] Ego-Stengel | Disruption of ripple-associated hippocampal activity during rest impairs spatial learning in the rat[END_REF].

These data reinforce the notion that neuronal activity patterns do correlate with the acquired information. More importantly, associating a rewarding stimulation of the medial forebrain bundle with an hippocampal place cell activity induced a place preference at the place cell location [START_REF] De Lavilléon | Explicit memory creation during sleep demonstrates a causal role of place cells in navigation[END_REF], demonstrating causal links between a particular place cell's firing rate and a specific location memory. In all these examples, the exact neural activity feature (its firing rate or its temporal spike patterns) correlated with the animal's location is unknown, but causal relationships do exist. Yet, causality is not enough to define a neural code.

Clearly, more caution is necessary when discussing the neural code as overstatements made [START_REF] Ferster | Cracking the neuronal code[END_REF][START_REF] Panzeri | Cracking the Neural Code for Sensory Perception by Combining Statistics, Intervention, and Behavior[END_REF] tend to generate the illusion that (1) the same code operate in any sensory and motor systems, which is obviously not the case and (2) the brain's cognitive functions consist of manipulating encoded representations of the world, a theory that is controversial. Does this mean that the concept of neural code should be abandoned or should it be used to describe studies linking neural activity to brain function? We believe that the neural code definition should be freed from the notion of representation and we should clarify what we refer to when investigating the neural mechanism of brain functions.


Martin Nudischer

Hansgeorg Binz

Ma Hias Bachmann

Stefan Recker

S Y Mpo

O Tatin

Experimental Investigation of Two Centrifugal Fans in a Serial Arrangement

Purely theoretical descriptions for reliably designing serially arranged fans are provided in the relevant literature. In order to investigate such serial arrangements with their individual fans, two geometrically identical fans were directly connected for performance tests. The resulting characteristic curve and its composition were measured on a test rig in accordance with ISO , the results of which are described in this paper. In addition, the bearing loads in an axial and radial direction were examined for both fans. The information collected on the serial arrangement was contrasted to the single fan. The results can be concluded as follows: The first fan in the serial arrangement (LPF) can be considered as a single fan. The head coefficient and bearing loads are in the same range as for a single fan. Despite the fact that the fans used are geometrically identical, the second fan (HPF) in the serial arrangement exhibits different behavior. At low volume flow rates, the head coefficient of the HPF is lower than the head coefficient of the LPF. With increasing volume flow rates, the HPF exhibits a higher head coefficient. Moreover, the axial loads on the HPF are lower. It can additionally be demonstrated that with serially arranged fans not only an increase in pressure is possible, but that the attainable volume flow rate also rises.

INTRODUCTION

Centrifugal fans are widely used for gas compression in a variety of industrial areas, especially for ventilation or as components in industrial processes. The demand for higher pressure requires the fans to have wider operating ranges. For example, if an existing exhaust system is expanded, this will create additional resistance in the exhaust system. To overcome this additional resistance, the performance of the exhaust fan must be increased. One way to increase the performance is by accelerating the impeller, leading to an increase in the pressure and volume flow rate. If the improved performance achieved through the acceleration of the impeller is not sufficient, another solution will be necessary.

There are two options: either a completely new and more powerful fan substituting the old one, or an additional fan to support the old fan by creating a serial arrangement. In this context, the serial arrangement of fans is defined as a combination of two or more single centrifugal fans and not one fan with two or more impellers on one shaft. Aside from the advantage of higher pressure, there is also potential to make significant energy savings in partial load conditions. Furthermore, the modularity of the system achieved by combining two or more centrifugal fans creates more variability and an advantage in production costs by using common parts. Additionally, for different system requirements a new design of the fan is not always necessary. However, the coupled fans influence each other and are interdependent. Despite the fact that only few scientific investigations have been conducted into serially arranged centrifugal fans, many performed serial arrangements can already be found incorporated into plants. A fan manufacturer's current product portfolio up to a volume flow rate of approximately 100 m 3 /min is shown in Figure [ ]. The chart illustrates the pressure increase over the volume flow rate for two types of turbo-machines with different working principles. Side channel blowers and centrifugal fans are the two types of fluid energy machines which are considered in the analysis. Therfore the focus of the product analysis is on the low cost sector and disregards expensive special solutions such as high-drive turbochargers or similar. The gap in the characteristic map for demands for high pressures levels and simultaneously high volume flow rates is obvious. By combining two or more centrifugal fans in series, the available pressure can be increased. Moreover, the inefficient side channel blowers with an energy efficiency of almost % [ ] can be replaced by multi-stage centrifugal fans in order to save energy. Centrifugal fans can help to reduce energy consumption with an efficiency of a single-stage up to 85 % [ ]. In addition, a serial system with two or more fans is able to disconnect the fans independently, for example under partial load conditions. The potential for serially arranged fans is enormous considering the attainable performances in terms of pressure and volume flow rate. Depending on the requirements for serially arranged fans, only two or more centrifugal fans and appropriate connecting pipes are required in general. This connection leads mechanically to a serial arrangement. However, the effects arising from this connection are still unknown and there is no scientifically based approach for dimensioning serially arranged fans. This system of two or more fans is also an oscillating system with two or more initiating vibration sources. The key question is thus: What are the changes incurred by directly connecting two or more fans compared to a single fan? This paper is an initial step toward identifying the correlations within the system in terms of the characteristic curve and bearing loads. [ ] and Banzhaf [ ] described failures in power plants in which more than one fan is installed within the pipe system. However, the focus in these publications was not on the fans and their influence on each other. This investigation was conducted into the adjustment settings for the single fans in order to prevent damages to the ventilation system. Carolus investigated the phenomena of the rotating stall and surge of fans in systems with two centrifugal fans [ ]. Rotating stall and surge are a kind of unsteady flow occurring in fans at low volume flow rates. The focus of this investigation was on the incidence and cause of these phenomena. Schulze-Dieckhoff outlined how defined control interventions influence the behavior of serially arranged fans at non-stationary operating points [ ]. Axial fans or gas turbines in a multi-stage arrangement are well understood and comprise the state of the art. There are numerous design methods and instructions for reliably designing multistage axial fans and turbines [ ]. Due to the different operating principles of those rotating machines, the knowledge cannot be transferred to centrifugal fans. The serial operating mode is widespread in pumps. However, the argument in favor of multi-stage pumps is not necessarily to gain more pressure: In pumps, the phenomenon of cavitation occurs when the local pressure drops too far, which leads to the pump being damaged due to the formation and subsequent collapse of vapor bubbles [ ]. To avoid this phenomenon pumps are set to serial operation mode. Consequently, neither the knowledge of pumps nor axial fans can be transferred to the design of serially arranged centrifugal fans. As already mentioned, the attainable ratios of volume flow rate and pressure increase are limited by the peripheral speed of the impeller. Moreover, the required input power is directly proportionate to the speed increase. Eight times as much power is required to operate the fan at twice the speed. According to literature, the performance of serially arranged fans can be approximated by multiplying the total pressure to n 2 , the pressure increases quadratically with the quotient (n 2 /n 1 ) 2 . The volume flow rate rises linearly. Generally, the fan forms part of a larger system, connected on either the suction or the pressure side or both. The connected system constitutes an obstruction for the fan which can be simplified as the system resistance curve. The intersection of the fan characteristic curve and system resistance curve is the operating point (op) of the fan. This volume flow rate flows through the system at the defined pressure. There is only a minor gain in pressure increase for lower system resistance curves through higher rotational speeds (op to op ). The higher the system resistance curve, the higher the increase in pressure (op to op ) [ ]. Nevertheless, two or more fans can be arranged in a serial operating mode to increase the pressure even further. Especially at high system resistances the increase in available pressure is enormous (change in op to op ). However, the changes coming up through the direct connection and the additional compression are unknown. If, for example, the summarized pressure increase is unequally distributed across the different stages (with two or more identical fans at the same rotational speed), this will result in a change in pressure distribution in the housing of each fan. These different pressure distributions lead to volatile bearing forces and thus altered service life expectations, not to mention higher maintenance rates. To summarize, there are only few scientific contributions dealing with more than one fan in a serial arrangement. However, the focus in these contributions is either on flow phenomena such as rotating stall or surge or the correct action in the event of malfunctions. For this reason, the direct connection of two centrifugal fans will be examined in detail in order to fill the gap in insufficient design guidelines for serially arranged fans. The focus of this investigation is on the performance of the serial arrangement and the changes resulting from direct connection. Bearing loads will be assessed in addition to the performance tests.

Pressure increase

. STATE OF THE ART

. TEST RIG Two centrifugal fans are required for the detailed experimental examination of the serial arrangement. All installed fans (low-pressure, high-pressure and single fan) are geometrically identical so as to avoid distortion due to design differences. The specifications of the fan are summed up in Table . The final performance of the arrangement is measured on a test rig in accordance with the international standard EN ISO : [ ]. In Figure , a sketch of the setup is shown above while the Cordier diagram is displayed below, in which the examined (single) fan is marked. The low-pressure fan (LPF) sucks the air from the atmosphere. The high-pressure fan (HPF) receives the already compressed air and compresses it further. The two fans are connected by a diffusor with its length of 400 mm. Downstream of the arrangement, the flow rate is determined via an orifice plate. A shutter is placed at the end of the outlet pipe to adjust the operating point. Static pressures and temperatures are measured between the two fans, at the outlet of the HPF before and after the orifice, in addition to the ambient pressure, ambient temperature and humidity. The Cordier diagram establishes the link between the volume flow rate V, the pressure increase ∆p and the diameter of the impeller D and the speed n at the optimum efficiency of the single fan. This link results in the dimensionless speed number and diameter number. Defined as follows and calculated for the examined single fan: Speed number:

σ opt = 2 • √ π • n • √ V (2 • ∆p/ρ) 3/2 = • • • = 0.239 ( )
Diameter number:

δ opt = √ π 2 • D • 4 2∆p/ρ V 2 = • • • = 3.811 ( )
With its high diameter number and low speed number, the examined fan is therefore typical for radial fans. In order to exclude any uncertainty regarding the performance of the fans implemented, both were measured separately. The

φ = 4 • V π • D 2 • u 2 ( )
In this instance D stands for the impeller Diameter and u for the circumferential speed of the impeller at its outlet. Differences in the flow coefficient between the LPF and the HPF are the result of the different average volume flow rates over the fans due to an increasing density within the system. All results regarding the pressure increase are displayed as the dimensionless head coefficient Ψ:

Ψ = ∆p t ρ/2 • u 2 2 ( )
The head coefficient establishes the link between the pressure increase and the circumferential speed at the impellers outlet. To compare the pressure increases of each stage, the dimensionless pressure ratio Π is introduced and defined as:

Π = ∆p outlet ∆p inlet ( )
The setup for measuring the bearing loads is shown in Figure in a sectional view. The base is a two-pole, asynchronous motor with several modifications. The red sections are the positions used for measuring the axial and radial forces. The drive shaft is movable in the z-direction. This axial displacement is measured by a load cell at the end of the shaft. The clutch decouples the rotation of the shaft from the axial displacement. The radial forces are measured by the radial load cell in the immediate vicinity of the impeller. In on the direction. For the presented results, the measurement direction of the load cell is oriented to the x-axis in order to disregard gravitation forces. In order to evaluate the axial forces independently of the impeller diameter D, the medium density ρ and the circumferential speed at the outer diameter of the impeller u 2 the dimensionless axial thrust coefficient α is introduced [ ]:

α = F ax ρ/2 • u 2 2 /4 • D 2 • π ( )
As with the axial thrust coefficient, the dimensionless radial thrust coefficient with regard to density ρ, the circumferential speed u 2 , the impeller diameter D and the impeller width (including the thickness of shroud and hub) b 2,tot is introduced for the radial forces [ ]:

k Ru = F r ρ/2 • u 2 2 • D • b 2,tot ( )
Fluctuations in the volume flow rate and pressure are caused by the rotating impeller with its finite blade number. These fluctuations influence on the force measurements. Radial and axial forces can thus be divided into a static or average load and a fluctuating load. One further point to note is that the radial forces occur simultaneously at the A-bearing and the B-bearing. Consequently, the presented loads acting in a radial direction only contain the loads on the A-bearing. Due to an equilibrium of forces in a radial direction, the other forces can be calculated using the distances of the bearings or impeller respectively. The relative total pressure equates to the total pressure increase at the different flow coefficients normalized by the (measured) maximum total pressure increase of the single fan:

∆p t,rel = ∆p t ∆p t,max,single f an ( ) So the maximum relative total pressure increase of the single fan is 1, 0. To calculate the pressure increase of a serial arrangement with two identical fans, the pressure increase at each volume flow rate is multiplied with the amount of .

RESULTS

The results are presented in two sections. The first part shows the performance, namely the pressure-to-volume characteristics of the serially arranged fans in contrast to the single fan.

The second section focuses on the resulting bearing loads in axial and radial direction.

.

Performance of the Serial Arrangement

According to the relevant literature the performance of the serial arrangement is calculated by doubling the total pressure at each volume flow rate. This is illustrated in Figure with the relative total pressure against the flow coefficient.

As can be seen from the graph, the maximum of the total pressure increase of the single fan is measured at a flow coefficient of approx. 0.04. The calculation for the serial arrangement based on the single fan results in the black dashed line with the maximum relative total pressure of 2.0. The maximum flow coefficient of the single fan and the The measured results of the serially arranged fans differ in some respects from the described calculated values. At low flow coefficients the measured relative total pressures are lower than the calculated ones. A possible reason could be, that at low flow coefficients the density of the air gets lower at the inlet of the high pressure fan (compared to the density at the inlet of the low pressure fan) through the high temperatures. The density could also be the reason for the higher relative total pressure at flow coefficients between 0.01 and 0.10.

Between a flow coefficient of 0.01 and 0.10 the measured pressure increases are higher until both curves are relatively equal at approx. 0.01. At flow coefficients higher than 0.01 the measured pressure increases are higher than the calculated ones and an increase in the maximum flow coefficient of over 0.16 is possible. The measurement indicates an increase in the maximum flow coefficient. Another noteworthy point is that the difference of the relative total pressure between the calculation and the measurement tends to increase with a rising flow coefficient at flow coefficients > 0.10. This leads to the conclusion that the maximum flow coefficient is higher than expected. The maximum attainable increase in the flow rate of the serial arrangement is thus also higher. To clarify these points, the two fans in the serial arrangement are considered separately. The chart in Figure shows at low flow coefficients, the head coefficient in the LPF is higher than in the HPF. At lower flow coefficients, the increase in pressure is therefore delivered by the LPF. With an increasing flow coefficient there is a reversal of the trend. At higher flow coefficients, the head coefficient of the HPF is higher. Thus, despite the geometrically identical structure of both fans, their behavior differs. These results indicate that both the LPF and the HPF only deliver the same head coefficient at one operating point. However, the blue line of the single fan and the red line of the LPF are almost identical. Therefore, the LPF in the serial arrangement has a comparable behavior to a single fan. It stands to reason that the first stage will not be affected by the downstream stages in any serial arrangement. Based on the presented data, the LPF of the serial arrangement has a behavior comparable to the single fan in terms of performance. Therefore, the state of the art for designing single fans can be transferred to the LPF. It is also important to note that the identified effects on the characteristic curve of the serial arrangement (see Figure

) are caused by the HPF. Even if the basic concept of the serially arranged fans is to increase the available pressure an increase in the maximum volume flow rate must be observed. A different head coefficient results in a different pressure increase at a particular operating point. These changes resulting in with the uneven pressure distribution and concerning the bearing loads are evaluated in the next section.

. Bearing Loads

In order to calculate the bearing loads, durability, shaft stresses and deflection, the radial and axial forces primarily generated by the impeller must be known. Simplified in the form of the non-uniform pressure distribution at the outlet of the rotating impeller, radial forces are created. As already mentioned, there is a fluctuation in the volume flow rate and the pressure increase caused by the rotating impeller with its finite blade number. These fluctuations influence force measurement. Therefore, the radial and axial forces can be divided into a static or average load and a fluctuating load. In this case, the mean over three revolutions of the impeller in terms of the average load is presented in of almost zero. With increasing flow coefficients, there is a nearly linear increase of the radial thrust coefficient. The LPF starts with a radial thrust coefficient of less than zero and exhibits a milder increase compared to the single fan. The radial load on the high-pressure fan is comparable to that of the single fan. In summary, the maximum radial thrust coefficients are at the highest flow coefficients. The amplitude of the radial thrust coefficient of the single fan, the LPF and the HPF in Figure fluctuates against the flow coefficient. It can be assumed that the amplitude of the radial thrust coefficient follows a parabolic trend with its apex in the range of the maximum efficiency. The changes incurred by the serial arrangement in terms of radial forces are unremarkable. In absolute terms, the radial loads of the single fan and the serial arrangement are not that high compared to the loads in axial direction.

Axial forces occur due to uneven pressure distributions on the surface of the hub and the shroud of the impeller and due to the momentum caused by the deflection of the flow. with an axial thrust coefficients of 0.25 and increases slightly until an φ of 0.07 is reached. After this, there is a jump to α = 0.28 at φ = 0.08 followed by a decrease. The characteristic of the LPF (red) is comparable to that of the single fan: a moderate increase in the axial thrust coefficient at small flow coefficients with a peak of α = 0.28 at φ = 0.08 and a subsequent reduction. The HPF deviates from this behavior. Only the trend of the curve is similar to the single fan and LPF. Starting with an axial thrust coefficient of 0.20 at low flow coefficients, there is a slight drop with increasing flow coefficients. What is remarkable is the sudden increase of the axial thrust coefficient with its peak in the same region of the flow coefficient like those of the single fan and the LPF, but with a lower maximum of α = 0.25. With further increasing flow coefficients, the axial thrust coefficient of the HPF converges with the others. There is a reduction in the axial forces in the HPF in serial arrangements. This reduction is almost 20 % at certain operating points. At flow coefficients higher than 0.10 the loads begin to converge with each other. The influence of the momentum is also visible with regard to the trend of the axial thrust coefficient at high flow coefficients. For the fluctuating loads presented in Figure , HPF, the static pressure ratio is examined against the flow coefficient in Figure . Here, the pressure ratio is defined as the quotient of the static pressure at the outlet to the static pressure of the inlet on the relevant fan. Obviously, at low flow coefficients, the pressure ratio by the LPF is higher. With increasing flow coefficients the static pressure increase of the LPF drops more rapidly than that of the HPF. There is an intersection of the curves at a flow coefficient of approximately 0.075. Up to this point, the pressure ratio of the HPF is higher. The lower pressure ratio at flow coefficients of 0.075 0 0.05 0.10 0.15 1.00 [-]

LPF HPF

Figure .

Pressure ratio over flow coefficient of LPF and HPF and below could contribute to the lower axial thrust coefficient in these flow coefficient ranges. There is a difference in the radial and axial forces in a serial arrangement depending on the fan being observed. The LPF has loads in a radial and axial direction comparable to the single fan. The loads that the HPF has in a radial direction are also comparable to the single fan. However, the axial loads in the HPF are lower compared to the single fan over a wide range of flow coefficients.

. CONCLUSION

In the current paper, the influence of two radial fans in a serial arrangement on the characteristic curve and the resulting effects on the bearing forces were examined. Both fans are geometrically identical centrifugal fans which were directly connected to each other and run at the same speed. The theoretically described behavior of the serial arrangement differs from the measured results. The measured pressure increase at defined volume flow rates is equivalent to or even higher than the calculated values (according to literature) which can be obtained by doubling the total pressure of the single fan. Another remarkable point is that the pressure not only increases, but there is also a higher flow coefficient or a higher volume flow rate respectively. The two fans in the serial arrangement, the low-pressure fan (LPF) and the high-pressure fan (HPF) show deviations in their characteristics as well as their bearing loads. The LPF displays behavior comparable to that of a single fan in its characteristics as well as in the occurring loads. This statement could be validated for the characteristic curve and the loads. It stands to reason that the first stage will not be affected by the downstream stages in any serial arrangement. The state of the art of designing single fans can therefore be transferred to the low-pressure stage of a multi-stage arrangement. Despite the identical geometric parameters of the HPF, its characteristic curve has a different shape. At low flow coefficients, the head coefficient is lower in comparison to the LPF. The pressure increase is thus unevenly composed. Only at one operating point is the pressure increase equally distributed across both fans. For lower flow coefficients, the main portion of the pressure increase is delivered by the LPF. With increasing flow coefficients, however, there is a reversal of the trend. At higher volume flow rates, the ratio of the pressure increase is higher for the HPF. So, despite the geometrically equal design of both fans, their behavior differs. A different pressure distribution also leads to a change in the loads of the bearings. The axial loads at the HPF in particular decrease in comparison to the single fan. Thus, for future design work on serial arrangements, the focus must be specifically on the HPF.

. OUTLOOK

These results show that the theoretically described characteristic curve differs from the measured one. The HPF in particular exhibits a significant change in behavior in terms of the characteristic curve and bearing loads. The focus will therefore be placed on the HPF to clarify these changes via computational fluid dynamics in combination with additional experiments. [ ] J.F. Guelich. Centrifugal pumps. Springer, edition, .

[ ] EN ISO -: Measurement of fluid flow by means of pressure differential devices inserted in circular crosssection conduits running full -Part : General principles and requirements.

. .

  Figure . Product portfolio analysis [ ]

Figure .

 . Figure . Characteristic curves of a single fan and serially arranged fans [ ]

Figure .

 . Figure . Test rig sketch and Cordier diagram

Figure .

 . Figure . Measurement of bearing loads in sectional view [ ]

Figure .

 . Figure . Relative total pressure against flow coefficient

Figure .

 . Figure . Head coefficient against flow coefficient of LPF and HPF

Figure .

 . Figure . Averaged radial thrust coefficient against flow coefficient of single fan, LPF and HPF

FigureFigure . Figure .

 .. Figure . Amplitude of radial thrust coefficient against flow coefficient of single fan, LPF and HPF

Figure .

 . Figure . Amplitude of axial thrust coefficient against flow coefficient of single fan, LPF and HPF

  [ ] EN ISO -:Measurement of fluid flow by means of pressure differential devices in circular cross-section conduits running full -Part : Orifice plates. . [ ] ISO : Measurement of fluid flow -Evaluation of uncertainties.

  Institute for Engineering Design and Industrial Design, University of Stu gart, Germany *Corresponding author: Martin.Nudischer@iktd.uni-stu gart.de

	NOMENCLATURE
	b 2,tot	Impeller width at the outlet	m
	D	Diameter Impeller	m
	F ax	Axial force	N
	F r	Radial force	N
	k Ru	Radial thrust coefficient	-
	n	Speed	s -1
	u 2	Circumferential speed at	m/s
		the outlet of the impeller
	V α	Volume flow rate Axial thrust coefficient	m 3 /s -
	∆p	Pressure increase	Pa
	δ opt	Diameter number	-
	Π st at. Pressure ratio	-
	ρ	Density	kg/m 3

  Over the decades, a lot of scientific research has been conducted into single centrifugal fans, see for example [ ], [ ], [ ] or [ ]. Due to the feasible strength of the materials used (usually steel or aluminum), the peripheral speed is limited

	to approximately	m/s for welded or riveted sheet-metal
	impellers [ ]. The available pressure increase and volume
	flow rate are thus limited with a single centrifugal fan. Sub-
	stituting the impeller material for high-strength materials
	enables higher peripheral speeds but also leads to excessive
	production costs.
	Only short theoretical descriptions are given in the literature
	for centrifugal fans with regard to the serial connection of
	two or more fans [ ], [ ]. Considerable changes incurred by
	coupling two or more fans cannot be identified. Eck theo-
	retically described the resulting characteristic curve of more
	than one centrifugal fan serially arranged through super-
	position of the single characteristic curves [ ]. The total
	pressure increase at a defined volume flow rate is therefore
	multiplied by the number of coupled fans. Trautmann et al.

Table .

 . Technical data of the fan

	Volume flow rate Total pressure increase 14800 Pa 38 m 3 /min
	Engine power	11 kW
	Rotational speed	11850 min -1
	Impeller diameter	240 mm

  Table the technical data of the used load cells are summed up. The measurement direction of the load cell is dependent Table . Technical data of the load cells

		Radial load cell
	Rated capacity	±1500 N
	Tolerance of sensitivity < ±0.2%
	Accuracy class	±0.3% (F nominal )
		Axial load cell
	Rated capacity	±2000 N
	Tolerance of sensitivity < ±0.2%
	Accuracy class	±0.1% (F nominal )

Table .

 . List of uncertainties

	Static pressure Volume flow rate ±2.0 % ±1.4 % Axial force ±2.6 % Radial force ±3.6 %

follows: Axial forces are positive in a positive z-direction and radial forces are positive in a negative x-axis.

  ] A. Krake and K. Fiedler, editors. Laser--Focus Measurements and Flow Visualisation within the Rotating Rotor and the Side Channel of a Side Channel Compressor, Hawaii, .[ ] A. Starck, A. Mühlbauer, and C. Kramer. Handbook of Thermoprocessing Technologies. Vulkan Verlag, edition, . Theoretische und experimentelle Untersuchung des Pumpens von lufttechnischen Anlagen mit Radialventilatoren. University of Karlsruhe, . Das instationäre Betriebsverhalten ungeregelter und geregelter lufttechnischer Anlagen mit einem, zwei in Reihe oder zwei parallelgeschalteten Ventilatoren, volume . . ] C.E. Brennen. Hydrodynamics and cavitation of pumps. In d'Agostino L. and Salvetti M. V., editors, Fluid Dynamics of Cavitation and Cavitating Turbopumps, pages -. Springer, .

	[ ] Elektror	[online].	Available	at:
	http://www.elektror.de/en/products. Accessed / /	.
	[ ] EN ISO	: Industrial Fans -Performance Testing Using
	Standardized Airways.	.

[ [ ] W. Bohl. Ventilatoren: Berechnung Konstruktion, Versuch, Betrieb. Vogel Fachbuch, edition, . [ ] L. Bommes, J. Fricke, and R. Grundmann. Ventilatoren. Vulkan Verlag, edition, . [ ] T. Carolus. Ventilatoren: Aerodynamischer Entwurf, Schallvorhersage, Konstruktion. Vieweg+Teubner Verlag, edition, . [ ] B. Eck. Ventilatoren: Entwurf und Betrieb der Radial-, Axial-und Querstromventilatoren. Springer, Berlin, edition, . [ ] T. Carolus. [ ] G. Trautmann and K. Jun. Leittechnische Maßnahmen zur Vermeidung unzulässiger Rauchgasdruecke bei Nachruestung von Dampferzeugern mit REA und DE-NOX. Jahrbuch der Dampferzeugungstechnik, : -, . [ ] H.-U. Banzhaf. Anlagespezifische Fragen bei Ventilatoren in Reihenschaltung. Brennstoff, Wärme, Kraft, : -, . [ ] B. Schulze Dieckhoff. [ [ ] S. Recker. Ganzheitliche Optimierung schnelldrehender Hochdruckventilatoren. Institutsbericht IKTD Universitaet Stuttgart, , .