Moustafa Nakechbandi
email: moustafa.nakechbandi@univ-lehavre.fr

Jean-Yves Colin Litis
email: jean-yves.colin@univ-lehavre.fr

AN ALGORITHM AND SOME NUMERICAL EXPERIMENTS FOR THE SCHEDULING OF TASKS WITH FAULT-TOLERANCY CONSTRAINTS ON HETEROGENEOUS SYSTEMS

Keywords: DAG, Scheduling with communication, Fault tolerant, Heterogeneous systems

In this paper, we propose an efficient scheduling algorithm for problems in which tasks with precedence constraints and communication delays have to be scheduled on an heterogeneous distributed system with an one fault hypothesis. Based on an extension of the Critical-Path Method CPM/PERT, our algorithm combines an optimal schedule with some additional tasks duplication, to provide fault-tolerance. Backup copies are not established for tasks that have already more than one original copy. The result is a schedule in polynomial time that is optimal when there is no failure, and is a good resilient schedule in the case of one server failure. We finally compare the optimal solutions with the resilient solutions found by this algorithm on several semi-random DAGs.

I. INTRODUCTION

Heterogeneous distributed systems have been increasingly used for scientific and commercial applications. Recent examples of such applications include Automated Document Factories (ADF) in banking environments where several hundred thousands documents are produced each day on networks of several multiprocessors servers. Or high performance Data Mining (DM) systems [START_REF] Palmerini | On performance of data mining: from algorithms to management systems for data exploration[END_REF]) that need to process very large data collections using very time-consuming algorithms. Or Grid Computing systems [START_REF] Ruffner | An Enhanced Application Model for Scheduling in Grid Environments[END_REF][START_REF] Venugopal | A Grid Task Broker for Scheduling Distributed Data-Oriented Applications on Global Grids[END_REF]) such as Computational Grids which focus primarily on very computationallyintensive operations, or Data Grids which control the sharing and management of large amounts of distributed data.

However, efficiently using these heterogeneous systems is a hard problem, because the general problem of optimally scheduling tasks is NP-complete, even when there are no communication delays [START_REF] Kwok | Static scheduling algorithms for allocating directed task graphs to multiprocessors[END_REF]Ahmad 1999, Garey and[START_REF] Garey | Computers and Intractability, a Guide to the Theory of NP-Completeness[END_REF]. When the application tasks can be represented by Directed Acyclic Graphs (DAGs), many dynamic scheduling algorithms have been devised. For some examples, see [START_REF] Maheswaran | A Dynamic matching and scheduling algorithm for heterogeneous computing systems[END_REF][START_REF] Iverson | Dynamic, Competitive Scheduling of Multible DAGs in a Distributes Heterogeneous Environment[END_REF][START_REF] Chen | Distributed Dynamic Scheduling of Composite Tasks on Grid Computing Systems[END_REF]. Also, several static algorithms for scheduling DAGs in meta-computing systems are described in [START_REF] Colin | Scheduling with Small Communication Delays and Task Duplication[END_REF][START_REF] Topcuoglu | Task scheduling algorithms for heterogeneous processors[END_REF][START_REF] Alhusaini | A Unified Resource Scheduling Framework for Heterogeneous, Computing Environments[END_REF][START_REF] Kwok | Static scheduling algorithms for allocating directed task graphs to multiprocessors[END_REF]. Most of them suppose that tasks compete for limited processor resources, and thus these algorithms are mostly heuristics. Problems with fault tolerant aspects are less studied. Reliable execution of a set of tasks is usually achieved by task duplication and backup copies [START_REF] Qin | A Novel Fault-tolerant Scheduling Algorithm for Precedence Constrained Tasks in Real-Time Heterogeneous Systems[END_REF][START_REF] Randell | System structure for software fault-tolerance[END_REF][START_REF] Chen | N-version programming: a fault tolerant approach to reliability of software operation[END_REF][START_REF] Girault | A scheduling heuristics for distributed real-time embedded systems tolerant to processor and communication media failures[END_REF].

A very classical and useful tool to study static scheduling problems with DAG is the Critical Path Method (also known as CPM, or PERT method, or CPM/PERT) [START_REF] Maheswaran | A Dynamic matching and scheduling algorithm for heterogeneous computing systems[END_REF]. Using a relaxation of the constraint on the number of available processors, this method gives results such as a lower bound on the execution time (or makespan) of the application and lower bounds on the execution dates of all tasks of the DAG. Because of the relaxation, tasks can be executed as soon as possible. Improvements and limits of this method to distributed systems with communications delays may be found in [START_REF] Colin | Scheduling with Small Communication Delays and Task Duplication[END_REF][START_REF] Colin | A multi-valued DAG model and an optimal PERT-like Algorithm for the Distribution of Applications on Heterogeneous, Computing Systems[END_REF][START_REF] Nakechbandi | Bounding the makespan of best pre-scheduling of task graphs with fixed communication delays and random execution times on a virtual distributed system[END_REF], for example. In [START_REF] Colin | A multi-valued DAG model and an optimal PERT-like Algorithm for the Distribution of Applications on Heterogeneous, Computing Systems[END_REF], we studied the problem of scheduling the tasks of a DAG on the servers of an heterogeneous system. There, the relaxation used in CPM/PERT was replaced by the dual relaxation that each server has no constraint on the number of tasks it can simultaneously process. That is, each server can simultaneously process a non limited number of tasks without loss of performances. Our goal was to compute a lower bound on the execution time of a realistic solution, and compute lower bounds on the execution dates of all tasks of the DAG. In [START_REF] Nakechbandi | An efficient fault-tolerant scheduling algorithm for precedence constrained tasks in heterogeneous distributed systems[END_REF], we further supposed that one server (and at most one) could suffer from a crash fault. The algorithm presented there improved on the one presented in [START_REF] Colin | A multi-valued DAG model and an optimal PERT-like Algorithm for the Distribution of Applications on Heterogeneous, Computing Systems[END_REF]) by adding backup copies to the optimal solution build.

The solution we propose now is simpler than the one presented in [START_REF] Nakechbandi | An efficient fault-tolerant scheduling algorithm for precedence constrained tasks in heterogeneous distributed systems[END_REF]. Additionally, we present some numerical experiments and simulation results. This rest of this paper is divided into four main parts. In the first one, we present the problem, and in the second one, we present a solution to the problem. In the third part, we make some numerical experiments using randomly generated tasks graphs, comparing the optimal solutions with the resilient solutions found by this algorithm. Finally, in the fourth part, we discuss the advantages and disadvantages of the proposed solution.

II. THE CENTRAL PROBLEM

The Distributed Servers System

We call Distributed Servers System (DSS) a virtual set of geographically distributed, multi-users, heterogeneous or not, servers. Therefore, a DSS has the following properties: first, the processing time of a task on a DSS may vary from a server to another. The processing time of each task on each server is supposed known. Second, although it may be possible that some servers of a DSS are potentially able to execute all the tasks of an application, it may also be possible in some applications that some tasks may not be executed by all servers. In a DSS problem, we suppose that the needs of each task of an application are known, and that at least one server of the DSS may process it.

The classical CPM/PERT relaxation of the number of processors, is replaced in the DSS problem with the dual relaxation that each server has no constraint on the number of tasks it can simultaneously process. Thus we suppose that the concurrent executions of some tasks of the application on a server have a negligible effect on the processing time of any other task of the application on the same server.

The transmission delay of a result between two tasks depends on the tasks and on their respective sites. The communication delay between two tasks executed on the same server is supposed equal to 0.

Directed Acyclic Graph

An application is decomposed into a set of indivisible tasks that have to be processed. A task may need data or results from other tasks to fulfil its function and then send its results to other tasks. The transfers of data between the tasks introduce dependencies between them. The resulting dependencies form a Directed Acyclic Graph. Because the servers are not necessarily identical, the processing time of a given task can vary from one server to the next. Furthermore, the duration of the transfer of a result on the network cannot be ignored. This communication delay is function of the size of the data to be transferred and of the transmission speed that the network can provide between the involved servers. Note that if two dependent tasks are processed themselves on the same server, this communication delay is considered to be 0.

The central scheduling problem P on a Distributed Server System, is represented therefore by the following parameters:

• a set of servers, noted Σ = {σ 1 , ..., σ s }, interconnected by a network, • a set of the tasks of the application, noted I = {1,..., n}, to be executed on Σ. The execution of task i, i ∈ I, on server σ r , σ r ∈ Σ, is noted i/σ r . The subset of the servers able to process task i is noted Σ i , and may be different from Σ,

• the processing times of each task i on a server σ r is a positive value noted r i σ π / . The set of processing times of a given task i on all servers of Σ is noted Π i (Σ). r i σ π / = ∞ means that the task i cannot be executed by the server σ r .

• a set of the transmissions between the tasks of the application, noted U. The transmission of a result of an task i, i ∈ I, toward a task j, j ∈ I, is noted (i, j). It is supposed in the following that the tasks are numbered so that if (i, j) ∈ U, then i < j, • the communication delays of the transmission of the result (i, j) for a task i processed by server σ r toward a task j processed by server σ p is a positive value noted

p r j i c σ σ / , /
. The set of all possible communication delays of the transmission of the result of task i, toward task j is noted ∆ i,j (Σ). Note that a zero in ∆ i,j (Σ) mean that i and j are on the same server, i.e. = ∞ means that either task i cannot be executed by server σ r , or task j cannot be executed by server σ p , or both.

Let Π (Σ) = U I i∈ Π i (Σ) be the set of all processing times of the tasks of P on Σ.

