# Using loop transformations for precision tuning in iterative programs

Youssef Fakhreddine, Guillaume Revy

► To cite this version:

Youssef Fakhreddine, Guillaume Revy. Using loop transformations for precision tuning in iterative programs. ARITH 2023 - IEEE 30th Symposium on Computer Arithmetic, Sep 2023, Portland, United States. pp.159-166, 10.1109/ARITH58626.2023.00031 . lirmm-04939356

HAL Id: lirmm-04939356

https://hal-lirmm.ccsd.cnrs.fr/lirmm-04939356v1

Submitted on 10 Feb 2025

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Using loop transformations for precision tuning in iterative programs

Youssef Fakhreddine and Guillaume Revy

Univ Perpignan Via Domitia, DALI, Perpignan, France

LIRMM, Univ Montpellier, CNRS (UMR 5506), Montpellier, France

***Abstract***—Many floating-point formats are available, providing all different levels of precision. By mixing several of these formats in the same program, it is possible to achieve good performance while maintaining an acceptable level of output accuracy. Therefore various tools have been designed to adapt the precision of computations in floating-point programs for performance and accuracy purposes. However most of them do not consider the iterative nature of these programs. This article presents a tool that enables to adapt the precision of floating-point computations in iterative routines, at the iteration level. This tool is based on multiple-precision computations to evaluate the impact of some format adaptations on the output accuracy, and it uses the delta-debugging to isolate the most relevant instruction set to be tuned. The originality of our approach is that it relies on static loop transformations to duplicate loop body instructions, and thus to increase the number of possible instructions that can be targeted. These transformations include especially the loop splitting and unrolling, which enable to allocate different precisions for different iterations, and thus to improve the tuning process. We show the advantages of this approach on a representative set of iterative programs.

**Keywords:** floating-point arithmetic, mixed precision, dynamic precision tuning tool, multiple-precision, loop splitting and unrolling, delta-debugging algorithm

# I. INTRODUCTION

The representation of real numbers in computer programs is essential for scientists in many fields such as numerical analysis, physics, and engineering. For example, real numbers are used to represent a wide range of physical quantities, such as temperature, distance, and velocity. However, representing such numbers in a computer program can be challenging due to the limited precision of digital computers.

Nowadays, floating-point arithmetic is widely used to represent real numbers on computer systems. For this purpose, the IEEE-754 2019 standard provides four floating-point formats, i.e. binary16, binary32, binary64, and binary128 [1]. Other non IEEE formats are also available, like bfloat16 developed by Google, TensorFloat-32 by NVidia, or again Posit proposed by Gustasfon. All these offer different levels of accuracy. Given these formats, the difficulty now for the developers is to choose the appropriate format for each floating-point element of a program to satisfy an accuracy constraint. Consequently, from an accuracy point of view, they use most of the time the binary64 floating-point format, which is the more accurate format directly available on most modern hardware architectures. Yet, from a performance point of view, the length of the type chosen for floating-point computations can have a great

impact on the execution time of the program. For example, loading binary32 data from memory is cheaper than loading binary64 data, and computing with binary32 data is cheaper, as well. This is particularly true on SIMD units since it is possible to pack twice as many binary32 than binary64 data for the same register size.

Several works have already shown that sometimes lowering the precision level of some variables, functions, or instructions in numerical programs does not cause too much loss of final accuracy, but it allows the program to gain speed. For example, Baboulin et al. showed in [2] that such mixed-precision programs can achieve similar accuracy to the original binary64 precision program while being significantly faster and reducing memory pressure. However, this work was done by hand, and guided by the complexity of the subroutines in use. And, in the general case, choosing the appropriate format for each data or each instruction in a program is a difficult task, which justifies the growing interest in auto-tuning tools.

In the last decade, many research projects are thus interested in auto-tuning for floating-point precision in numerical programs. In this context, two approaches coexist. On the first side, static approaches: such tools analyze the program without executing it to determine a precision requirement of each variable in order to satisfy an accuracy constraint, and therefore it works at the source code level. Among these tools, we find for example Daisy [3], [4]. It combines affine and interval arithmetics together with rewriting optimizations and SMT-solving to produce finite precision programs that fulfill accuracy requirements. In addition, it outputs correctness certificates that can be checked by Coq or HOL4. FPTuner [5] produces also mixed-precision programs, and it determines rigorous errors bounds by using Symbolic Taylor Expansions [6]. Another alternative is proposed by Martel in [7]. This tool relies on a 2-step static analysis (forward and backward) done by abstract interpretation to determine the precision of each variable. It outputs mixed-precision programs, as well, and constraints that can be afterward checked by an SMT solver. These techniques have also been used in the SMT-based version of POP [5]. These approaches give good results but remain limited for large programs and complex structures like loops. On the other side, dynamic approaches: Most of them rely on a trial-and-error strategy, where a step is to compare the results of both the original program and the transformed one to estimate the accuracy loss. This is for example implemented in Craft [8], [9] or Precimonious [10]. The former uses a
```c
int main(void) {
    double S = 0;
    #pragmaclang loop split_optimization(enable)
    for (int i = 1; i <= 1000; i++)
        S = S + 0.01;
    printf("S = %.20e", S);
    // ... |S - Optim|/|S| < 1e-6 ?
    check_relnum_error(S, 1e-6);
    return 0;
}
```

<caption>Listing 1: C program computing the sum S in (1).</caption>

breadth-first search strategy to explore various mixed-precision configurations, while the latter relies on the delta-debugging algorithm [11] to isolate a configuration that fulfills the accuracy requirement and produces a mixed-precision program faster than the original one. One of the limitations of these approaches is the combinatorics linked to the number of possible transformations. To tackle this issue, two extensions of Precimoniou s have been proposed. Blame [12] reduces its search space by quickly exhibiting variables whose precision does not impact the final result, and HiFPTuner [13] groups dependent data requiring the same level of precision in communities and it applies transformations directly to these communities. Other alternatives also exist: for example, Promise [14] and Verrou [15] combine both Discrete Stochastic Arithmetic through the CADNA library [16] and delta-debugging to produce mixed-precision programs, while ADAPT [17] relies on automatic differentiation to determine precision requirements of inputs and intermediate variables. Recently, Pherbie [18] has been developed on top of Herbie [19], which uses rewriting rules and multiple-precision together with local error analyses to adapt the precision of sub-expressions and to produce a set of mixed-precision candidate implementations optimized for accuracy and speed. These tools allow to tackle larger problem than static ones, but they do not consider the iterative nature of the program.

Generally, less attention has indeed been paid to the tuning of iterative programs. This is one of the limitations of Precimoniou s mentioned by its authors, and only a few attempts have been proposed in [20], [21], [22]. Hence, in this article, we focus on dynamic auto-tuning tools targeting iterative routines. The originality of our approach is to use static loop transformations to increase the number of instructions appearing in the program, improving consequently the tuning process by allowing the allocation of different precisions at different iterations.

The main contributions of this article are the following:

* A dynamic auto-tuning tool targeting iterative routines, and that relies on three modules: static loop transformations to increase the combinatorics of possible transformations, the delta-debugging to find the most relevant ones, and the fp2mp module to evaluate the impact on numerical accuracy of tested configurations,

* And a set of examples, to show the interest of our approach, and particularly to show how it allows us to improve the result of the tuning process to the detriment of an increase in the tuning time.

<figure><img src="image_2.png" /><figcaption>Figure 1: Relative error for the computation of S in (1), according to the number of binary32 iterations.</figcaption></figure>

The article is organized as follows: It starts by a motivating example in Section II to show the interest of our approach. Then Section III presents our tool, and in particular the different modules it relies on. Section IV gives some experimental results before a conclusion in Section V.

## II. MOTIVATING EXAMPLE

To demonstrate the usefulness of our approach, let us consider the program in Listing 1,<fnref n="1" /> that computes the sum

$$S = \sum_{i=1}^{1000} 0.01 = 10. \qquad (1)$$

and where all the computations are performed in binary64 arithmetic. It returns $S = 9.99999999999983124610^2$. If the program was implemented using binary32 arithmetic, the result would be $S = 10.00013351440429687500^3$, meaning with a relative error of about $10^{-5}$ compared to the result in binary64.

Now, when dealing with auto-tuning tools, the question is the following: What data or what instructions can be used in lower precision (i.e. binary32) to achieve a result as accurate as the one obtained using only the higher precision (i.e. binary64) for a given threshold? For example, let us consider relative error and a threshold of $10^{-6}$. The program contains two data to be tuned (the variable S and the constant 0.01), and existing approaches will give no solution because they will detect only three possible transformations, which all give a result with a relative error higher than $10^{-6}$ compared to binary64 result. In our case, we aim at adapting the precision of instructions: it is worth it since only one possible transformation is detected resulting in the binary32 result.

However, some of these 1000 iterations, whose body contains only one instruction each, could be performed in binary32 arithmetic, while keeping an acceptable output accuracy. For example, let us see what happens if the program would perform the first iterations in binary32, and the last ones in binary64, only. Figure 1 shows the relative error of such a result compared to the full binary64 result, according to the number of binary32 iterations. We can observe that at most 458 iterations could be performed in binary32 before the relative error exceeds the threshold of $10^{-6}$.

<footnote n="1">¹The program is here decorated to be used by our tool.</footnote>

<footnote n="2">²Exact value in binary64 is 5629499534213025 · 2⁻⁴⁹.</footnote>

<footnote n="3">³Exact value in binary32 is 20971800 · 2⁻²¹.</footnote>
<figure><img src="image_3.png" /><figcaption>Figure 2: Workflow of our auto-tuning tool.</figcaption></figure>

Now to push this idea to the limit, we can fully unroll this loop, thus resulting in a program of 1000 distinct instructions. Then by using the delta-debugging algorithm presented in Section III-D, we can isolate 526 instructions out of these 1000 that can be transformed while maintaining the relative error below the threshold of $10^{-6}$.

Both of these results offer a huge improvement, but they remain untrackable by existing tools. Therefore this justifies the interest of auto-tuning tools that consider the iterative nature of programs, like the one presented in this article.

### III. AUTO-TUNING TOOL FOR ITERATIVE PROGRAMS

Our goal is to develop a dynamic auto-tuning tool that considers the iterative nature of scientific computing routines. The originality of our approach is that it relies on compilation techniques, including loop splitting and unrolling, to increase the number of instructions to be targeted in the program, without modifying its semantics, and consequently to increase the number of possible configurations. Then these configurations are exposed to the delta-debugging algorithm, to isolate a more relevant one, meaning one that allows to keep an acceptable output accuracy. And each time a configuration is considered, its impact on output accuracy is evaluated using multiple-precision through the `fp2mp` tool.

This section starts with an overview of its general workflow. Then it explains the different components of this flow, before going back to motivating example.

#### A. General workflow

Our tool is a command-line tool written in Python, and that requires a program together with a constraint. In this program, the loop to be tuned is pointed out by a new pragma we have introduced in the LLVM front-end, as shown in Listing 1. Moreover, the command-line gives also to the user the capability of selecting the formats including standard floating-points formats (i.e. binary16, binary32, and binary64) and bfloat16, as well as the strategy (i.e. splitting or unrolling) used in the tuning process and detailed in Section III-B. Then it works iteratively, by reducing the format of each instruction to the format immediately below it in the list of selected formats. It then repeats this step, as long as format adaptations are possible. Thus, it allows to mix several formats, and not only two, in the same program.

The general flow of our tool is given in Figure 2. As most dynamic tool, ours uses a trial-and-error strategy. Precisely, it proceeds in 3 main steps:

1. First, it starts by applying loop transformations on the input. These transformations are guided by hints given by the user through command-line parameters, and they currently include loop splitting and unrolling. This step outputs a modified LLVM intermediate representation.

2. Second, for each floating-point instruction of the programs, it inserts in this LLVM IR the equivalent computation in multiple-precision. This is done using the `fp2mp` tool, as explained in Section III-E. This tool outputs also a list of MPFR objects, whose precision can be further modified. This is used to build the configuration sets as explained in Section III-C.

3. Third, using the delta-debugging algorithm, we determine a maximal subset of transformations that allows to still obtain an acceptable output. This is explained in Section III-D.

At the end of this process, we do not generate a new optimized program, but we give some hints to the user on how to produce it through a list of acceptable transformations.

#### B. Static loop transformations

The originality of our approach comes from the static loop transformations applied on the original program. Their objective is to increase the number of possible transformations. For doing this, our tool leverages the LLVM capabilities of transforming programs.

LLVM is a set of tools and libraries that are used for compiler construction. One of the key features of LLVM is its intermediate representation, which serves as a high-level assembly language that is portable across different platforms and architectures. The LLVM IR is designed to be easily optimized by various transformations during multiple passes. These passes allow developers to perform various optimizations on the program, such as unrolling loops, inlining functions, and vectorizing programs. These optimizations can help to improve the performance of a program, particularly for programs that have a high degree of iterative or repetitive behavior.

In this project, we use the capability of unrolling loops available in LLVM, and the one of splitting loops developed by Revy [23]. These are both loop transformations that make instructions appeared in the program, but without modifying its semantics. Here we called *strategy* a combination of a loop transformation and a factor. For example, on the program of Listing 1, the strategy of unrolling loop by a factor of 2 gives the following program,

```c
for (int i = 1; i <= 500; i++) {
    S = S + 0.01;
    S = S + 0.01;
}
```

while the one of splitting the loop by a factor of 2 gives the following one.

```c
for (int i = 1; i <= 500; i++)
    S = S + 0.01;
for (int i = 501; i <= 1000; i++)
    S = S + 0.01;
```
This strategy is given to LLVM core by inserting appropriate metadata in the LLVM IR.

These transformations allow for determining two different patterns of transformations. Indeed the unrolling allows the detection of some transformations appearing regularly throughout the iterations. On the other hand, the splitting allows detecting transformations that are relevant on a subset of iterations, only (e.g. at the beginning, at the end, ...).

Note that this approach is a bit antagonist to the existing ones. Indeed the current trend is to reduce the combinatorics of targeted instructions to speedup the process: the fewer the possible instructions, the less configurations to test, and the faster the tuning process. This is actually done in Blame or HiFPTuner that reduce the combinatorics compared to Precimonious, for example. Conversely, in our approach, we increase this combinatorics, but in a reasonable way. Doing this way, we certainly increase the time necessary for the process to converge, but we significantly improve the quality of this process, as illustrated in Section IV.

### C. Building transformation sets

This section deals with the way used to build the set of possible transformations. This set is given to the delta-debugging to determine the final configuration.

Once the loop transformation is performed, and the `fp2mp` is run on the modified LLVM IR, we get a list of MPFR objects for each function. This list corresponds to all the floating-point elements of the input program (i.e. variables, constants, and operations). It is given in a JSON file, whose entries are like the one below.

```json
{
  "name": "add.mpfr",
  "fp_name": "add",
  "fp_type": "double",
  "precision": "53",
  "tuned": "true",
  "location": {"line": "12", "column": "11"},
  "op_code": "fadd",
  "block": "for.body",
  "loops": [
    {"name": "for.cond", "depth": "1"}
  ]
}
```

Hence for each of them, it contains especially its floating-point type and its precision, its location, and possibly the loop it belongs to in the program.

Using this information, it is possible to build a list of possible transformations. In our case, we target floating-point instructions. Hence a transformation is a couple composed of an instruction and a list of formats, all lower than the format of the instruction itself. They correspond to the formats to which the tool will try to reduce the instructions. For example, if the user ask for using `binary16`, `binary32`, and `binary64`, in this case, we will have the following possible transformations for the addition `add` below.

`add -> [b32, b16]`

This means that the tool will start by trying to reduce the format of the `add` instruction to the `binary32` format. If it succeeds, it will then try to reduce it to the `binary16` format.

<figure><img src="image_4.png" /><figcaption>Figure 3: Delta-Debugging workflow.</figcaption></figure>

For performance purposes, when the splitting strategy is used, it can make sense to use the same format for all the instructions in a loop. In this sense, we implemented another way to define transformation sets, by grouping all the instructions of a loop body in the same set and by assigning a unique set of possible target formats to this group. This is used for example in Sections IV-B and IV-C.

### D. The delta-debugging algorithm

The delta-debugging algorithm, developed by Zeller and Hildebrandt, was initially designed for the simplification of test cases. This is a process of reducing the complexity of a test case while still maintaining the ability to expose a bug. This is useful in debugging because it can help to isolate the cause of a bug, making it easier to fix. Its workflow is given in Figure 3. It is implemented in a tool called `ddmin` which is widely used and it has been proven to be effective in reducing the size of test cases for many types of bugs. It works like a binary search to determine a 1-minimal set $\Delta^+$ of changes that reproduces the bug. This means that removing only one change from $\Delta^+$ makes the bug go away.

Several works have been published to adapt it to specific use cases. In the context of auto-tuning, it has been adapted to be used by other tools, like for example Precimonious or Promise. In our case, we use `ddmax`, a variant of the delta-debugging `ddmin` routine [24]. Given a set $\Delta$ of possible transformations, its goal is to maximize the subset $\Delta^+$ of the input set $\Delta$ that makes the initial program still computing an acceptable output. At the end of the algorithm, $\Delta^+$ is 1-maximal, which means that adding one transformation in $\Delta^+$ will make the accuracy test failed. Actually `ddmax` uses the same techniques as `ddmin`: maximizing $\Delta^+$ can be done by systematically minimizing $\Delta - \Delta^+$, starting with large differences and then smaller and smaller ones, until every remaining differences would cause $\Delta^+$ to produce unacceptable output.

From an implementation point of view, our tool relies on the DD Python module.<fnref n="4" /> This provides a class that implements the delta-debugging algorithm. In order to adapt it to our needs, we have overloaded the `test()` method. It is in charge of testing whether a configuration of transformations allows to satisfy the input constraint or not. For doing this, it modifies the precision of the corresponding MPFR objects in the LLVM IR using the `update` module of `fp2mp`.

<footnote n="4">⁴See https://github.com/grimm-co/delta-debugging/.</footnote>
## E. The fp2mp verification tool

Dynamic auto-tuning tools often rely on a trial-and-error strategy, for a dataset for which a reference result is known. Particularly, they iteratively apply transformations on the input until the accuracy test fails. For doing this, this test procedure uses a module to evaluate the impact on the output accuracy of some transformations on the input program.

In our case, we use multiple-precision computations. Particularly, our auto-tuning process relies on the `fp2mp` tool. This tool is designed for experimentation with floating-point precision in computer programs. It works by duplicating in the LLVM IR each floating-point instruction used in a program with its equivalent in multiple-precision through the MPFR library. This allows the user to define the precision they desire, and to evaluate the impact of these precision modifications on the behavior of the program.

The interest of `fp2mp` in our approach is that this module offers the capability of adapting the precision of certain computations only. Even if it works on `LLVM IR`, for a sake a clarity, let us see how it would work at a C level. For example, after having unrolled the loop by a factor of 2, Listing 1 would be transformed as follows where `Smpfr` holds the equivalent value of `S` but in `MPFR`.

```c
double S = 0;
mpfr_t SMPfr;
mpfr_init2(Smpfr, 53);
mpfr_set_d(Smpfr, 0., MPFR_RNDN);
for(int i = 0; i <= 1000; i++) {
    S = S + 0.01;
    mpfr_add_d(Smpfr, SMPfr, 0.01, MPFR_RNDN);
    S = S + 0.01;
    mpfr_add_d(Smpfr, SMPfr, 0.01, MPFR_RNDN);
}
mpfr_clears(Smpfr, (mpfr_ptr) 0);
```

Hence adapting the precision of only the second addition to the binary32 format, for example, would result in the following piece of C code, where `C1`, `C2` and `C3` are three temporary variables representing values in binary32, that is, in precision 24, to ensure the computation is done in precision 24, as well.

```c
double S = 0;
mpfr_t SMPfr, C1, C2, C3;
mpfr_init2(Smpfr, 53);
mpfr_inits2(24, C1, C2, C3, (mpfr_ptr) 0);
mpfr_set_d(Smpfr, 0., MPFR_RNDN);
for(int i = 0; i <= 1000; i++) {
    S = S + 0.01;
    mpfr_add_d(Smpfr, SMPfr, 0.01, MPFR_RNDN);
    S = S + 0.01;
    mpfr_set(C1, SMPfr, MPFR_RNDN);
    mpfr_set(C2, 0.01f, MPFR_RNDN);
    mpfr_add(C3, C1, C2, MPFR_RNDN);
    mpfr_set(Smpfr, C3, MPFR_RNDN);
}
mpfr_clears(Smpfr, C1, C2, C3, (mpfr_ptr) 0);
```

## F. Back to the motivating example

To reinforce the interest of our tool, let us now go back to the motivating example of Section II. Using a strategy of splitting the loop by a factor of 20 would result in the following piece of program.

```c
int main(void) {
    double S = 0;
```