Let ∆ (Σ) = U U j i ∈ ) , (
∆ i,j (Σ) be the set of all communication delays of transmissions (i, j) on Σ.

The central scheduling problem P on a distributed servers system DSS can be modelled by a multi-valued DAG G = {I, U, Π(Σ), ∆(Σ)}. In this case we note P={G,

Σ}.

Example 1 : Figure 2 presents an example of DAG.

Fig. 2. Example of DAG : the Πi vector on a node is the vector of the processing time of task i on the various servers, and ∆i,j on an arc is the communication delays matrix between the two tasks depending on the

servers that process them.

On this example, if we have 4 servers {σ 1 , σ 2 , σ 3 , σ 4 }

and if Π 1 = ( 3 3 3 3, ∞, 2, ∞), then 1 / 1 σ π =3. And 2 / 1 σ π = ∞ ,
meaning that server σ 2 cannot execute task 2 etc.

On the same example, communications from task 1 to task 2 are given by matrix ∆ 1,2 in Fig. 3.

σ1 σ σ σ σ2 2 2 2 σ3 σ4 σ1 0 3 2 ∞ σ2 ∞ ∞ ∞ ∞ σ σ σ σ3 3 3 3 2 3 3 3 3 0 ∞ σ4 ∞ ∞ ∞ ∞ Fig. 3. Example of communication delays matrix ∆1,2
between task 1 and task 2. In the matrix of Fig. 3, one can see that if task 1 is processed on server σ 3 and task 2 is processed on server σ 2 , then c 1/σ3, 2/σ2 = 3.

Definition of a feasible solution

We note PRED(i), the set of the predecessors of task i in

G: { } ) , ( et / ) PRED( U i k I k k i ∈ ∈ =
And we note SUCC(i), the set of the successors of task

i in G: { } ) , ( et / ) SUCC( U j i I j j i ∈ ∈ =
A feasible solution S for the problem P is a subset of executions { i/σ r , i∈I } with the following properties:

• each task i of the application is executed at least once on at least one server σ r of Σ i ,

• to each task i of the application executed by a server σ r of Σ i , is associated one positive execution date

r i t σ / ,
• for each execution of a task i on a server σ r , such that PRED(i) ≠ ∅, there is at least an execution of a task k, k ∈PRED(i), on a server σ p , σ p ∈ Σ κ , that can transmit its result to server σ r before the execution date

r i t σ / .
The last condition, also known as the Generalized Precedence Constraint (GPC) [START_REF] Colin | A multi-valued DAG model and an optimal PERT-like Algorithm for the Distribution of Applications on Heterogeneous, Computing Systems[END_REF], can be expressed more formally as:

     + + ≥ Σ ∈ ∃ ∈ ∀ ≥ ∈ ∀ r p p p r r i k k k i k p i r c t t i k t S i σ σ σ σ σ σ π σ σ / , / / / / / / ), PRED( 0 / else ) PRED( if ∅ = i
It means that if a communication must be done between two scheduled tasks, there is at least one execution of the first task on a server with enough delay between the end of this task and the beginning of the second one for the communication to take place. A feasible solution S for the problem P is therefore a set of executions i/σ r of all i tasks, i ∈ I, scheduled at their dates r i t σ / , and verifying the Generalised Precedence Constraints GPC. Note that, in a feasible solution, several servers may simultaneously or not execute the same task. This may be useful to generate less communications.

All the executed tasks in this feasible solution, however, must respect the Generalized Dependence Constraints.

Optimality Condition

Let T be the total processing time of an application (also known as the makespan of the application) in a feasible solution S, with T defined as:

) ( max / / r r r i i S i/ t T σ σ σ π + = ∈ A feasible solution S* of the problem P modelled by a DAG G = {I, U, Π(Σ), ∆(Σ)} is optimal if its total
processing time T* is minimal. That is, it does not exist any feasible solution S with a total processing time T such that T < T*.

III. THE DSS_1FAULT ALGORITHM

The algorithm proposed here, named DSS_1FAULT, has two phases: the first one is for the scheduling of original copies where we use the DSS-OPT algorithm [START_REF] Colin | A multi-valued DAG model and an optimal PERT-like Algorithm for the Distribution of Applications on Heterogeneous, Computing Systems[END_REF]) and the second one is for adding and scheduling additional backups copies when necessary.

Scheduling the original copies

We schedule original copies of tasks in our algorithm with the DSS-OPT algorithm [START_REF] Colin | A multi-valued DAG model and an optimal PERT-like Algorithm for the Distribution of Applications on Heterogeneous, Computing Systems[END_REF]. The DSS-OPT algorithm is an extension of CPM/PERT algorithms type to the distributed servers problem. In its first phase, it computes the earliest feasible execution date of each task on every server, and in its second phase it builds a feasible solution (without server fault) starting from the end of the graph with the help of the earliest dates computed in the first phase.

Let P be a DSS scheduling problem, and let G = {I, U, Π(Σ), ∆(Σ)} be its DAG.

One can first note that there is an optimal trivial solution to this DSS scheduling problem. In this trivial solution, all possible tasks are executed on all possible servers, as soon as possible, and their results are then broadcasted to all others servers. This is an obvious waste of processing power and communication resources, however, and something as optimal, but less wasteful in terms of used resources, is usually needed.

The first phase of the DSS_OPT routine, DSS_LWB(), goes from the initial tasks to the final ones, computing along the way the earliest feasible execution dates The second phase of the DSS_OPT routine determines, for every task i that does not have any successor in G, i.e. task i is a "leaf" or final task, the execution i/σ r ending at the earliest possible date r / i r σ . If several executions of task i end at the same smallest date r / i b σ , one is chosen, arbitrarily or using other criteria of convenience, and kept in the solution. Then, for each kept execution i/σ r that has at least one predecessor in the application, the subset L i of the executions of its predecessors that satisfy GPC(i/σ r ) is established. This subset of executions of predecessors of i contains at least an execution of each of its predecessors in G. One execution k/σ p of every predecessor task k of task i is chosen in the subset, arbitrarily or using other criteria of convenience, and kept in the solution. It is executed at its earliest possible date p / k b σ . The examination of the predecessors is pursued in a recursive manner until the studied tasks do not present any predecessors in G.

Adding backup copies

The ADD_BACKUP_COPIES routine starts from tasks without any predecessors, similarly to DSS_LWB(), and proceed from there to the end of the DAG. First, if there is currently only one copy of a given task, it determines what is the worst possible delay it may encounter if a failure occurs on another server, while satisfying its GPC. It also determines the fastest server (not considering the server executing the only current copy of this task in the current solution) able to execute this task, and adds a backup copy on this server to the solution, again considering the worst possible delay resulting from this failure, while satisfying the GPC of this copy. Else the task has already several copies in the optimal solution, and the routine determines for each original copy of this task, what is the worst possible delay it may encounter if a failure occurs on another server, while satisfying its GPC.

The complete DSS_1FAULT algorithm is the following: By adding backup copies using ADD_BACKUP_COPIES we get the following fault-tolerance scheduling (Fig. 5.):