```c
    for (int i = 1; i <= 50; i++)
        S = S + 0.01;
    for (int i = 51; i <= 100; i++)
        S = S + 0.01;
    // ...
    for (int i = 951; i <= 1000; i++)
        S = S + 0.01;
    // ... |S - Optim|/|S| < 1e-6 ?
    checkreverse_rel_error(S, 1e-6);
    return 0;
}
```

Hence we no longer have a single instruction, but 20 separate ones. Assuming two formats (i.e. `binary64` and `binary32`), the combinatorics amounts to $2^{20} = 1048576$, meaning 1048576 possible configurations.

By running our tool, in this case, we found in about 10 sec. 9 instructions out of these 20 ones that could be transformed so that the result of this transformed program remains at the required threshold of $10^{-6}$. This means 45% of the instructions are performed in `binary32` arithmetic, which is actually a significative improvement compared to existing methods. And this result is not far from the one obtained by fully unrolling the program in Section II. However, if the unroll strategy by a factor of 20 would have been used, only 1 instruction out of the 20 would have been isolated to be relevant, that is, only 5% of the instructions.

Therefore, the quality of the result directly depends on the strategy used. However, finding the best strategy for a given problem is a difficult task that we do not discuss in this article, we leave this choice to the user.

# IV. EXPERIMENTAL RESULTS

In this section, we present some experimental results to show the interest of our analysis techniques.

## A. Impact of different strategies

In this first experiment, we use 5 iterative programs: the `sum` problem from Listing 1 (1000 iterations), the `riemann` and `simpson` methods that compute the integral of a function (400 and 1000 iterations, respectively), the `arclength` algorithm inspired from [10] and that computes the arc length of a function (1000 iterations), and finally the `nbody` problem inspired from [25] that predicts the individual motions of a group of celestial objects interacting with each other gravitationally (100 iterations). All these programs are implemented in `binary64` floating-point arithmetic, and we consider reducing the precision of instructions to the binary32. For doing this, for each problem, we ran our auto-tuning tool with different strategies and different thresholds. For all the problems, the threshold considered is the relative error of the result computed by the transformed program compared to the one in `binary64` (which is considered as the "exact" result), but for `nbody` where the threshold is the absolute distance between the 3D coordinates computed by the transformed program and the ones by the original one. Tables I and II show the impact of using our approach. In these tables, "delta-debugging" means the delta-debugging is applied without loop transformations, like other tools would have done. Moreover, for each problem,
<figure><img src="image_5.png" /><figcaption>Figure 4: Iteration patterns for the *riemann* method and different splitting factors.</figcaption></figure>

<figure><img src="image_6.png" /></figure>

these give the number of binary64 and binary32 instructions, together with the transformation percentage.

We observe that in all cases, applying static loop transformations allow to improve the quality of the tuning result. For example, let us consider the *arclength* problem and the *splitting* strategy. For a threshold of $10^{-8}$, the delta-debugging alone found 5 out of the 11 instructions of the main loop that can be transformed, i.e. 45.5%. Splitting the loop into 50 subloops increases the combinatorics: now 550 instructions are available in the program, i.e. $2^{50}$ possible configurations. Applying the delta-debugging on it allows to achieve 86.9% of transformations, i.e. 478 instructions among these 550 ones, which is a huge improvement. And the improvement is even better for sharper thresholds. Still on *arclength* problem and *splitting* strategy by a factor of 50, for a threshold of $10^{-11}$, no solution is found by applying only the delta-debugging, while our solution found 195 instructions out of the 550 to be tuned, i.e. 35.5%. The same holds for other problems, and the improvement can be even better, like for *riemann* method and a threshold of $10^{-8}$, where the *unrolling* strategy is used with a factor of 50. In this case, almost all the instructions are transformed (i.e. 99.2%), while the delta-debugging alone isolates only 40% of the instructions.

However, whatever the strategy, the higher the factor, the greater the number of instructions in the transformed program, and therefore the more important the combinatorics. In this case, the tuning time also increases. For example, still on the example of *arclength* problem and *splitting* strategy, the tuning time goes from a few seconds to several tens of minutes for a factor of 50 (i.e. 22 min. for a threshold of $10^{-8}$ to 41 min. for $10^{-11}$). But this can be contained: for example, for a factor of 20, the tuning time does not exceed 10 min.

## B. Same transformations for the whole loop body

From a performance point of view, it is better to have instructions of the same format in a given loop body. For example, it helps in activating some compiler optimizations like vectorization [13]. Hence, in this experiment, we ran our auto-tuning tool by forcing all the instructions of a same loop to be in the same format. (This does not make sense for *unrolling* strategy, but only for *splitting* strategy.) We did

this experiment for the *riemann* problem, for a splitting factor ranging from 1 to 30. Figure 4 shows the results for a threshold of $10^{-8}$ and $10^{-9}$, where the iterations in binary32 are in blue, and the ones in binary64 are in red.

As expected, we observe that the value of the splitting factor impacts the number of instructions transformed. For example, in Figure 4a, for a threshold of $10^{-8}$, and for a factor in $\{1, 2\}$, no transformation is found and all the iterations are in binary64 arithmetic. For a factor of 3, the first third of iterations (i.e. 133 iterations) are transformed in binary32, while the two last thirds remain in binary64. And this ratio increases with the value of the splitting factor. The same holds for a threshold of $10^{-9}$ in Figure 4b, but since the threshold is sharper, the number of binary32 iterations is lower.

Note in this experiment that, as mentioned in Section III-B, the *splitting* strategy helps to find patterns where the transformations appear at the beginning of all the iterations.

## C. Auto-tuning of unbounded loops

In the previous experiments, the constraint was on the accuracy of the output. The constraint can also relate to the number of iterations, for example, which is particularly interesting in the case of unbounded loops. In this third experiment, we consider the *conjugate gradient* given in Algorithm 1. This is a well-known algorithm to solve the linear system $A \cdot x = b$ in the case of symmetric definite-positive matrices. It iterates as long as the residual $\|r_k\|$ is greater than a given threshold. And it is actually known to be extremely sensitive to the computation precision, which impacts directly this number of iterations.

**Algorithm 1** Conjugate Gradient method.

1: $r_0 := p_0 := b - Ax_0$, and $k = 0$

2: **while** $\|r_k\| \ge \epsilon$ and $k < \maxiter$ **do**

3: &nbsp;&nbsp;&nbsp;&nbsp;α_k := $\frac{r_k^T r_k}{p_k^T A p_k}$

4: &nbsp;&nbsp;&nbsp;&nbsp;x_{k+1} := x_k + α_k p_k

5: &nbsp;&nbsp;&nbsp;&nbsp;r_{k+1} := r_k - α_k A p_k

6: &nbsp;&nbsp;&nbsp;&nbsp;β_k := $\frac{r_{k+1}^T r_k}{r_k^T r_k}$

7: &nbsp;&nbsp;&nbsp;&nbsp;p_{k+1} := r_{k+1} + β_k p_k

8: &nbsp;&nbsp;&nbsp;&nbsp;k = k + 1

9: **end while**
<table><caption>Table I: Auto-tuning results for splitting strategy.</caption><thead><tr><th rowspan="2">Program</th><th rowspan="2">Threshold</th><th colspan="4">delta-debugging</th><th colspan="4">split / factor = 10</th><th colspan="4">split / factor = 20</th><th colspan="4">split / factor = 50</th></tr><tr><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th></tr></thead><tbody><tr><td rowspan="4">sum</td><td>1e-5</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>5</td><td>5</td><td>50.0</td><td>0:00:07</td><td>9</td><td>11</td><td>55.0</td><td>0:00:13</td><td>22</td><td>28</td><td>56.0</td><td>0:00:39</td></tr><tr><td>1e-6</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>6</td><td>4</td><td>40.0</td><td>0:00:06</td><td>11</td><td>9</td><td>45.0</td><td>0:00:13</td><td>27</td><td>23</td><td>46.0</td><td>0:00:43</td></tr><tr><td>1e-7</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>9</td><td>1</td><td>10.0</td><td>0:00:09</td><td>17</td><td>3</td><td>15.0</td><td>0:00:18</td><td>43</td><td>7</td><td>14.0</td><td>0:00:53</td></tr><tr><td>1e-8</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>10</td><td>0</td><td>0.0</td><td>0:00:08</td><td>20</td><td>0</td><td>0.0</td><td>0:00:18</td><td>48</td><td>2</td><td>4.0</td><td>0:01:01</td></tr><tr><td rowspan="4">riemann</td><td>1e-8</td><td>3</td><td>2</td><td>40.0</td><td>0:00:04</td><td>9</td><td>41</td><td>82.0</td><td>0:00:25</td><td>1</td><td>99</td><td>99.0</td><td>0:00:43</td><td>9</td><td>241</td><td>96.4</td><td>0:01:42</td></tr><tr><td>1e-9</td><td>3</td><td>2</td><td>40.0</td><td>0:00:04</td><td>21</td><td>29</td><td>58.0</td><td>0:00:48</td><td>32</td><td>68</td><td>68.0</td><td>0:01:27</td><td>38</td><td>212</td><td>84.8</td><td>0:03:23</td></tr><tr><td>1e-10</td><td>4</td><td>1</td><td>20.0</td><td>0:00:04</td><td>29</td><td>21</td><td>42.0</td><td>0:01:04</td><td>66</td><td>34</td><td>34.0</td><td>0:02:08</td><td>55</td><td>195</td><td>78.0</td><td>0:04:10</td></tr><tr><td>1e-11</td><td>5</td><td>0</td><td>0.0</td><td>0:00:04</td><td>44</td><td>6</td><td>12.0</td><td>0:00:55</td><td>73</td><td>27</td><td>27.0</td><td>0:03:56</td><td>59</td><td>191</td><td>76.4</td><td>0:04:27</td></tr><tr><td rowspan="4">arclength</td><td>1e-8</td><td>6</td><td>5</td><td>45.5</td><td>0:00:10</td><td>28</td><td>82</td><td>74.5</td><td>0:02:09</td><td>34</td><td>186</td><td>84.5</td><td>0:04:19</td><td>72</td><td>478</td><td>86.9</td><td>0:22:28</td></tr><tr><td>1e-9</td><td>8</td><td>3</td><td>27.3</td><td>0:00:10</td><td>59</td><td>51</td><td>46.4</td><td>0:03:36</td><td>108</td><td>112</td><td>50.9</td><td>0:06:24</td><td>121</td><td>429</td><td>78.0</td><td>0:38:03</td></tr><tr><td>1e-10</td><td>9</td><td>2</td><td>18.2</td><td>0:00:10</td><td>89</td><td>21</td><td>19.1</td><td>0:02:45</td><td>161</td><td>59</td><td>26.8</td><td>0:09:42</td><td>284</td><td>266</td><td>48.4</td><td>0:42:00</td></tr><tr><td>1e-11</td><td>11</td><td>0</td><td>0.0</td><td>0:00:08</td><td>106</td><td>4</td><td>3.6</td><td>0:02:22</td><td>203</td><td>17</td><td>7.7</td><td>0:08:52</td><td>355</td><td>195</td><td>35.5</td><td>0:41:16</td></tr><tr><td rowspan="4">simpson</td><td>1e-7</td><td>4</td><td>6</td><td>60.0</td><td>0:00:05</td><td>18</td><td>82</td><td>82.0</td><td>0:01:31</td><td>4</td><td>196</td><td>98.0</td><td>0:00:32</td><td>14</td><td>486</td><td>97.2</td><td>0:02:15</td></tr><tr><td>1e-8</td><td>5</td><td>5</td><td>50.0</td><td>0:00:05</td><td>36</td><td>64</td><td>64.0</td><td>0:02:05</td><td>45</td><td>155</td><td>77.5</td><td>0:05:16</td><td>92</td><td>408</td><td>81.6</td><td>0:21:44</td></tr><tr><td>1e-9</td><td>6</td><td>4</td><td>40.0</td><td>0:00:09</td><td>56</td><td>44</td><td>44.0</td><td>0:02:22</td><td>105</td><td>95</td><td>47.5</td><td>0:06:17</td><td>281</td><td>219</td><td>43.8</td><td>0:26:33</td></tr><tr><td>1e-10</td><td>9</td><td>1</td><td>10.0</td><td>0:00:09</td><td>63</td><td>37</td><td>37.0</td><td>0:02:47</td><td>116</td><td>84</td><td>42.0</td><td>0:08:23</td><td>295</td><td>205</td><td>41.0</td><td>0:33:23</td></tr><tr><td rowspan="4">nbody</td><td>1e-7</td><td>9</td><td>15</td><td>62.5</td><td>0:00:17</td><td>6</td><td>234</td><td>97.5</td><td>0:02:37</td><td>9</td><td>471</td><td>98.1</td><td>0:04:18</td><td>21</td><td>1179</td><td>98.2</td><td>0:16:04</td></tr><tr><td>1e-8</td><td>16</td><td>8</td><td>33.3</td><td>0:00:25</td><td>57</td><td>183</td><td>76.2</td><td>0:07:40</td><td>62</td><td>418</td><td>87.1</td><td>0:18:52</td><td>40</td><td>1160</td><td>96.7</td><td>0:28:19</td></tr><tr><td>1e-9</td><td>20</td><td>4</td><td>16.7</td><td>0:00:29</td><td>163</td><td>77</td><td>32.1</td><td>0:10:19</td><td>176</td><td>304</td><td>63.3</td><td>0:43:10</td><td>172</td><td>1028</td><td>85.7</td><td>0:30:19</td></tr><tr><td>1e-10</td><td>21</td><td>3</td><td>12.5</td><td>0:00:29</td><td>197</td><td>43</td><td>17.9</td><td>0:12:07</td><td>388</td><td>92</td><td>19.2</td><td>0:33:18</td><td>770</td><td>430</td><td>35.8</td><td>5:30:57</td></tr><tr><td rowspan="4">arclength</td><td>1e-8</td><td>6</td><td>5</td><td>45.5</td><td>0:00:10</td><td>29</td><td>81</td><td>73.6</td><td>0:01:13</td><td>44</td><td>176</td><td>80.0</td><td>0:03:20</td><td>86</td><td>464</td><td>84.4</td><td>0:11:22</td></tr><tr><td>1e-9</td><td>8</td><td>3</td><td>27.3</td><td>0:00:10</td><td>36</td><td>74</td><td>67.3</td><td>0:02:09</td><td>119</td><td>101</td><td>45.9</td><td>0:07:07</td><td>80</td><td>470</td><td>85.5</td><td>0:15:54</td></tr><tr><td>1e-10</td><td>9</td><td>2</td><td>18.2</td><td>0:00:10</td><td>95</td><td>15</td><td>13.6</td><td>0:02:25</td><td>168</td><td>52</td><td>23.6</td><td>0:05:51</td><td>391</td><td>159</td><td>28.9</td><td>0:53:48</td></tr><tr><td>1e-11</td><td>11</td><td>0</td><td>0.0</td><td>0:00:08</td><td>110</td><td>0</td><td>0.0</td><td>0:01:52</td><td>210</td><td>10</td><td>4.5</td><td>0:07:01</td><td>487</td><td>63</td><td>11.5</td><td>0:30:51</td></tr><tr><td rowspan="4">simpson</td><td>1e-7</td><td>4</td><td>6</td><td>60.0</td><td>0:00:05</td><td>16</td><td>84</td><td>84.0</td><td>0:00:51</td><td>20</td><td>180</td><td>90.0</td><td>0:01:04</td><td>41</td><td>459</td><td>91.8</td><td>0:03:53</td></tr><tr><td>1e-8</td><td>5</td><td>5</td><td>50.0</td><td>0:00:05</td><td>26</td><td>74</td><td>74.0</td><td>0:01:14</td><td>26</td><td>174</td><td>87.0</td><td>0:02:29</td><td>52</td><td>448</td><td>89.6</td><td>0:09:43</td></tr><tr><td>1e-9</td><td>6</td><td>4</td><td>40.0</td><td>0:00:09</td><td>56</td><td>44</td><td>44.0</td><td>0:01:17</td><td>74</td><td>126</td><td>63.0</td><td>0:03:38</td><td>165</td><td>335</td><td>67.0</td><td>0:20:59</td></tr><tr><td>1e-10</td><td>9</td><td>1</td><td>10.0</td><td>0:00:09</td><td>65</td><td>35</td><td>35.0</td><td>0:02:01</td><td>119</td><td>81</td><td>40.5</td><td>0:06:36</td><td>269</td><td>231</td><td>46.2</td><td>0:20:14</td></tr><tr><td rowspan="4">nbody</td><td>1e-7</td><td>9</td><td>15</td><td>62.5</td><td>0:00:17</td><td>7</td><td>233</td><td>97.1</td><td>0:01:31</td><td>18</td><td>462</td><td>96.2</td><td>0:06:26</td><td>16</td><td>1184</td><td>98.7</td><td>0:08:19</td></tr><tr><td>1e-8</td><td>16</td><td>8</td><td>33.3</td><td>0:00:25</td><td>136</td><td>104</td><td>43.3</td><td>0:06:46</td><td>77</td><td>403</td><td>84.0</td><td>0:22:23</td><td>143</td><td>1057</td><td>88.1</td><td>0:55:47</td></tr><tr><td>1e-9</td><td>20</td><td>4</td><td>16.7</td><td>0:00:29</td><td>177</td><td>63</td><td>26.2</td><td>0:10:57</td><td>285</td><td>195</td><td>40.6</td><td>0:46:33</td><td>301</td><td>899</td><td>74.9</td><td>2:45:53</td></tr><tr><td>1e-10</td><td>21</td><td>3</td><td>12.5</td><td>0:00:29</td><td>208</td><td>32</td><td>13.3</td><td>0:09:28</td><td>385</td><td>95</td><td>19.8</td><td>0:38:00</td><td>900</td><td>300</td><td>25.0</td><td>3:30:26</td></tr><tr><td rowspan="4">arclength</td><td>1e-11</td><td>24</td><td>0</td><td>0.0</td><td>0:00:29</td><td>225</td><td>15</td><td>6.2</td><td>0:09:32</td><td>447</td><td>33</td><td>6.9</td><td>0:36:41</td><td>1057</td><td>143</td><td>11.9</td><td>2:45:13</td></tr><tr><td>1e-12</td><td>24</td><td>0</td><td>0.0</td><td>0:00:30</td><td>238</td><td>2</td><td>0.8</td><td>0:06:59</td><td>468</td><td>12</td><td>2.5</td><td>0:33:17</td><td>1154</td><td>46</td><td>3.8</td><td>2:09:48</td></tr></tbody></table>