Input: G = {I, U, Π(Σ), ∆(Σ)
) ( min max / ) ( SUCC / r i r i i i r T σ σ Σ ∈ ∀ ∅ = ∀ = for all tasks i such that SUCC(i) = ∅ do i L ← { i/σ r / σ r ∈ Σ ι and T r r i ≤ σ / } i/σ r ← keepOnefrom( i L ) schedule (i/σ
(i) ≠ ∅ then for each task k such that k ∈ PRED(i) do r i k L σ / ← { k/σ q / σ p ∈ Σ κ and r r p p p i i k k k b c b σ σ σ σ σ π / / , / / / ≤ + + } k/σ q ← keepOneFrom( r i k L σ / ) schedule (k/σ q ) end

Fig. 5. Gantt chart given by DSS_1FAULT

Because the computed execution time of each task on each server is its earliest execution time on this server, and because only the copy with the earliest ending time, of each task without any successor, is used in the solution calculated by DSS_OPT() , and finally because all other copies are used only if they ensure that the final copies receives their data in time else they are not used, it follows that the feasible solution computed by DSS_OPT is optimal in execution time for the problem without server failure.

Lemma 1: The feasible solution calculated by the DSS_OPT algorithm is optimal if there is no server failure.

Because the copies in the DSS_1FAULT solution coming from the DSS_OPT solution will not be delayed if there is no server failure, and because additional backup will not be used in this case, then we have: Theorem 1: The solution calculated by DSS_1FAULT is optimal if there is no server failure. Also, in the final solution computed by DSS_1FAULT(), each task of the DAG has at least two copies (coming from the DSS_OPT() routine), or one copy (coming from the DSS_OPT() routine) and one backup copy (build by the ADD_BACKUP_COPY() routine) , always executed on different servers.

Furthermore, the execution date of each backup copy and the delayed execution date of each original copy coming from DSS_OPT is always evaluated by ADD_BACKUP_COPIES() taking into account the delayed execution dates of the copies and the execution dates of the backups copies of each predecessor, using the worst possible case of failure of a predecessor, we have: Theorem 2: The solution calculated by DSS_1FAULT is feasible if there is at most one server failure.

The most computationally intensive part of DSS_OPT() is the first part DSS_LWB(). In this part, for each task i, for each server executing i, for each predecessor j of i, for each server executing j, a small computation is done. Thus the complexity of DSS_LWB() is Ο(n 2 s 2 ), where n is the number of tasks in P, and s is the number of servers in DSS. Thus, the complexity of the DSS_OPT() algorithm is Ο(n 2 s 2 ).

Similarly, in ADD_BACKUP_COPIES(), for each task i, for each copy of task i (at most one copy per server), for each predecessor j of i, for each copy of j (at most one per server), one small computation is done. Thus the complexity of ADD_BACKUP_COPIES() is bounded by Ο(n 2 s 2 ), where n is the number of tasks in P, and s is the number of servers in DSS. Thus we have: Theorem 3: The complexity of the DSS_1FAULT algorithm is Ο(n 2 s 2 ).

IV. PERFORMANCE EVALUATION

To evaluate DSS_1FAULT, we have compared the fault tolerant solutions it generated on some classical problems and DAG to optimal solutions without fault tolerancy. These numerical experiments were done using simulations on three different kinds of graphs. The first one is a simple, semi-random, one level 'fork-join' DAG (see Fig. 6. a.), with limited parallelism. The second one is a regular simple two-dimensional grid DAG (see Fig. 6. b.), exhibited by some numerical applications, with lot of parallelism and very local communications. The last one is the "butterfly" DAG (see Fig6. c.) present in applications such as the FFT or shuffles algorithms, again with lot of parallelism, but a more complex communication pattern. The servers performances are independent random values for each task of the DAG, and so is each communication delay.

Fork-Join DAG

As expected, this kind of DAG does show a very limited parallelism. On the Gantt chart example in Fig. 7, one can distinguish the original copies of the tasks on the left part of each server's simulated activity chart, and the added backup copies on the right part.

2-Dimensional grid DAG

This highly parallel DAG is much more efficiently executed on the servers. On the Gantt chart example in Fig. 8, the original copies of the tasks are grouped in the left part of each server's simulated activity chart, with the added backup copies spread more widely on the rights part. Although this is not clearly visible on the black and white figure, some added backup copies of the earliest tasks of the DG are present in the mist of the original copies.

Butterfly DAG

The Butterfly is a highly parallel DAG too. On the Gantt chart example of Fig. 9, the original copies are clearly distinguishable one the left, as bands. Because of the random nature of the server's performances, these bands tend to become fuzzier as time passes, however. The backup copies are scheduled later on the right part of the chart.

Makespan with and without backup copies

In all three kinds of DAGs, it is found that the makespan average with backup copies is between 1.5 (usually) and 2 (at most) times the makespan without backup copies. For example, in the Butterfly DAGs, we obtained the following figure (Fig . 10). In this simulation the number of tasks varies from 10 to 1200 tasks and we have the average over 50 random DAGs.

V. ANALYSIS

This algorithm has two advantages: • when no server fails, the DSS-1FAULT's solution is optimal as it uses the optimal solution computed by DSS-OPT. • when there is a failure of one server, the DSS-1FAULT's solution is certain to finish correctly, because every tasks has two or more scheduled copies on different servers in the final solution. If more than one fault occur, the solution may still finish, but there is no guaranty there. We do not establish backup copies for tasks which have already two or more original copies from the DSS-OPT algorithm scheduling to limit tasks duplication and processor. It also gives indications on the sensibility of an application to one server failure when compared to the solution without any server failure, because the makespan in the presence of one failure is a worst case analysis.

The model of failure, as it features at most 1 crash, may seem poor. However, if the probability of any failure is very low, and the probabilities of failure are independent, then the probability of two failures will be much smaller indeed. Furthermore, the algorithm may be extended to 2 or more failures, by using two or more backup copies per task. The efficiency of this kind of solution to the "kfailures" problem is not investigated, yet.

Finally, the solution solved by this new algorithm uses the classical CPM/PERT relaxation, namely that an unbounded number of tasks may be processed on each server in parallel without any effect on the tasks' processing time, in the same sense that the CPM/PERT method do not consider resources constraints in order to get earliest execution dates. This relaxation is not far from the reality, if each server is a multiprocessors architecture. Or if each server is a time-shared, multi-users system with a permanent heavy load coming from other applications, and the tasks of an application on each server represent a negligible additional load. In other cases, the same way these CPM/PERT results are used in some real-life systems as the priority values of tasks in some list-scheduling algorithms, the result found by our algorithm may be used as the first step of a list scheduling algorithm, in which the earliest execution dates of primary and backup copies are used as priority values to schedule these copies on the servers of a real-life system.

VI. CONCLUSION AND FUTURE WORKS

In this paper, we have proposed a polynomial scheduling algorithm in which tasks with precedence constraints and communication delays have to be scheduled on an heterogeneous distributed system environment with one fault hypothesis. To provide a faulttolerant capability, we employed primary and backup copies. But no backup copies were established for tasks which have more than one primary copy.

The result have been a schedule in polynomial time that gives earliest execution dates to copies of tasks when there is no failure, and is a good resilient schedule in the case of one failure. Performance evaluation on some DAGs gave an increase in case of one server failure in makespan of 1.5 to 2 times the optimal makespan without server failure.

The execution dates of the original and backup copies may be used as priority values for list scheduling algorithm in cases of real-life, limited resources, and systems.

In our future work, we intend to study the same problem with sub-networks failures. Also, we intend to consider the problem of non permanent failures of servers. Finally, we want to consider the problem of the partial failure of one server, in which one server is not completely down but loses the ability to execute some tasks and keeps the ability to execute at least one other task.

Fig. 1 :

 1 Fig. 1: Example of Distributed Servers System with the list of the executable services for each server.

  for all possible executions i/σ r of each task i of problem P.

  each task i where PRED(i) = ∅ do for each server σ r such that σ r ∈ Σ i do there is a non marked task i such that all its predecessors k in G are marked do for each server σ r such that σ r ∈ Σ i do

Find

  Fig. 4. Gantt chart given by DSS_OPT. The fact that task 3 is executed at the same time that task 2 on server σ 3 comes from the CPM/PERT relaxation.

Fig. 6 .

 6 Fig. 6. Three different kind of graphs

Fig. 7 .

 7 Fig. 7. Gantt chart for Fork-Join DAG A. These lines represent the execution of the originals copies

Fig. 8 .

 8 Fig. 8. Gantt chart for grid DAG A. These lines represent the execution of the originals copies B. These lines represent the execution of the backups copies.

Fig. 9 .

 9 Fig.9. Gantt chart for Butterfly DAG A. These lines represent the execution of the originals copies B. These lines represent the execution of the backups copies.

Fig. 10 .

 10 Fig.10. Makespan average for Butterfly DAGs A. Makespan without backup, B. Makespan with backup,

  }

	Output: A feasible solution with backup copies
	DSS_1FAULT ()
	DSS_OPT()	// first phase
	ADD_BACKUP_COPIES()	// second phase
	end DSS_1FAULT
	DSS_OPT()
	DSS_LWB ()

for end if end schedule keepOneFrom(L i

  ≠ σ i be the fastest server able to execute i Execute a new backup copy of i on σ r at date 0 Let σ r ≠ σ i be the fastest server able to execute i Execute a backup copy of i on σ r taking only into account the delayed execution dates of the copies and backups of each predecessor of i to verify the GPC else // i has at least two copies scheduled, on // separate servers, of course // compute the delayed execution date of the // copy of task i on each server, if the failure is

	end if
	mark (i)
	end for
	while there is a non marked task i such that all its
	predecessors k in G are marked do
	if i has only one copy scheduled then
	Let σ i be the server executing the copy of i
	// First compute the delayed execution date of
	// task i on this server, if the failure is on an
	// another server
	find the delayed execution date of the copy of i
	on σ i taking only into account the delayed
	execution dates of the copies and backups of
	each predecessor of i to verify the GPC
	// Second compute one backup copy on the
	// fastest server left, if failure is on server of
	// primary

) return an execution i/σ r of task i in the list of the executions L i . end keepOneFrom.

ADD_BACKUP_COPIES()

for each task i such that PRED(i) = ∅ do if i has only one copy scheduled then //compute one backup on the fastest server left, if // failure is on server of this copy

Let σ r // on an another server for each server σ i executing a copy of i do


Valerii Sopin
email: vvs@myself.com

Existence of a quadratic polynomial, which represents infinitely many prime numbers

Keywords: Bunyakovsky's conjecture, Landau's problems, complete and subcomplete sequences, prime numbers, p-adic analysis, Fermat's theorem on sums of two squares, Goldbach'

No single case of Bunyakovsky's conjecture for degree greater than one has been proved, although numerical evidence in higher degree is consistent with the conjecture. In this paper we overcome such misfortune (using Friedlander-Iwaniec theorem, Fermat's theorem on sums of two squares and Brahmagupta-Fibonacci Identity, Bezout's lemma and a connection to SL(2, Z) and Hyperbolic Prime Number Theorem).

Introduction

Since the number of all primes is infinite we have the question of whether the number of primes, which for example are contained in the form n 2 + 1 (the 4th Landau problem, the 5th Hardy-Littlewood conjecture), is also infinite. As the 4th Landau problem (it goes back to Euler, who observed in a letter to Goldbach in 1752 that n 2 + 1 is often prime for n up to 1500) demonstrates, particular problems of this type have long been of interest; however, the first formulation of a precise conjecture had to wait until Viktor Bunyakovsky in 1857 [START_REF] Bunyakovsky | Sur les diviseurs numeriques invariables des fonctions rationnelles entieres[END_REF]: Suppose f (x) is a non-constant polynomial with integer coefficients and positive leading coefficient. Moreover, suppose that f (x) is irreducible in Z[x] (Z is for integers) and that there is no prime q which divides f (n) for every positive integer value of n (the numbers f (n) should be relatively prime). Then f (n) is prime for infinitely many positive integer values of n.

Bunyakovsky's conjecture was later generalized to a finite family of polynomials by Schinzel, who with Sierpinski [START_REF] Schinzel | Sur certaines hypotheses concernant les nombres premiers[END_REF] gave several applications to elementary number theory.

In the case of a single linear polynomial, Bunyakovsky's conjecture amounts to Dirichlet's 1837 theorem on primes in progressions. But Dirichlet's theorem has remained the only proven case of it.

Probably the main problem with the conjecture is the lack of good reformulations of its conditions in case of degree higher than 1. This leads to the idea of consideration not one polynomial, but aggregation of polynomials in the following way:

Conjecture. If the leading coefficient of a polynomial f (x) with integer coefficients is positive, then there exists integer c such that f (N) + c contains infinitely many primes.

It is helpful to keep in mind next picture: every integer point (x, y) on coordinate plane represents tuple {x, f (x) + y}. Notice that for any fixed n f (n) + c (c is any integer) contains all prime numbers, as it covers range of arithmetic progression x + 1. Moreover, Hilbert's irreducibility theorem guarantees that the polynomial f (x) + c is irreducible for almost every c. In addition, it is worthy to highlight that the Conjecture is obvious for a linear polynomial and it can be used to give a simple proof of Dirichlet's theorem on arithmetic progressions.

In this paper we show existence of a quadratic polynomial, range of which contains infinitely many prime numbers. On the whole, we prove the first case of Bunyakovsky's conjecture for degree greater than 1.

Main result

To obtain the pronouncement we are going to use Fermat's Theorem on sums of two squares.

Fermat's Theorem on sums of two squares

Fermat's Theorem on sums of two squares states that if p = 4k + 1 is a prime number, it can be expressed as the sum of two squares. Euler succeeded in proving Fermat's theorem on sums of two squares in 1749 with such propositions:

(i) if a number which is a sum of two squares is divisible by a prime which is a sum of two squares, then the quotient is a sum of two squares;

(ii) if a number which can be written as a sum of two squares is divisible by a number which is not a sum of two squares, then the quotient has a factor which is not a sum of two squares.

In addition, formula (Sum of Squares Function, see [START_REF] Beiler | Recreations in the Theory of Numbers: The Queen of Mathematics Entertains[END_REF]) for number of representations of a natural t = 2 a0 q 2a1 1 . . . q 2ar r p b1 1 . . . p bs s ,

where the q i are primes of the form 4k + 3 and the p j are primes of the form 4k + 1, as the sum of two squares, ignoring order and signs, is

r 2 (t) =      0 if any a i is a half-integer 1 2 B
if all a i are integer and B is even

1 2 (B -(-1) a0
) if all a i are integer and B is odd

B = (b 1 + 1)(b 2 + 1) . . . (b s + 1)
Accordingly, representation as sum of two squares is unique for any prime 4k + 1 and any prime 4k + 3 is not sum of two squares.

Remark 1. B. M. Bredihin proved the infinitude of primes of the form x 2 + y 2 + 1.

Theorem

Theorem. There exists a quadratic polynomial, which represents infinitely many prime numbers.

Proof. If range of x 2 + 1 or x 2 + (x -1) 2 contains infinitely many primes, we have done. Otherwise, assume the contrary: there exists natural N that for any natural number n > N none of the numbers n 2 + 1, n 2 + (n -1) 2 are primes.

If p = 4k+1 is a prime number, then there must be natural m such that m 2 +1 is divisible by p (we can see this by Euler's criterion or via Lagrange's approach with quadratic forms). Thus, for p = 4k + 1 > N we have that

p = t 2 + l 2 , m 2 + 1 = (v 2 + w 2 )(t 2 + l 2 ),
where t, l > 1, t -l > 1 and v 2 + w 2 > 1 by the assumption. Diophantus/Brahmagupta-Fibonacci Identity says that the product of two sums of two squares is itself a sum of two squares. Namely,

(a 2 + b 2 )(c 2 + d 2 ) = (ac -bd) 2 + (ad + bc) 2 = (ac + bd) 2 + (ad -bc) 2 ,
where, without loss of generality, a ≥ b ≥ 0, c ≥ d ≥ 0.

Observe:

(1) If c, d > 1, then ad + bc, ac + bd > 1.

( It is well-known that there exist infinitely many integers m such that m 2 + 1 is either prime or the product of two primes. Thus, by the assumption it is the product of two primes (Sum of two squares theorem is notable here and the case (4k + 3) 2 = m 2 + 1 is impossible). We notice that conditions a, b, c, d > 1 with c -d > 1 and a -b > 1 for p 1 = a 2 + b 2 and p 2 = c 2 + d 2 , where m 2 + 1 = p 1 p 2 , mean that, according to the observations on the previous page, ad-bc = 1. Note here that r 2 (m 2 +1 = p 1 p 2 ) = 2 then and Brahmagupta-Fibonacci Identity yields all cases for Sum of Squares Function.

If one of the primes is always bounded, then from Bezout's lemma we know that all solutions of ad -bc = 1 for fixed natural a, b can be represented in the form (d + kb, c + ka), where k is an arbitrary integer, i.e. infinite number of primes under consideration satisfy finitely many certain linear patterns. In particular, there exist a, b with corresponding c, d that

(d + kb) 2 + (c + ka) 2 = (c 2 + d 2 ) + (2ac + 2db)k + (a 2 + b 2 )k 2
are primes for infinitely many k.

Otherwise, according to the first lines of the page we always have relation ad -bc = 1 between two primes p 1 = a 2 + b 2 and p 2 = c 2 + d 2 in all sufficiently big m 2 + 1 = p 1 p 2 , i.e. we have come to an analogy of the so-called Hyperbolic Prime Number Theorem, which established upper and lower bounds for the number of primes p = a 2 + b 2 + c 2 + d 2 up to x with the hyperbolic condition ad -bc = 1.

Remark 2. Semisubgroups of SL 2 (Z)? Notice that the Semigroup Membership problem is decidable for 2-by-2 integer matrices with non-zero determinant (I. Potapov and P. Semukhin).

SL 2 (Z) = {2 -by -2 integer matrices with determinant 1} 1 0 1 1 n = 1 0 n 1 1 0 1 1 n-1 1 1 0 1 = 1 1 n -1 n M atrices 1
Remark 3. Euler's 6k + 1 theorem states that every prime of the form 6k + 1 can be written in the form x 2 + 3y 2 with x and y positive integers (every prime number other than 2 and 3 is of the form 6k ± 1). A prime dividing n 2 + n + 1 can only be 3 or of the form 6k + 1.

Remark 4. All primes appear as factor on the three polynomials: x 2 + 1, 2x 2 + 1, 2x 2 -1. All primes of 2x 2 -1 are of the form 8k + 1 or 8k + 7.

Generalization

Waring's problem

Waring's problem asks whether each natural number k has an associated positive integer s such that every natural number is the sum of at most s natural numbers raised to the power k. For every k let g(k) denote the minimum number s of kth powers of naturals needed to represent all positive integers. It has been determined that g(k) < ∞. For example, every natural is sum of at most 4 squares, 9 cubes, or 19 fourth powers. Remark 5. From [START_REF] Panaitopol | On the Representation of Natural Numbers as Sums of Squares[END_REF] it is known that for each odd natural n there exist non-negative integers x, y, z such that n = x 2 + y 2 + 2z 2 . Moreover, it is not possible to represent any natural number in the manner ax 2 + by 2 with 1 ≤ a ≤ b. Remark 6. Notice that all numbers of the form 2(4 m ), 6(4 m ) or 14(4 m ) cannot be represented as a sum of four non-zero squares. Moreover, 2(4 m ), 6(4 m ) and 14(4 m ) have only one representation as a sum of four squares. Is 2(4 m ) + 2 a sum of two primes in the form 4k + 1?

Goldbach-Linnik type and Waring-Goldbach type problems

Goldbach's conjecture is one of the oldest and best-known unsolved problems in number theory and all of mathematics. It states that every even natural number greater than 2 is the sum of two prime numbers, as Goldbach was following the now-abandoned convention of considering 1 to be a prime number.

It is known that the set of even integers that are not the sum of two primes has density zero. Olivier Ramare showed that every even number is sum of at most six primes and Terence Tao indicated that every odd integer is a sum of at most five primes. Moreover, it was exhibited by Yuri Linnik that every sufficiently large integer can be represented as the sum of two primes and K powers of two, where K is an absolute number (K = 13 is acceptable).

Furthermore, all sufficiently large odd integers can be represented in the form n = p 1 +p 2 2 +p 3 3 +p 4 4 +p 5 5 , where p i , i = 1, . . . , 5 are primes. Additionally, every sufficiently large integer congruent to 14 modulo 240 may be written as the sum of 14 fourth powers of prime numbers. And every sufficiently large odd integer may be written as the sum of 21 fifth powers of prime numbers. Besides, every sufficient large odd integer is a sum of one prime, two squares of primes and 31 powers of two.

Remark 7. It is worth bringing up here Yitang Zhang and Polymath Project collaborative result: there are infinitely many prime gaps (a prime gap is the difference between two successive prime numbers) that do not exceed some constant (246 can be taken).

Cubes and Fourth Powers

Which primes are sums of two cubes? The prime 2 and primes of the form 3x 2 -3x + 1 for some integer x. However, there are infinitely many primes of the form x 3 + 2y 3 . Moreover, Ω(n 0.9 ) of the numbers from 1 to n have representations as sums of three cubes of non-negative integers. And it is conjectured that every sufficiently large natural number can be represented as a sum of four cubes of natural numbers. So, is the Conjecture true for f (x) = x 3 ?

It is important problem to determine whether there are infinitely many prime numbers, which are represented in the form p = a 4 + b 4 . However, we know that Ω(n 0.8 ) of the numbers from 1 to n have representations as sums of four fourth powers.

Remark 8. Notice that some probabilistic models suggest that the sums of four fourth powers, and more generally sums of k perfect kth powers for k ≥ 3, should have positive natural density. In particular the gaps between these numbers are conjectured to have bounded average size.

3 Useful observations

Complete and subcomplete sequences

Given a sequence S of positive integers, let P (S) be the set of numbers which can be represented as the sum of a finite subsequence of S. Then S is complete if N \ P (S) is finite, and subcomplete if there are positive integers m and n such that {mk + n, k ∈ N} ⊂ P (S).

Roth & Szekeres [START_REF] Roth | Some asymptotic formulae in the theory of partitions[END_REF] and Graham [START_REF] Graham | Complete sequences of polynomial values[END_REF] showed that if a polynomial f with real coefficients maps integers to integers, a necessary and sufficient condition for the range of f to be complete is that the leading coefficient is positive and for any prime p there exist an integer n such that p does not divide f (n).

Any integer-valued polynomial f (x) = a 0 + a 1 x + • • • + a d x d can be written in the basis of binomial coefficient polynomials f (x) = â0 + â1

x 1 + • • • + âd x d
, where each âi is an integer. As claimed in [START_REF] Roth | Some asymptotic formulae in the theory of partitions[END_REF] [START_REF] Graham | Complete sequences of polynomial values[END_REF] the range of a polynomial f is complete if and only if âd > 0 and gcd(â 0 , â1 , . . . , âd ) = 1.

Additionally, if a (rational) polynomial maps positive integers to positive integers, its image is subcomplete [START_REF] Roth | Some asymptotic formulae in the theory of partitions[END_REF][6] [START_REF] Burr | On the completeness of sequences of perturbed polynomial values[END_REF]. Szemeredi and Vu [START_REF] Szemeredi | Finite and infinite arithmetic progressions in sumsets[END_REF] showed that there is a constant c such that the following holds: any

increasing sequence A = a 1 < a 2 < a 3 < . . . satisfying A(n) ≥ cn 1/2 (with A(n) = |{a i ∈ A : a i ≤ n}|)
for all n is subcomplete. A subcomplete sequence S is complete if and only if for all m P (S) has an element in all residue classes mod m [5].

The p-adic analysis

Every polynomial with integer coefficients is 1-lipschitz mapping over the ring of p-adic integers Z p . It turns out that any 1-lipschitz measure-preserving transformation on Z p is an isometry which induces permutations on all residue rings Z/p k Z, k = 1, 2, . . . , and vice versa. That's why measure-preserving polynomials with integer coefficients on Z p are called permutation polynomials. See for details [START_REF] Anashin | Applied algebraic dynamics[END_REF] [START_REF] Rivest | Permutation Polynomials Modulo 2 w[END_REF]. Even more, a polynomial f with integer coefficients induces a measure-preserving transformation on the ring of p-adic integer Z p if and only if the mapping x → f (x) mod p 2 is a bijective transformation on residue ring Z/p 2 Z [START_REF] Anashin | Applied algebraic dynamics[END_REF].

A permutation polynomial f is called ergodic polynomial on Z p if the mapping x → f (x) mod p k is transitive modulo p k for all k = 1, 2, . . . , where transitivity modulo p k means that this permutation x → f (x) mod p k on residue ring Z/p k Z has only one cycle of length p k . It is known that f is ergodic if and only if it is transitive modulo p 3 , see [START_REF] Anashin | Applied algebraic dynamics[END_REF].

Obviously, if a polynomial f is a permutation polynomial on Z p , then polynomial f + c is also a permutation polynomial on Z p for any integer c. However, a linear polynomial from the Dirichlet's theorem is not always measure-preserving, see [START_REF] Anashin | Applied algebraic dynamics[END_REF]. And gcd(â 0 , â1 , . . . , âd ) = 1 is not always true for permutation polynomials, see [START_REF] Anashin | Applied algebraic dynamics[END_REF].

Additionally, it seems that general case, an ordinary polynomials, can be considered as deviation from permutation polynomials by meaning that not whole Z p , but some ball of non-zero radius in it is in the range. Paper [START_REF] Zieve | Polynomial Mappings mod p n[END_REF] studies the topological structure of polynomial mappings over the ring of p-adic integers.

Friedlander-Iwaniec theorem

It states that there are infinitely many prime numbers of the form a 2 + b 4 , see [START_REF] Friedlander | Using a parity-sensitive sieve to count prime values of a polynomial[END_REF]. Moreover, [START_REF] Iwaniec | Almost-primes represented by quadratic polynomials[END_REF] proves that there exist infinitely many integers n such that n 2 + 1 is either prime or the product of two primes. The last implies that one prime factor p of m 2 + 1 is strictly smaller than m, and therefore also divisor of (the usually much smaller) m2 + 1, where m = m mod p.

For an arbitrary polynomial f that is irreducible and doesn't have a fixed prime divisor, one can say that f represents infinitely many elements with at most 1 + deg(f ) factors.

  ) If c, d > 1 and b = 0, then ac -bd = ac, ad -bc = ad > 1. (3) If b = 0, c, d > 1 and (c -d) > 1, then ac -bd = b(c -d) + (a -b)c > 1. (4) If b = 0, then |ad -bc| = |b(d -c) + (a -b)d| can possibly be one. Notice that c, d is relatively prime. Given relatively prime integers c, d, there are integers a, b such that ad -bc = 1, as there are integers s, t that 1 = sc + td (the Euclidean algorithm). Thus, a := t, b := -s.

  For any p = c 2 + d 2 there exist pair of integers a, b such that ad -bc = 1, i.e.On the other hand, a matrix a b c d corresponds to (ac + bd) 2 + (ad -bc) 2 = (ac + bd) 2 + 1.


Davide Beneventi
email: davide.beneventi@pagora.grenoble-inp.fr

Andrea Graziano
email: arch.a.graziano@gmail.com

Gioia Furia
email: gioia.furia@cea.fr

Mathieu Lionel Charpin
email: charpin@alprobotic.fr

Mathieu Jean-Louis
email: mathieu.jean-louis@rtone.fr

Mathieu Jean

Optimization of trajectory processing algorithms to print 3D circuit boards using piezo ink

Keywords: Six-axis robot, Jetting, Printed electronics, Lines morphology, Conductivity

Introduction

Over the last decade, additive manufacturing (AM) progressively evolved from a laboratory curiosity to a ubiquitous manufacturing approach able to process a wide variety of materials (i.e. from metal powders to cell-inoculated hydrogels) for the production of both structural and functional objects [START_REF] Beaman | Additive Manufacturing Review: Early Past to Current Practice[END_REF]. Moreover, the straightforward link between the digital model and the final physical product allows including additive manufacturing among the production processes driving the Industry 4.0 revolution [START_REF] Dilberoglu | The Role of Additive Manufacturing in the Era of Industry 4.0[END_REF]. Despite the wide variety of technologies implemented on 3D printers, most of commercial processes are based on the use of cartesian robots and on the layer by layer construction of the 3D object. The intrinsic limitations of this approach, such as the generation of anisotropic objects and the impossibility to deposit material along out-ofplane freeform paths [START_REF] Bhatt | Expanding of additive manufacturing through use of robotics technologies: A survey[END_REF], led to an ever-increasing interest for poly articulated industrial robots.

Owing to i) their ability to follow freeform trajectories while keeping a controlled orientation of the dispensing tool, ii) the ease of scale-up, and iii) high payload/positioning and repeatability, 6-axis arms have been recently identified as ideal manipulators to carry dispensing tools and have been integrated in both conventional AM technologies [START_REF] Bhatt | Expanding of additive manufacturing through use of robotics technologies: A survey[END_REF]- [START_REF] Shah | Emerging trends in robotic aided additive manufacturing[END_REF] and hybrid additive-subtractive processes [START_REF] Li | A novel 6-axis hybrid additive-subtractive manufacturing process: Design and case studies[END_REF].