<table><caption>Table II: Auto-tuning results for unrolling strategy.</caption><thead><tr><th rowspan="2">Program</th><th rowspan="2">Threshold</th><th colspan="4">delta-debugging</th><th colspan="4">unroll / factor = 10</th><th colspan="4">unroll / factor = 20</th><th colspan="4">unroll / factor = 50</th></tr><tr><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th><th>b64</th><th>b32</th><th>%</th><th>h:m:s</th></tr></thead><tbody><tr><td rowspan="4">sum</td><td>1e-5</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>6</td><td>4</td><td>40.0</td><td>0:00:09</td><td>11</td><td>9</td><td>45.0</td><td>0:00:13</td><td>25</td><td>25</td><td>50.0</td><td>0:00:49</td></tr><tr><td>1e-6</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>10</td><td>0</td><td>0.0</td><td>0:00:08</td><td>19</td><td>1</td><td>5.0</td><td>0:00:17</td><td>45</td><td>5</td><td>10.0</td><td>0:00:52</td></tr><tr><td>1e-7</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>10</td><td>0</td><td>0.0</td><td>0:00:07</td><td>20</td><td>0</td><td>0.0</td><td>0:00:16</td><td>48</td><td>2</td><td>4.0</td><td>0:00:53</td></tr><tr><td>1e-8</td><td>1</td><td>0</td><td>0.0</td><td>0:00:01</td><td>10</td><td>0</td><td>0.0</td><td>0:00:08</td><td>20</td><td>0</td><td>0.0</td><td>0:00:17</td><td>48</td><td>2</td><td>4.0</td><td>0:00:57</td></tr><tr><td rowspan="4">riemann</td><td>1e-8</td><td>3</td><td>2</td><td>40.0</td><td>0:00:04</td><td>6</td><td>44</td><td>88.0</td><td>0:00:22</td><td>6</td><td>94</td><td>94.0</td><td>0:00:58</td><td>2</td><td>248</td><td>99.2</td><td>0:01:32</td></tr><tr><td>1e-9</td><td>3</td><td>2</td><td>40.0</td><td>0:00:04</td><td>26</td><td>24</td><td>48.0</td><td>0:00:29</td><td>48</td><td>52</td><td>52.0</td><td>0:01:16</td><td>42</td><td>208</td><td>83.2</td><td>0:03:58</td></tr><tr><td>1e-10</td><td>4</td><td>1</td><td>20.0</td><td>0:00:04</td><td>44</td><td>6</td><td>12.0</td><td>0:00:52</td><td>88</td><td>12</td><td>12.0</td><td>0:02:19</td><td>179</td><td>71</td><td>28.4</td><td>0:09:45</td></tr><tr><td>1e-11</td><td>5</td><td>0</td><td>0.0</td><td>0:00:04</td><td>44</td><td>6</td><td>12.0</td><td>0:00:54</td><td>88</td><td>12</td><td>12.0</td><td>0:02:13</td><td>104</td><td>146</td><td>58.4</td><td>0:08:53</td></tr><tr><td rowspan="4">arclength</td><td>1e-8</td><td>6</td><td>5</td><td>45.5</td><td>0:00:10</td><td>29</td><td>81</td><td>73.6</td><td>0:01:13</td><td>44</td><td>176</td><td>80.0</td><td>0:03:20</td><td>86</td><td>464</td><td>84.4</td><td>0:11:22</td></tr><tr><td>1e-9</td><td>8</td><td>3</td><td>27.3</td><td>0:00:10</td><td>36</td><td>74</td><td>67.3</td><td>0:02:09</td><td>119</td><td>101</td><td>45.9</td><td>0:07:07</td><td>80</td><td>470</td><td>85.5</td><td>0:15:54</td></tr><tr><td>1e-10</td><td>9</td><td>2</td><td>18.2</td><td>0:00:10</td><td>95</td><td>15</td><td>13.6</td><td>0:02:25</td><td>168</td><td>52</td><td>23.6</td><td>0:05:51</td><td>391</td><td>159</td><td>28.9</td><td>0:53:48</td></tr><tr><td>1e-11</td><td>11</td><td>0</td><td>0.0</td><td>0:00:08</td><td>110</td><td>0</td><td>0.0</td><td>0:01:52</td><td>210</td><td>10</td><td>4.5</td><td>0:07:01</td><td>487</td><td>63</td><td>11.5</td><td>0:30:51</td></tr><tr><td rowspan="4">simpson</td><td>1e-7</td><td>4</td><td>6</td><td>60.0</td><td>0:00:05</td><td>16</td><td>84</td><td>84.0</td><td>0:00:51</td><td>20</td><td>180</td><td>90.0</td><td>0:01:04</td><td>41</td><td>459</td><td>91.8</td><td>0:03:53</td></tr><tr><td>1e-8</td><td>5</td><td>5</td><td>50.0</td><td>0:00:05</td><td>26</td><td>74</td><td>74.0</td><td>0:01:14</td><td>26</td><td>174</td><td>87.0</td><td>0:02:29</td><td>52</td><td>448</td><td>89.6</td><td>0:09:43</td></tr><tr><td>1e-9</td><td>6</td><td>4</td><td>40.0</td><td>0:00:09</td><td>56</td><td>44</td><td>44.0</td><td>0:01:17</td><td>74</td><td>126</td><td>63.0</td><td>0:03:38</td><td>165</td><td>335</td><td>67.0</td><td>0:20:59</td></tr><tr><td>1e-10</td><td>9</td><td>1</td><td>10.0</td><td>0:00:09</td><td>65</td><td>35</td><td>35.0</td><td>0:02:01</td><td>119</td><td>81</td><td>40.5</td><td>0:06:36</td><td>269</td><td>231</td><td>46.2</td><td>0:20:14</td></tr><tr><td rowspan="4">nbody</td><td>1e-7</td><td>9</td><td>15</td><td>62.5</td><td>0:00:17</td><td>7</td><td>233</td><td>97.1</td><td>0:01:31</td><td>18</td><td>462</td><td>96.2</td><td>0:06:26</td><td>16</td><td>1184</td><td>98.7</td><td>0:08:19</td></tr><tr><td>1e-8</td><td>16</td><td>8</td><td>33.3</td><td>0:00:25</td><td>136</td><td>104</td><td>43.3</td><td>0:06:46</td><td>77</td><td>403</td><td>84.0</td><td>0:22:23</td><td>143</td><td>1057</td><td>88.1</td><td>0:55:47</td></tr><tr><td>1e-9</td><td>20</td><td>4</td><td>16.7</td><td>0:00:29</td><td>177</td><td>63</td><td>26.2</td><td>0:10:57</td><td>285</td><td>195</td><td>40.6</td><td>0:46:33</td><td>301</td><td>899</td><td>74.9</td><td>2:45:53</td></tr><tr><td>1e-10</td><td>21</td><td>3</td><td>12.5</td><td>0:00:29</td><td>208</td><td>32</td><td>13.3</td><td>0:09:28</td><td>385</td><td>95</td><td>19.8</td><td>0:38:00</td><td>900</td><td>300</td><td>25.0</td><td>3:30:26</td></tr><tr><td rowspan="4">arclength</td><td>1e-11</td><td>24</td><td>0</td><td>0.0</td><td>0:00:29</td><td>225</td><td>15</td><td>6.2</td><td>0:09:32</td><td>447</td><td>33</td><td>6.9</td><td>0:36:41</td><td>1057</td><td>143</td><td>11.9</td><td>2:45:13</td></tr><tr><td>1e-12</td><td>24</td><td>0</td><td>0.0</td><td>0:00:30</td><td>238</td><td>2</td><td>0.8</td><td>0:06:59</td><td>468</td><td>12</td><td>2.5</td><td>0:33:17</td><td>1154</td><td>46</td><td>3.8</td><td>2:09:48</td></tr></tbody></table>

For this last experiment, we use the 494\_bus matrix from the Suite Sparse Matrix Collection [26], where *b* is so that the solution *x* is a vector of ones, and we ran the algorithm with a threshold of $10^{-6}$ on the residual. Implemented in binary64, it requires 1315 iterations to converge, while in binary32 this number of iterations scales to 2494. The objective was to find in the binary64 version of the program some iterations that could be transformed into binary32 so that the total number of iterations remains no larger than a user-defined bound.

For doing this, we ran our tool by forcing all the instructions of the main loop to be in the same format, as in Section IV-B. And we repeated this for a splitting factor ranging from 2 to 20. Figure 5 shows the results for a number of iterations bounded by 1500 and 1800.

Our tool allows to find automatically some iterations to be transformed in binary32 arithmetic. For example, for a splitting factor of 2, we obtained two loops in binary64 arithmetic, of 657 and 658 iterations, respectively. Then if we bound the

total number of iterations to 1500 (Figure 5a), no iteration can be lowered in binary32, and all the 1315 iterations remain in binary64 (in red in the figures). However, if we bound the number of iterations to 1800 (Figure 5b), the tool found that the second loop can be performed in binary32 (in blue in the figures). However, since the computation precision impacts the number of iterations, lowering the precision of instructions in the second loop makes its number of iterations increased. It scales to 1046, for a total of 1703 iterations, which is still below the bound.

## V. CONCLUSION