Recent studies on cold spray [START_REF] Nault | Multi-axis tool path optimization and deposition modeling for cold spray additive manufacturing[END_REF] and hydrogel [START_REF] Wulle | Multi-axis 3D printing of gelatin methacryloyl hydrogels on a non-planar surface obtained from magnetic resonance imaging[END_REF] deposition on non-planar substrates highlighted that the redundant degree of freedom of 6 axis robots imposes to tightly correlate the material dispensing tool characteristics, the geometry of the printing substrate and the tool head velocity in order to depose high-precision patterns. In order to develop a general method for the generation of freeform printing trajectories, algorithms to optimize the generation of printing paths accounting for robot joints kinematics [START_REF] Xie | A potential field based multi-axis printing path generation algorithm[END_REF] and to fit 2D trajectories on regular 3D substrates [START_REF] Hosseini Jafari | Surface parameterization and trajectory generation on regular surfaces with application in robot-guided deposition printing[END_REF] with minimum deformation have been recently proposed demonstrating that, despite the inherent complexity of 6 axis manipulators, a simplified approach can be envisaged.

As observed earlier for structural and large-scale applications [START_REF] Urhal | Robot assisted additive manufacturing: A review[END_REF], the use of poly articulated manipulators is triggering an ever-increasing interest for the AM of 3D electronics.