This article addresses the tuning of precision in iterative routines. We have designed an auto-tuning tool that relies on compiler capabilities and not on the user action to improve the tuning process. Particularly, it applies static loop transformations on the original program to increase the number of instructions appearing in the program and thus to yield
<figure><img src="image_7.png" /><figcaption>(a) number of iterations bounded by 1500</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(b) number of iterations bounded by 1800</figcaption></figure>

<caption>Figure 5: Increase of iteration numbers for the *conjugate gradient* and different splitting factors.</caption>

better solutions. We showed on some examples, that it helps to improve significantly the existing solutions, even in the case of unbounded loops, to the detriment of an increase in the tuning time. A direct extension would be to give the capability to user of combining both these loop transformations.

Now the research direction is twofold: First, other static loop transformations exist, and it could be interesting to investigate their interest. For example, the loop peeling and versioning lead also to a duplication of the instructions of the loop, and they could be useful for our auto-tuning approach. Second, we have applied our techniques on small routines embracing simple loops or, in the case of the *conjugate gradient*, on the main loop. It could be interesting to study how this approach scales with the size of the loop body since instructions are duplicated, and how it behaves on nested loops.

## ACKNOWLEDGMENT

This work was supported by the PADOC ANR project under grant n°ANR-18-CE25-0004.

## REFERENCES

[1] IEEE Computer Society, *IEEE Standard for Floating-Point Arithmetic*, 2019.

[2] M. Baboulin, A. Buttari, J. Dongarra, J. Kurzak, J. Langou, J. Langou, P. Luszczek, and S. Tomov, "Accelerating scientific computations with mixed precision algorithms," *Computer Physics Communications*, vol. 180, no. 12, pp. 2526–2533, 2009.

[3] A. Izycheva and E. Darulova, "On Sound Relative Error Bounds for Floating-Point Arithmetic," in *17th Conference on Formal Methods in Computer-Aided Design*, 2017, p. 15–22.

[4] E. Darulova, A. Izycheva, F. Nasir, F. Ritter, H. Becker, and R. Bastian, "Daisy - Framework for Analysis and Optimization of Numerical Programs," in *Tools and Algorithms for the Construction and Analysis of Systems*, D. Beyer and M. Huisman, Eds., 2018, pp. 270–287.

[5] W.-F. Chiang, M. Baranowski, I. Briggs, A. Solovyev, G. Gopalakrishnan, and Z. Rakamarić, "Rigorous Floating-Point Mixed-Precision Tuning," in *POPL 2017*, 2017, pp. 300–315.

[6] A. Solovyev, C. Jacobsen, Z. Rakamarić, and G. Gopalakrishnan, "Rigorous Estimation of Floating-Point Round-off Errors with Symbolic Taylor Expansions," in *FM 2015*, 2015, pp. 532–550.

[7] M. Martel, "Floating-Point Format Inference in Mixed-Precision," in *9th NASA Formal Methods Symposium*, 2017.

[8] M. O. Lam, J. K. Hollingsworth, B. R. de Supinski, and M. P. Legendre, "Automatically adapting programs for mixed-precision floating-point computation," in *ICS'13*, 2013, pp. 369–378.

[9] M. O. Lam and J. K. Hollingsworth, "Fine-grained floating-point precision analysis," vol. 32, no. 2, 2018, pp. 231–245.

[10] C. Rubio-González, C. Nguyen, H. D. Nguyen, J. Demmel, W. Kahan, K. Sen, D. H. Bailey, C. Iancu, and D. Hough, "Precimonious: tuning assistant for floating-point precision," in *SC'13*, 2013, pp. 1–12.

[11] A. Zeller and R. Hildebrandt, "Simplifying and Isolating Failure-Inducing Input," *IEEE Transactions on Software Engineering*, pp. 183–200, 2002.

[12] C. Rubio-González, C. Nguyen, B. Mehne, K. Sen, J. Demmel, W. Kahan, C. Iancu, W. Lavrijsen, D. H. Bailey, and D. Hough, "Floating-point precision tuning using Blame analysis," in *ICSE 2016*, 2016, pp. 1074–1085.

[13] H. Guo and C. Rubio-González, "Exploiting Community Structure for Floating-Point Precision Tuning," in *27th ACM SIGSOFT International Symposium on Software Testing and Analysis*. Association for Computing Machinery, 2018, p. 333–343.

[14] S. Graillat, F. Jézéquel, R. Picot, F. Févotte, and B. Lathuilière, "Autotuning for floating-point precision with Discrete Stochastic Arithmetic," *Journal of Computational Science*, vol. 36, p. 101017, 2019.

[15] F. Févotte and B. Lathuilière, "Debugging and optimization of HPC programs with the Verrou tool," in *International Workshop on Software Correctness for HPC Applications (Correctness)*, 2019, pp. 1–10.

[16] F. Jézéquel and J.-M. Chesneaux, "CADNA: a library for estimating round-off error propagation," *Computer Physics Communications*, pp. 933–955, 2008.

[17] H. Menon, M. O. Lam, D. Osei-Kuffuor, M. Schordan, S. Lloyd, K. Mohror, and J. Hittinger, "ADAPT: Algorithmic Differentiation Applied to Floating-Point Precision Tuning," in *International Conference for High Performance Computing, Networking, Storage, and Analysis*. IEEE Press, 2018.

[18] B. Saiki, O. Flatt, C. Nandi, P. Panchekha, and Z. Tatlock, "Combining Precision Tuning and Rewriting," in *28th IEEE Symposium on Computer Arithmetic (ARITH 2021)*, 2021, pp. 1–8.

[19] P. Panchekha, A. Sanchez-Stern, J. R. Wilcox, and Z. Tatlock, "Automatically improving accuracy for floating point expressions," in *PLDI'2015*, 2015, pp. 1–11.

[20] Y. Chatelain, E. Petit, P. De Oliveira Castro, G. Lartigue, and D. Defour, "Automatic Exploration of Reduced Floating-Point Representations in Iterative Methods," in *25th International Conference Euro-Par 2019 Parallel Processing*, 2019, pp. 481–494.

[21] H. Brunie, C. Iancu, K. Z. Ibrahim, P. Brisk, and B. Cook, "Tuning Floating-Point Precision Using Dynamic Program Information and Temporal Locality," in *SC20: International Conference for High Performance Computing, Networking, Storage and Analysis*, 2020, pp. 1–14.

[22] R. Gu, P. Beata, and M. Becchi, "A Loop-Aware Autotuner for High-Precision Floating-Point Applications," in *International Symposium on Performance Analysis of Systems and Software (ISPASS 2020)*, 2020, pp. 285–295.

[23] G. Revy, "Analyzing the impact of floating-point precision adaptation in iterative programs," in *2021 IEEE 28th Symposium on Computer Arithmetic (ARITH)*, 2021, pp. 25–32.

[24] L. Kirschner, E. Soremekun, and A. Zeller, "Debugging Inputs," in *2020 IEEE/ACM 42nd International Conference on Software Engineering (ICSE)*, 2020, pp. 75–86.

[25] D. Ben Khalifa and M. Martel, "A study of the floating-point tuning behaviour on the n-body problem," in *Computational Science and Its*