Leaving aside the manipulator technology, most of additive manufacturing processes have been developed to process one-two materials with the aim of manufacturing structural objects. Only after the pioneering works of Macdonald [START_REF] Macdonald | 3D Printing for the Rapid Prototyping of Structural Electronics[END_REF] and Lewis [START_REF] Lewis | Device fabrication: Three-dimensional printed electronics[END_REF], where conventional FDM 3D printing was coupled to needle-dispensing of conductive inks in order to embed conductive paths into 3D plastic objects, a wide variety of multi materials AM processes have been developed [START_REF] Lu | Additive manufacturing frontier: 3D printing electronics[END_REF]- [START_REF] Tan | 3D printed electronics: Processes, materials and future trends[END_REF] allowing the simultaneous fabrication of both structural and electronic components. Electronic circuits surface integration has been also proposed as alternative to their bulk integration into 3D objects and high line density circuits on complex 3D shapes have been manufactured by lased direct structuring (LDS) of conductive paths on injection-moulded plastics [START_REF] Heininger | Manufacturing of molded interconnect devices from prototyping to mass production with direct structuring[END_REF]- [START_REF] Huttunen | Multilayer plastic substrate for electronics[END_REF]. Regardless of the high precision and high conductivity of surface circuits obtained by LDS, this manufacturing technology relies on the use of a specific thermoplastic composite and a sequence of laser etching and metallization steps [START_REF] Heininger | Manufacturing of molded interconnect devices from prototyping to mass production with direct structuring[END_REF] which preclude its use to generic 3D objects. The direct printing of conductive and dielectric inks by piezo jetting [START_REF] Ledesma-Fernandez | High viscosity jetting of conductive and dielectric pastes for printed electronics[END_REF], ink jet and aerosol spray [START_REF] Seifert | Additive manufacturing technologies compared: morphology of deposits of silver ink using inkjet and aerosol jet printing[END_REF] or needle dispensing [START_REF] Tricot | Fabrication of 3D conductive circuits: print quality evaluation of a direct ink writing process[END_REF] appear as a viable alternative to LDS. Indeed, recent studies demonstrated that the integration of piezo or pneumatic jetting printheads on a 6-axis robot can be effectively used to print silver tracks on a 3D paper substrate [START_REF] Furia | Use of a 6-axis robot and ink piezo-jetting to print conductive paths on 3D objects. Printed circuit geometry, and conductivity predictive model[END_REF], RF antennas on PC-ABS objects [START_REF] Delfaut | Study of feasibility for 2D and 3D innovative jet printed RF devices[END_REF], [START_REF] Le Dam | Study of a 3D FSS-Based Printed Radome Demonstrator for a Patch Antenna at 3.6 GHz[END_REF] and large-scale composites [START_REF] Riviere | Les enjeux de l'impression directe de fonctions radiofréquence sur pièces 3D, Mecaplastronic Connection[END_REF]. Whereas, pneumatic needle dispensing has been used to print conductive tracks during the 3D printing of a PLA helicoidal structure on a double robot bench [START_REF] Bao | 3D Structural electronics via multidirectional robot 3D printing[END_REF].

Despite first works demonstrate that 6 axis robots can be effectively used to print electronic circuits on the surface of generic 3D objects, the control of the full process (i.e. from the definition of the 3D object, circuit design and the identification/optimization of trajectory processing algorithms to high precision circuit board printing) and its impact on the final print quality remains rather unexplored. The aim of this work is to investigate a 3D circuits printing process based on the 3D scan of a generic object followed by the transfer of 2D circuits on the object surface and the generation of smooth 3D printing trajectories using the combination of mesh-and normals-smoothing algorithms. In order to validate the viability of the proposed process, the trajectory generation and printing conditions optimized on model geometries (semi spherical substrate and simple circuit design) were used to design and print a hybrid 3D circuit board including surface mount discrete components, a USB PCB and a NFC antenna.

Materials and methods

Printing cell

Cell setup. A inhouse designed robotic cell [START_REF] Furia | Use of a 6-axis robot and ink piezo-jetting to print conductive paths on 3D objects. Printed circuit geometry, and conductivity predictive model[END_REF], composed of a 6-axis robot (Staubli TX2 60, 670 mm range and 0.02 mm repeatability) equipped with a high-speed laser profilometer (Keyence, LJ-V7080/LJ-X8000A) and a piezo jetting dispensing valve (Vermes MDV 3200A) mounting a 0.1 mm nozzle was wired and mounted on a low vibration bench (AlpRobotic, France). The laser profilometer was used to scan the surface of 3D objects at speeds of 40 or 50 mm/s and sampling frequencies of 200 and 1000 Hz, corresponding to sampling periods of 200 and 50 µm, respectively. The planar scan trajectory was made of 10 cm length parallel lines with a step of 3 cm. Height data were retrieved as .csv files and combined. The scan tool centre position (TCP) was used to obtain the coordinates of each sampled point in the robot reference cartesian system (World). Sampling periods (d) of 0.2 and 0.05 mm were used to obtain low or high definition 3D scans which were processed in order to generate the printing trajectory and to analyse the morphology of printed paths, respectively.

The piezo jetting non-contact technology was selected for this study since: i) it allows processing high density inks with limited/no flowing when deposed on inclined/vertical surfaces, ii) the printhead can be operated at a distance of 1-3 mm from the substrate thus preventing the collision of the printhead with the substrate in case of inaccurate calibration of the robot-printhead kinematic chain. The jetting valve has a fully programmable duty cyle. In this study the tappet-rising, -opening and -falling time had constant values (of 0.5, 2 and 0.2 ms, respectively), whereas, the delay time between each dispensing cycle was varied between 0.1 and 1000 ms in order to adapt the overall jetting frequency to the printing speed (i.e. the jetting valve TCP speed) and the desired path width using a semi-empirical model developed in a previous study [START_REF] Furia | Use of a 6-axis robot and ink piezo-jetting to print conductive paths on 3D objects. Printed circuit geometry, and conductivity predictive model[END_REF]. The protocol to generate an accurate and smooth printhead trajectory and to synchronize the ink jetting frequency with the printing speed was based on four sequential operations: i) Meshing and smoothing points obtained with the laser profilometer. In order to remove the measurement noise, points were meshed using rectangular cells and the Paul Bourke (PB)

smoothing method [START_REF] Bourke | Surface Relaxation and Smoothing, Geometry, Surfaces, Curves, Polyhedra: Notes on Polygon and Meshes[END_REF], i.e.:

𝑃 = 𝑃 + ∑ 𝑃 -𝑃 (1)
where, 𝑃 is the smoothed point, n is the number of points 𝑃 surrounding the pristine point 𝑃 (i.e. n = 4 for rectangular cells). Two major constrains dictated the number of loops when applying the PB smoothing method:

the offset of smoothed points with respect to the pristine ones, had to be < 200 µm in order to limit the scattering of the smoothed mesh from the real geometry.

the computational load, ie. the number of basic operations necessary to obtain a sufficient noise reduction both in terms of point position and orientation of normals to the 3D mesh, which for the PB method is 𝑁𝑏𝑜 = 30 .

ii) Transferring 2D paths on the 3D mesh with minimum geometric deformation. Planar circuits were transferred on 3D substrates, using quasiisothermic mesh parametrization [START_REF] Sechelmann | Quasiisothermic Mesh Layout[END_REF] and planarization (GrassHopper, Kangaroo). 2D circuit curves are discretized with 0.2 mm resolution, positioned, and then projected on the planarized mesh in order to extract their points' positions in terms of mesh face number and mesh parameter (position on the mesh face). Curves are then reconstructed as polylines onto the original 3D mesh using the mesh face number and points parameters. Length variations between 2D and 3D homologue curves are then calculated to evaluate the deformation induced by the transformation.

iii) Smoothing normals to the 3D mesh of points in the 3D circuit in order to remove the noise generated by the discrete mesh. The normal orientation alignment algorithm (NA) aims at smoothing the robotic toolpath in terms of cinematic rotations, allowing a more homogeneous traveling speed. The iterative algorithm for each trajectory point i checks the adjacent vectors and does the average, 𝜗 , of the three (the previous one 𝜃 , the considered one 𝜃 , and the next one 𝜃 ):

𝜗 = (2).
In this way, the direction of each vector tends to align to the contiguous vectors reducing normals oscillations due to the scan mesh geometry and noise. As for mesh smoothing, the angular scattering among normals of neighbouring points and the computational load, which for the NA is 𝑁𝑏𝑜 = 9

, where used to constrain the number of loops. iv) To synchronize the jetting frequency with the effective printhead speed/position along the trajectory. Since the complex kinematic and the presence of singularities does not allow poly articulated robots to keep a constant speed along 3D paths, the TCP speed was retrieved during a blank run along the printing trajectory. According to a semi-empirical model described in a previous work [START_REF] Furia | Use of a 6-axis robot and ink piezo-jetting to print conductive paths on 3D objects. Printed circuit geometry, and conductivity predictive model[END_REF], the average speed (𝑉 ), the diameter of individual ink drops ejected by the printhead (0 = 0.280 mm in the tested conditions) and the target line width (w) were used to calculate the jetting frequency (𝐹 ) and the drop number (𝑛 ) necessary to print the path using the equations:

𝑤 =  • •( ) - ( ) -1 -(1 -𝐶) (3)
and

𝑛 = 𝐿 (4),
where

𝐶 =   (5)
is the linear coverage ratio between adjacent drops, L is the length of the printed line and 𝑆 = .

. • .

is the experimental drop spreading coefficient for the tested Ag ink as determined measuring the diameter of 2-10 superimposed drops printed on the ABS substrate (Fig. 2.).

Printing of conductive inks.

Effect of smoothing and circuit design on the printing speed. The effect of circuit design and smoothing on the printing speed was quantified using an ABS 3D printed semi sphere (45 mmradius, 0.1 mm layer) polished with a P600 sandpaper and the two planar circuits shown in Fig. 3.a and b. Right angles in the circuits were filleted with radii ranging between 0 and 4 mm and they were printed with a jetting frequency of 89 Hz corresponding to a drop coverage ratio of 66%. The circuit shown in Fig. 3.c was used to evaluate the print accuracy with a reduced pitch and to limit the line width it was printed using a jetting frequency lowered to 75 Hz (i.e. a drop coverage ratio of 60%). Owing to the absence of tiny details on the 3D substrate and to limit the number of points to process, the 3D scan sampling period (d) of 0.2 mm was under sampled to 0.4 mm and the region of interest was restricted to the print area, i.e. a spherical dome with a base diameter of 50 mm. Fig. 4.a shows that he planar scan with a fixed sampling period d generates rectangular cells with size (𝑙 ) and normal angle deviation between neighbouring cells (𝜗 ) depending on the point distance from the centre of the half-sphere:

𝑙 = 𝑅 𝐶𝑜𝑠 -𝐶𝑜𝑠 (6)
𝜗 = (7)
where 𝐷 = 𝑖 𝑑 is the distance from the half sphere centre and i the number of considered periods.

According to eq. ( 6) and ( 7) the minimum angle variation with a sampling step of 0.4 mm is 0.51 deg (Fig. 4.b) and, in the region of interest, it has an average value of 0.53 deg which was considered as the target value in order to obtain a noise-free trajectory after PB smoothing.

According to the lower computational load of the normal smoothing algorithm NA, 3D circuits were discretized with a period of 0.2 mm, corresponding to a target angle variation of 0.25 deg for a noise-free trajectory. The effect of the PB and NA smoothing algorithms on the motion of the printhead TCP along the 3D trajectory was evaluated by varying the number of loops (i.e from 0 to 40 and 0 to 300 for PB and NA, respectively) and retrieving from the robot controller the effective TCP speed and orientation with a sampling period of 4 ms.

3D circuits printing. After the optimization of the printing trajectories, a commercial silver ink (Henkel, Loctite Edag 418SS) was used to print the 3D circuits whose geometry, before and after dying in an air circulation oven (30 min at 90 °C), was recorded with high definition 3D scans (0.05 mm sampling period). Dry circuits were also characterized with an USB optical microscope (Dino-Lite).

Printed circuits geometry (i.e. length and interline pitch) and line section were compared with the theoretical ones, i.e.: the pristine 2D circuit geometry and the line sectional area (A) and height (h) calculated as

𝐴 = (8)
and

ℎ = • • (9)
where Vd can be either the wet ink volume

𝑉 = 𝑉 • 𝑛 (10)
or the dry ink volume

𝑉 = 𝑉 • 𝑛 • (11),
w is the line width and nd is the number of drops deposed for a line length L, V0 is the drop unit volume (i.e. 11.5 10 -3 mm 3 in the tested jetting conditions), l and d are the liquid-and the dry-ink densities, i.e. 1.98 and 3.73 g/cm 3 , respectively.

The electric resistance of printed lines was measures using a digital multimeter (Fluke 113)

and it was compared to the corresponding theoretical value calculated using the Ohm law and equations 4, 8, 11

𝑅 = Ω L ( 12
)
where  is the nominal resistivity of the Ag ink used in this study (i.e. ≤ 7.5 10 -7 Ohm m), L is the line length (i.e. 35 mm), 𝑉 is the average TCP speed and 𝐹 is the jetting frequency and 𝑉 is the dry ink volume calculated with equation [START_REF] Macdonald | 3D Printing for the Rapid Prototyping of Structural Electronics[END_REF].

The shape fidelity and electric resistance of 3D printed circuits and the computational load and processing time were used to evaluate the overall reliability of the developed process.

3D printed circuit board.

A NFC card reader composed by a hybrid circuit board including a USB PCB, surface mount discrete (SMD) components (resistors, capacitors, transistors, leds), printed interconnections and a NFC antenna, was designed using a commercial PCB designer software (Altium Designer), Fig. 5a. The 3D substrate was constituted by a structural part of the card reader (Fig. 5b) and included a flat region for the USB PCB mount and a connection and a large area concave regions with10 and 90 mm bending radius, respectively.

Optimized printing conditions and geometrical constrains imposed by the 6-axis robot were used to re-design and print the 3D circuit. SMD components were manually placed on the 3D circuit before ink curing (30 min at 90 °C). Line resistance and components characteristics were measured on the final 3D printed board.

Results and discussion

3.1 Effect of PB and NA algorithms on normals alignment.

Noise generated by the laser profilometer during height measurement induced slight deviations in the overall shape of the reconstructed dome. However, the impact on the orientation of mesh cells was not negligible and, as shown in Fig. 6, the average angle variation between neighbouring cells attained 1.8 deg with peaks of 4-5 deg instead of the theoretical 0.53 deg.

With the print configuration used in this study (i.e. piezo jetting with a printhead aligned to the normal of the substrate and vertical distance of 2 mm) these uneven angle variations can induce offsets up to 0.16 mm in the target ink deposition point and irregular movements of the 6-axis robot.

Fig. 6 shows that PB smoothing of the pristine mesh efficiently removes the noise and after 20 loops attains the minimum theoretical average angle variation that can be obtained on a 90 mm diameter dome and 0.4 mm sampling period with a minimum dispersion.

Similarly, the application of the NA algorithm to the normals of points lying on the 3D trajectory produced a drop in the average angle which reached the theoretical value of 0.25 deg and minimum dispersion after 20-30 loops (Fig. 7). In the conditions tested in this study, the NA algorithm quickly converged and had a limited computational load compared to PB mesh smoothing (Fig. Owing to the heavy computational load for the relatively small mesh, the large sampling period used in this study (i.e. 3800 mm 2 and 0.4 mm) and the constrain to limit the offset of sampled points position below 0.2 mm, the PB smoothing algorithm was limited to 10 loops. In these conditions the average angle variation was close to the theoretical limit, however individual angles had quite scattered values with standard deviation of 0.2 deg and maximal amplitude of 0.8 deg (Fig. 6.a). The NA algorithm was then used to attain the minimal angle variation of 0.25 deg applying 30 smoothing loops. 3.2 Effect of smoothing and circuit design on TCP speed.

TCP speeds retrieved during a blank run along the circuit shown in Fig. 3.a projected on the ABS dome (Fig. 9) show that the sampling noise of the bare mesh affects the TCP trajectory and robot motion which appear discontinuous with abrupt changes in the printhead orientation both along straight lines and right angles. Robot joints are subjected to continuous acceleration/deceleration and, as shown in Fig. 9, the ensuing TCP speed fluctuates between 4 and 14 m/s, well below the targeted 15 mm/s. Filleting right angles with 1 mm radius does not improve TCP speed and oscillations, however a radical improvement in robot motion was obtained applying the optimized PB and NA smoothing conditions. Fig. 9 shows that with a minimum fillet radius of 1 mm the average TCP speed was close to 15 mm/s with minimal fluctuations between 15. the average TCP speed matches with the nominal speed up to 30 mm/s. Above 30 mm/s the average speed starts to diverge with a pronounced increase of speed fluctuations from 10 to 70% when the nominal speed increases from 15 to 50 mm/s, respectively. According to equation [START_REF] Lewis | Device fabrication: Three-dimensional printed electronics[END_REF] and the proportional correlation between the printing speed and line resistance, a maximum speed fluctuation of ca. 10% was used as constraint to limit line resistance variations in the printed circuit. In the present study, nominal speeds up to 30 mm/s can be used to print linear segments, whereas, in the presence of 1 mm filleted right angles or U turns, the nominal speed should not exceed 15 mm/s.

Geometry and shape fidelity of printed circuits

Fig. 11 shows that circuits printed using the optimized protocol (ie. mesh smoothing with10 PB loops, 2D circuits projection on the 3D mesh, discretization with 0.2 mm segments and normals alignment with 30 NA loops) displayed excellent shape fidelity with no detectable fluctuations in the geometry of printed patterns and a maximum deviation of 2.5% in the length of 3D printed circuits with respect to pristine 2D circuits. 3D scans of the circuit just after ink deposition (Fig.

11.b and 12) showed that printed lines had a 0.265 x 0.20 mm height semi-elliptical section with relatively small variations in width and height. A 10-15% increase in line width was measured in filleted right angles and reflected the local speed decrease from ca. 15 to 13.5 mm/s, as recorded during blank runs. As summarized in Table 1, ink drying led to a decrease of the line height from 0.20 to 0.11 mm. Line width did not change and the resistance of printed lines on a control length of 35 mm was consistent with values calculated with equations ( 11) and ( 12), i.e. 0.8 and 0.76 Ohm, respectively. Line width and pitch determined by optical microscopy were ca.10% smaller than those determined using the laser profilometer (i.e 0.46 and 2 mm against 0.53 and 2.1 mm), however this discrepancy was associated with the profilometer sampling period of 0.05 mm and in line with the measurement accuracy and the slight geometrical deformation on ca. 2.5% generated by the planar circuit projection on the half sphere. Complex patterns (Fig. Despite the good match between the 2D circuit and the printed one, printed lines displayed periodical width fluctuations when compared to lines printed on a smooth PET films under similar conditions (Fig. 14.a). Those shape defects were ascribed to the typical paths generated during the FDM 3D printing process which, despite an intensive surface polishing, were not removed from the ABS dome thus inducing a deformation of the ink/substrate contact line and the formation of kinks (Fig. 14.b). Despite the good control of the line pitch during the printing process, line fluctuations can attain amplitudes of 0.1 mm thus inducing significant variations in the interline distance which, in the tested conditions, ranged from 0.50 to 0.35 mm indicating that a decrease of line pitch to ca. 0.70 mm can lead to lines contact.

Despite the high precision in printhead positioning obtained with the process developed in this study, the surface morphology of the 3D printed substrate represents a major limitation to attain a high print quality and line density. Substrate surface treatment by polishing, sand blasting or by the application of a primer varnish/solvent smoothing will be therefore explored in future work to improve the print quality. All SMD components were tested by placing tester electrodes on printed lines at ca. 5 mm distance from component connections. Nominal resistance and capacity values were obtained for 200k-20 Ohms resistors and capacitors, green and red leds worked correctly and neighbouring pins of the USB controller did not display short circuit.

3D printed hybrid circuit board

Conclusions

Results presented in this work demonstrate that 6 axis industrial robots coupled to high accuracy 3D scan and adequate smoothing algorithms, ie:

-Paul Bourke smoothing of the 3D mesh with an iteration number limited to 10

-Normals alignment of points lying on the trajectory with an iteration number limited to 30 can be effectively used to move a piezo jetting printhead along smooth and high precision trajectories with a nearly constant speed and to print conductive paths on 3D objects.

Moreover, the direct projection of 2D circuits on the planarized mesh and the subsequent 3D development led to a deformation of the pristine 2D circuit (ca. 2% in the tested conditions) and allowed the quick generation of 3D trajectories using standard pre-existing circuit designs.

Trajectories generated on a 3D printed ABS dome and optimized printing conditions (i.e. speed and circuit design) were used to print silver ink circuits which displayed geometry and resistance in line with the semi empirical model developed to predict line width and resistance and the geometry of the pristine 2D circuit. Only the width of printed tracks displayed not negligible wave-like fluctuations which were associated to patterns generated during FDM 3D

printing of the ABS dome, thus highlighting the role of substrate surface morphology in affecting the print accuracy.

The optimized circuit design and trajectory processing conditions were successfully used to print a 3D circuit board mimicking a NFC card reader and including a USB PCB, SMD components and a NFC antenna thus showing the viability of the developed process to manufacture low-density 3D PCB.

Overall, this work demonstrates that 6 axis robots coupled to a sequence of unit operations integrated in a control software, i.e.:

i) 3D scan-mesh smoothing,

ii) 2D circuit projection on the planarized mesh and 3D trajectory generation, iii) normals alignment of points in the 3D trajectory, iv) retrieve of the effective toolhead TCP translation speed along the 3D trajectory with a blank run, v) synchronization of the ink jetting frequency with the effective average speed, vi) generation of the print code and print, allows to obtain a duty cycle, from the 2D circuit to the 3D printed one, of less than 30 min.

Fig. 1 .

 1 Fig. 1. Detail of the robotic cell developed in this study and 3D printed ABS half sphere used as model substrate a). Screenshot of the offline programming environment b).

Fig. 2 .

 2 Fig. 2. Influence of the coverage ratio and of the number of superimposed unit drops, 𝑛 = , on the diameter of the printed dot. Data obtained on a smooth PET film were added to compare ABS with a substrate with good ink spreading.

Fig. 3 .

 3 Fig. 3. Planar circuits design used in this study. a) Large interdigitated path, b) antenna-type three-windings coil, c) patterns with 1 mm pitch.

Fig. 4 .

 4 Fig. 4. a) 2D scheme of the effect of planar scan on the cell size and normal angle deviation between neighbouring points/cells. b) Cell size (arc length) and normal angle variations calculated with eqs. (6, 7) and plotted as a function of the distance from the centre of the half-sphere.

Fig. 5 .

 5 Fig. 5. a) design of the NFC reader hybrid PCB including SMD components and b) 3D structural element of the card reader used as printing substrate.

Fig. 6 .Fig. 7 .

 67 Fig. 6. a) Average normal angle variation and standard deviation along printing trajectory plotted as a function of the number of loops of the PB smoothing algorithm. Inset shows raw angle variation data. b) Images of the meshed surface and normals before and after 10 loops smoothing.

Fig. 8 .Fig. 8 .

 88 Fig. 8.b shows that mesh processing with 10 loops PB smoothing, further increases the convergence of the NA algorithm which attains stable angle variations and standard deviation

  7 and 14 mm/s. The progressive increase of the nominal TCP speed slightly affected the translation speed up to 30 mm/s (Fig. 10.a), when abrupt speed variations start to appear in 4 mm filleted right angles indicating that robot joints are no longer able to provide the angular acceleration needed to prevent TCP deceleration. Fig. 10.a shows that the further increase of the nominal to 50 mm/s led to large fluctuations in the TCP speed. The printhead attains the nominal speed in short segments along the circuit trajectory and most of the time is in a transient acceleration/deceleration regime since in 1 mm filleted right angles the TCP speed drops to ca. 23 mm/s.

Fig. 9 .

 9 Fig. 9. Effect of normal smoothing and angle fillet on the speed variation of the printhead TCP along the large interdigitated circuit (Fig. 3.a) projected on a 3D mesh smoothed with 10 PB loops. The robot nominal printing speed is 15 mm/s.

Fig. 10 .

 10 Fig. 10.b shows that, despite a slight increase in speed scattering (i.e. the standard deviation),

Fig. 10 .

 10 Fig. 10. a) Effect of the nominal printing speed on the effective speed along a trajectory with 4 mm filleted angles and processed with 10 PB and 30 NA loops. b) Plot of the average TCP speed along the circuit trajectory and of the maximum speed decrease as a function of the nominal printing speed.

Fig. 11 .Fig. 12 .

 1112 Fig. 11. a) Image of the antenna-type circuit printed on the ABS semi sphere (nominal speed 15 mm/s, frequency 88 Hz). b) Screenshot of the high definition scan used to characterize the circuit geometry after printing (scan sampling period 0.05 mm).

Fig. 13 .

 13 Fig. 13. Patterns with 1 mm pitch shown in Fig. 3.c printed on the ABS dome a), b). Parallel lines with 1 mm pitch and, highlighted in red, forms with bending radius of 0.5 mm c). axial lines with 1 mm gap d) and U turn with 1 mm radius e). Speed fluctuations in forms with bending radius of 0.5 mm, f).

Fig.Fig. 14

 14 Fig. 15a-e summarized the protocol used to project the optimized circuit (Fig. 15.b) on the substrate and to generate 3D trajectories which had a maximum length deviation from the pristine 2D circuit of 1.4%. The 4-turn antenna coil had a cumulative length (including segments to connect with the USB controller) and resistance of 953 mm and 21.1 Ohms, respectively, which were in line with: i) the resistance of 35 mm control line expected for the printing conditions (i.e. 0.77 Ohm/35 mm vs. the 0.8-0.85 Ohms/35 mm given in Tables1 and 2), ii) the low coil resistance required in NFC antennas to avoid signal attenuation.

Fig. 15

 15 Fig. 15 Design and printing on a 3D substrate of a circuit for badge identification via communication. a) CAD model of the 3D substrate. b) CAD model of the 2D circuit adapted for 6-axis robot printing. c) Scan 3D of the 3D-printed ABS substrate and planarization of the printing region. d) Orthogonal projection of the 2D circuit on the planarized surface. e) Development of the 3D shape and generation of the 3D circuit. f) Final 3D circuit printed on the ABS substrate. g) Detail of the USB controller region before and after components placement. h) Detail of the NFC antenna region before and after components placement.

Table 1 .

 1 Characteristics of circuits shown in Fig. 3a and 3b printed using a nominal speed of 15 mm/s and a jetting frequency of 88 Hz. *The index w refers to the wet ink and d to the dry ink.

			Straight line			Curve		Coil
		Width	Height*	Resistance	Width	Height*	Resistance	Pitch
		(mm)	mm)	(Ohm)	(mm)	(mm)	(Ohm)	(mm)
	Theoretical	0.44	0.20w/0.11d	≤ 0.76	0.46	0.22w/0.12d	≤ 0.64	2
	Measured	0.46±0.05	//	0.80±0.05	0.51±0.03	//	//	2±0.03
	Measured	0.53±0.05 0.20w±0.02	//	0.61±0.05 0.21w±0.02	//	2.1±0.05
	from 3D scan		0.12d±0.03			0.13d±0.02

Table 2 .

 2 Characteristics of the circuit shown in Fig. 3c printed using a nominal speed of 15 mm/s and a jetting frequency of 75 Hz. *The index w refers to the wet ink and d to the dry ink.

			Straight line		Curve (1 mm radius)	Line
			Height*	Resistance	Width	Height*	Resistance	Pitch
		(mm)	mm)	(Ohm)	(mm)	(mm)	(Ohm)	(mm)
	Theoretical	0.42	0.18w/0.10d	< 0.90	0.44	0.2w/0.11d	< 0.76	1
	Measured	0.50±0.05	//	0.85±0.05	0.54±0.03	//	//	1±0.02
	Measured	0.54±0.05 0.17w±0.03		0.58±0.05 0.21w±0.03
	from 3D scan		0.08d±0.02			0.10d±0.02

Acknowledgments

This work was supported by the Region Auvergne Rhone Alpes (project R&D Booster 2020 IOT 3D FLEXFAB, grant n° 20 010403 01 -22437).