*R* := [*PARAMS*] → {[*input tuple*] → [*output tuple*] | *constraints*}, where [*PARAMS*] is the list of relation parameters, [*input tuple*] represents dependence sources, [*output tuple*] represents dependence destinations, and *constraints* are the constraints—a system of affine equalities and inequalities on parameters and tuple variables.

In the case of dependencies when the dimensions of the left and right tuples of the relation are the same, the distance vector is the difference between the iteration vector of a dependence target and the iteration vector of the corresponding dependence source.

A dependence vector is uniform if all its elements are constants.

To extract dependencies present in the loop nest, we use the polyhedral model, which is returned by PET [5], and we use the iscc calculator [6], which performs calculations on polyhedral sets and relations [7]. The iscc is an interactive interface to the PET library and barvinok library that will let you count points in polytopes. Barvinok is available online at https://repo.or.cz/barvinok.git (accessed on 24 August 2021). We also used the iscc calculator to calculate distance vectors and generate target code.

To parallelize and tile the loop nest, a time partitioning constraint should be created [1] that states that if iteration *I* of statement *S*1 depends on iteration *J* of statement *S*2, then *I* must be assigned to a time partition that is executed no earlier than the partition containing *J*, i.e., schedule(*I*) ≤ schedule(*J*), where schedule(*I*) and schedule(*J*) denote the discrete execution time of iterations *I* and *J*, respectively.

Linear independent solutions to time partition constraints are needed to create schedules for each occurence of a single instruction of the loop nest allowing for parallelization and tiling of code. The schedule defines a strict partial order, i.e., an irreflexive and transitive relation on the statement instances that determines the order in which they are or should be executed. Details of use of linear independent schedules for generating parallel tiled code can be found in a number of articles, for example, in article [2].

We should extract as many linear independent solutions to time partition constraints as possible. The degrees of parallelism of the target code and the dimension of the tile are higher when more independent solutions are extracted [1]. When there is a single solution to the time partition constraints, parallelization and tiling of the corresponding loop nest using affine transformations is not possible [1].

#### **3. Methods. Parallel Tiled Code Generation**

Using PET and the iscc calculator, we extract dependencies available in the code in Listing 1; they are presented with the following relation.

*R* := (*loop*, *n*) → { (*l*, *i*, *k*) → (*l* - , *i*, *k*- ) | *l* > 0 ∧ *i* < *n* ∧ 0 ≤ *k* < *i* ∧ *l* < *l* - ≤ *loop* ∧ 0 ≤ *k*- < *i* } ∪ (*loop*, *n*) → { (*l*, *i*, *k*) → (*l* - , −1 + *i* − *k*, *k*- ) | *l* > 0 ∧ *i* < *n* ∧ *k* ≥ 0 ∧ *l* < *l* - ≤ *loop* ∧ 0 ≤ *k*- ≤ −2 + *i* − *k* } ∪ (*loop*, *n*) → { (*l*, *i*, *k*) → (*l* - , *i* - , −1 − *i* + *i* - ) | *l* > 0 ∧ 0 ≤ *k* < *i* ∧ *l* ≤ *l* - ≤ *loop* ∧ *i* < *i* - < *n* } ∪ (*loop*, *n*) → { (*l*, *i*, *k*) → (*l*, *i*, *k*- ) | 0 < *l* ≤ *loop* ∧ *i* < *n* ∧ *k* ≥ 0 ∧ *k* < *k*-< *i* },

where *R* is the relation name; *loop* and *n* are parameters; ∪ is the union operation of sets (relation *R* is composed as the set union of simpler relations); the tuple before the sign → of each simpler relation is the left tuple of this relation, for example, for the first simpler relation, the left tuple is presented with variables (*l*, *i*, *k*); the tuple after the sign → of each simpler relation is the right tuple of this relation, for example, for the first simpler relation, the right tuple is presented with variables (*l* - , *i*, *k*- ); the expressions after the sign | are the constraints of each simpler relation, each constraint is represented with the conjunctions of inequalities built on tuple variables and parameters; and ∧ is the logical AND operator.

The left tuple of each simpler relation represents dependence sources, whereas the right represents dependence destinations.

Applying the *deltas* operator of the iscc calculator to relation *R*, we obtain the three distance vectors presented with set *D* below.

$$D := (loop, \stackrel{\circ}{n}) \to \{ (l, i, k) \mid 0 \le l < loop \land ((l > 0 \land i < 0 \land i < k < n + 2i) \lor 1) \}$$

(*i* > 0 ∧ −*n* + 2*i* < *k* < *i*)) } ∪

(*loop*, *n*) → { (0, 0, *k*) | *loop* > 0 ∧ 0 < *k* ≤ −2 + *n* } ∪

(*loop*, *n*) → { (*l*, 0, *k*) | 0 < *l* < *loop* ∧ 2 − *n* ≤ *k* ≤ −2 + *n* }.

where the notations used are the same as for relation *R* above except from the set is represented with a single tuple.

Each conjunct in the set above represents a particular distance vector.

Taking into account the constraints of those distance vectors, we simplify them to the following form.

$$\begin{array}{lcl} D := \{(a\_1, a\_2, a\_3) \mid a\_1 \ge 0 \land -\infty \le a\_2 \le \infty \land -\infty \le a\_3 \le \infty\} \\ (0, 0, b\_3) \mid b\_3 > 0; \\ (1, -0, a\_2) \mid a\_1 \ge 0 \land \dots < a\_2 \le \infty. \end{array}$$

(*c*1, 0, *c*3) | *c*1 > 0 ∧ −∞ ≤ *c*3 ≤ ∞}.

The time partition constraints created from the resulting distance vectors according to article [1] are as follows.

$$h\_1 \* a\_1 + h\_2 \* a\_2 + h\_3 \* a\_3 \ge 0,\tag{1}$$

$$h\_3 \* b\_{\hat{3}} \ge 0,\tag{2}$$

$$h\_1 \* c\_1 + h\_3 \* c\_3 \ge 0,\tag{3}$$

where *h*1, *h*2, *h*3 are the unknowns.

Taking into consideration that −∞ ≤ *a*2 ≤ ∞, −∞ ≤ *a*3 ≤ ∞, −∞ ≤ *c*3 ≤ ∞, we can suppose that to satisfy all the above constraints, *h*2 and *h*3 should be 0, i.e., *h*2 = *h*3 = 0. Thus, the above constraints can be rewritten as follows.

$$h\_1 \* a\_1 \ge 0,\tag{4}$$

$$h\_1 \* c\_1 \ge 0.\tag{5}$$

Hence, we may cease that there exists a single solution to constraints (1), (2), and (3), namely (1, 0, 0)*T*. This means that all the three loops in the code in Listing 1 cannot be parallelized and tiled by means of affine transformations.

Next, we try to parallelize and tile only two inner loops *i* and *k* in the loop nest in Listing 1. For this purpose, we make the outermost loop *l* to be serial and extract dependencies for inner loops *i* and *k* described with the relation below.

*R* := (*n*) → { (*i*, *k*) → (*i*, *k*- ) | *i* < *n* ∧ *k* ≥ 0 ∧ *k* < *k*- < *i* } ∪ *n* → { (*i*, *k*) → (*i* - , −1 − *i* + *i* - ) | 0 ≤ *k* < *i* ∧ *i* < *i* -< *n* }.

Applying the *deltas* operator of the iscc calculator to relation *R*, we obtain the two distance vectors presented with set *D* below.

*D* := (*n*) → { (*i*, *k*) | *i* > 0 ∧ −*n* + 2*i* < *k* < *i* } ∪

(*n*) → { (0, *k*) | 0 < *k* ≤ −2 + *n* }.

Next, we simplify the representation of the distance vectors above to the form.

*D* := {(*a*1, *a*2) | *a*1 > 0 ∧ −∞ ≤ *a*2 ≤ ∞;

$$\left(\underline{0,b\_2}\right)\mid b\_2 > \underline{0}\}.$$

The time partition constraints formed on the basis of the distance vectors above are the following.

$$\begin{array}{c} h\_1 \ast a\_1 + h\_2 \ast a\_2 \geq 0, \\ h\_2 \ast b\_2 \geq 0, \end{array} \tag{6}$$

where *h*1, *h*2 are the unknowns.

Taking into consideration that −∞ ≤ *a*2 ≤ ∞, we can deduce that to satisfy constraints (6) and (7), *h*2 should be 0, i.e., *h*2 = 0.

So, there exists a single solution to constraints (6) and (7), namely (1, 0)*T*, and we conclude that provided the outermost loop *l* is serial, the two inner loops *i* and *k* cannot be parallelized and tiled by means of affine transformations.

To cope with that problem, we transform the code in Listing 1 to improve dependence properties. With this goal, we apply the following schedule to each iteration of the code in Listing 1:

(*l*, *i*, *k*)*T* → (*l*, *t* = *i* − *k*)*T*.

This schedule implies that each iteration of the code in Listing 1, represented with iteration vector (*l*, *i*, *k*)*T*, is mapped to the two-dimensional time (*l*, *t* = *i* − *k*)*T*. It means that iterations of loop *l* should be executed serially, while for a given value of iterator *l*, iteration (*i*, *k*)*T* should be executed at time *t* = *i* − *k*. This time guarantees that each iteration (*i*, *k*)*T* is executed when all its operands are ready. To justify that fact, let us noting that for iteration (*i*, *k*)*T*, operand *b*[*k*][*i*] is input data; hence, its value is ready at time 0, and operand *w*[(*i* − *k*) − 1] is ready at time (*i* − *k*) − 1 when an actual value of this operand is already calculated and written in memory. Thus, iteration (*i*, *k*)*T* can be executed at time *t* = *i* − *k*, i.e., at time, which is one more than the time when operand *w*[(*i* − *k*) − 1] is ready.

In other words, the schedule above is based on data flow software paradigm [8].

To generate target serial code, we form the following relation, which maps each statement instance within the iteration space of the code in Listing 1 to the two-dimensional schedule below.

*CODE* := (*loop*, *n*) → { (*l*, *i*, *k*) → (*l*, *t* = *i* − *k*) | *loop* > 0 ∧ 0 < *l* ≤ *loop* ∧ 0 < *i* < *n* ∧ 0 ≤ *k* < *i* },

where the constraints

*loop* > 0 ∧ 0 < *l* ≤ *loop* ∧ 0 < *i* < *n* ∧ 0 ≤ *k* < *i*

define the iteration space of the code in Listing 1. Applying the *iscc* codegen operator to the relation above, we get the pseudocode shown in Listing 2.

Listing 2: Target serial pseudocode.

```
1
2 for ( int counter = 1; counter <= loop ; counter += 1 )
3 for ( int var1 = 1; var1 < n ; var1 += 1 )
4 for ( int var2 = var1 ; var2 < n ; var2 += 1 )
5 ( counter , var2 , −var1 + var2 ) ; // pseudostatement
```
We transform the pseudocode code in Listing 2 to C code, taking into account that in that pseudocode, variables *counter*, *var*1, and *var*2 correspond to variables *l*, *t*, and *i*, respectively, in the tuple of set *CODE*; the second variable *var*2 in the pseudostatement relates to variable *i*, while the third expression −*var*1 + *var*2 corresponds to variable *k* in the tuple of set *CODE*. Thus, we replace the pseudostatement in the code in Listing 2 with the statement

*w*[*i*]+ = *b*[*k*][*i*] ∗ *w*[(*i* − *k*) − 1];

from Listing 1 changing variables *i* and *k* with variable *var*2 and the expression −*var*1 + *var*2, respectively. As a result, we obtain the compilable program fragment presented in Listing 3.

Listing 3: Target sequential compilable program fragment.

```
1 for ( int counter = 1; counter <= loop ; counter += 1 )
2 for ( int var1 = 1; var1 < n ; var1 += 1 )
3 for ( int var2 = var1 ; var2 < n ; var2 += 1 )
4 w[ var2 ] += b[ − var1 + var2 ] [ var2 ] * w[ ( var2 − (− var1 +
             var2 ) ) − 1];
```
The target serial code in Listing 3 is in the scope of reordered transformations. It performs the same computations as those performed with the initial code in Listing 1 but in a different order. It is well-known that a reordered transformation of a code is correct if it executes the same computations as those executed with the initial one (1) and respects all the dependencies that appear in that code (2) [1]. The transformed code is correct as it performs the same computations as those executed with the initial one (1) and it respects all the dependencies available in the initial one as explained below (2).

There exist three kinds of dependencies in the code presented in Listing 3: data flow dependencies (some statement instance first generates a result, then that result is used with another statement instance, those instances belong to different time units represented with the value of iterator counter), antidependencies (some statement instance first reads a result, then that result is updated with another statement instance, and output dependencies (two statement instances write their results to the same memory location).

Data flow dependencies are respected due to the fact that in the target code, the execution of a statement instance being the target of each data dependence starts only when all the arguments (data) of this operation are prepared, i.e., the processing of all the instruction instances generating these arguments has already finished, and the operand values are stored in the shared part of memory. This is guaranteed because the source of each data dependence is executed at a time unit defined with the value of iterator counter that is less than the one when the corresponding target is executed.

Anti- and output dependencies are honored due to the lexicographical order of the execution of dependent statement instances within each time partition represented with the value of iterator var1.

We also experimentally confirmed that the both loop nests presented in Listing 1 and Listing 3 generate correct results. The experiments used for input data prepared deterministically and randomly.

Dependencies available in the code in Listing 3 are represented with the following relation.

*R* := (*loop*, *n*) → { (*counter*, *var*1, *var*2) → (*counter*- , 1 + *var*2, *var*2- ) | *counter* > 0∧ *var*1 > 0∧ *var*2 ≥ *var*1∧ *counter* ≤ *counter*- ≤ *loop* ∧ *var*2 < *var*2- < *n* } ∪ (*loop*, *n*) → { (*counter*, *var*1, *var*2) → (*counter*- , *var*1- , *var*2) | *counter* > 0 ∧ *var*1 > 0 ∧ *var*1 ≤ *var*2 < *n*∧*counter* < *counter*- ≤ *loop*∧0 < *var*1- ≤ *var*2 } ∪(*loop*, *n*) → { (*counter*, *var*1, *var*2) → (*counter*- , *var*1- , −1 + *var*1) | *counter* > 0 ∧ *var*1 ≤ *var*2 < *n* ∧ *counter* < *counter*- ≤ *loop* ∧ 0 < *var*1- < *var*1 } ∪ (*loop*, *n*) → { (*counter*, *var*1, *var*2) → (*counter*, *var*1- , *var*2) | 0 < *counter* ≤ *loop* ∧ *var*1 > 0 ∧ *var*2 < *n* ∧ *var*1 < *var*1-≤ *var*2 },

where *loop* and *n* are parameters.

Applying the *deltas* operator of the iscc calculator to relation *R*, we obtain the three distance vectors presented below.

*D* := (*loop*, *n*) → { (*counter*, *var*1, *var*2) | 0 ≤ *counter* < *loop* ∧ ((*var*1 > 0 ∧ 0 < *var*2 < *n* − *var*1) ∨

(*counter* > 0 ∧ *var*1 < 0 ∧ −*n* − *var*1 < *var*2 < 0)) } ∪

(*loop*, *n*) → { (0, *var*1, 0) | *loop* > 0 ∧ 0 < *var*1 ≤ −2 + *n* } ∪

(*loop*, *n*) → { (*counter*, *var*1, 0) | 0 < *counter* < *loop* ∧ 2 − *n* ≤ *var*1 ≤ −2 + *n* }.

Taking into account the constraints of those distance vectors, we simplify their representation to the following form.

*D* := {(*a*1, *a*2, *a*3) | *a*1 ≥ 0 ∧ −∞ ≤ *a*2 ≤ ∞ ∧ −∞ ≤ *a*3 ≤ ∞;

(0, *b*2, 0) | *b*2 > 0;

$$\left(\underline{c\_1}, \underline{c\_2}, 0\right) \mid c\_1 > 0 \land -\infty \le c\_2 \le \infty\\\underline{\qquad}$$

The time partition constraints constructed according to article [1] are as follows.

$$h\_1 \* a\_1 + h\_2 \* a\_2 + h\_3 \* a\_3 \ge 0,\tag{8}$$

$$h\_2 \* b\_2 \ge 0,\tag{9}$$

$$h\_1 \* c\_1 + h\_2 \* c\_2 \ge 0,\tag{10}$$

where *h*1, *h*2, *h*3 are the unknowns.

Taking into account that −∞ ≤ *a*2 ≤ ∞ , −∞ ≤ *a*3 ≤ ∞, −∞ ≤ *c*2 ≤ ∞, we can deduce that *h*2 and *h*3 should be 0, i.e., *h*2 = *h*3 = 0 for the constraints (8), (9), and (10) to be compatible. Therefore, these constraints can be written with the following formulas.

$$h\_1 \* a\_1 \ge 0\_{\text{'}}\tag{11}$$

$$h\_1 \* c\_1 \ge 0.\tag{12}$$

Thus, we may conclude that there exists a single solution to constraints (8), (9), and (10), namely (1, 0, 0)*T*. This means that all thee loops in the code in Listing 3 cannot be parallelized and tiled by means of affine transformations.

Next, we try to parallelize and tile only two inner loops *var*1 and *var*2 in the loop nest presented in Listing 3. For this purpose, we make the outermost loop *counter* to be serial and extract dependencies for inner loops *var*1 and *var*2. They are expressed with the relation below.

*R* := (*n*) → { (*var*1, *var*2) → (1 + *var*2, *var*2- ) | *var*1 > 0 ∧ *var*2 ≥ *var*1 ∧ *var*2 < *var*2- < *n* } ∪ *n* → { (*var*1, *var*2) → (*var*1- , *var*2) | *var*1 > 0 ∧ *var*2 < *n* ∧ *var*1 < *var*1- ≤ *var*2 }.

Applying the *deltas* operator of the iscc calculator to relation *R*, we obtain the two distance vectors presented below.

*D* := (*n*) → { (*var*1, *var*2) | *var*1 > 0 ∧ 0 < *var*2 < *n* − *var*1 } ∪

*n* → { (*var*1, 0) | 0 < *var*1 ≤ −2 + *n* }.

After the simplification of the representation of the distance vector above, we obtain the following vectors.

$$D := \{ (a\_1, a\_2) \mid a\_1 > 0 \land a\_2 > 0; \}$$

(*b*1, 0) | *b*1 > 0}

The time partition constraints created from the distance vectors above are the following:

$$h\_1 \* a\_1 + h\_2 \* a\_2 \ge 0,\tag{13}$$

$$h\_1 \* b\_1 \ge 0,\tag{14}$$

where *h*1, *h*2 are the unknowns. There are two linear independent solutions to the constraints above: (1, 0)*T* and (0, 1)*T*. Applying those solutions, we are able to parallelize and tile the two inner loops of the code in Listing 3 using the technique presented in paper [2].

The target parallel tiled code presented by means of the OpenMP C/C++ API is shown in Listing 4. It is generated for the best tile size equal to 24 × 54; choosing the best tile size is explained in Section 5.

Listing 4: Transformed parallel loop nest.

1 2 **#define** min ( lhs , rhs ) ( ( lh s ) < ( rhs ) ? ( lh s ) : ( rhs ) ) 3 **#define** max ( lhs , rhs ) ( ( l h s ) > ( rhs ) ? ( l h s ) : ( rhs ) ) 4 **#define** floord ( val ,d) ( ( ( val ) <0) ? −(( −( val ) +(d) −1)/(d) ) : ( val ) /(d) ) 5 **#define** ceild ( val ,d) ceil ( ( ( **double** ) ( val ) ) /(( **double** ) (d) ) ) 6 7 **for** ( **int** i 0 = 1; i 0 <= loop ; i 0 += 1 ) { 8 **for** ( **int** w0 = 0; w0 <= floord (26 \*n−26, 675) ; w0+=1) { 9 #pragma omp p a r a l l e l **for** 10 **for** ( **int** h0 = max(0 , w0 − (n + 49) / 50 + 1) ; h0 <= min ( ( n − 1) / 54, (25 \* w0 + 24) / 52) ; h0 += 1 ) { 11 **for** ( **int** i 1 = max(1 , 54 \* h0 ) ; i 1 <= min(min(n − 1 , 50 \* w0 − 50 \* h0 + 49) , 54 \* h0 + 53) ; i 1 += 1 ) { 12 **for** ( **int** i 2 = max(50 \* w0 − 50 \* h0 , i 1 ) ; i 2 <= min( n − 1 , 50 \* w0 − 50 \* h0 + 49) ; i 2 += 1 ) { 13 w[ i 2 ] += ( b[ − i 1 + i 2 ] [ i 2 ] \* w[ i 1 − 1]) ; 14 } 15 }

```
16 }
17 }
18 }
```
In that code, outermost loop *i*0 is serial nontiled, and loops *w*0 and *h*0 enumerate tile identifiers, while loops *i*1 and *i*2 enumerate iterations within each tile. Parallelism is extracted with the wavefront technique [2] and presented with the OpenMP directive #*pragma omp parallel f or* inserted before loop *h*0 that means that this loop is parallel.

#### **4. Related Work**

Related techniques can be divided into the following two classes: the class of reordering transformations and the one of nonreordering transformations. There are numerous publications concerned with both of the classes. Approaches based on affine transformations [1,2,9–11] and those based on the transitive closure of dependence graphs [12–16] belong to reordering transformations. Reordering techniques are code-independent and are used in optimizing compilers, for example [17–19], which automatically generate optimized target code for source code.

Nonreordering transformations are code-dependent, i.e., for a given code, a transformation is fulfilled manually. The following publications within nonreordering transformations concern the problem similar to that implemented with the code in Listing 1 but not exactly the same problem [3,4,8,20–26].

Both classes allow for generating the target program that is semantically identical to the original one. However, there are the following differences in target code generated using techniques of those classes.

Reordering transformations do not introduce any additional computations to generated target code in comparison with those of original code; they only reorder loop nest iterations of the original code allowing for tiling and/or parallelism. They are codeindependent and are aimed at automatic code generation.

Nonreordering transformations allow us to parallelize and tile code similar to some extent to the code in Listing 1 but not exactly the same. All of those techniques introduce additional computations to generated target code in comparison with those in the original code. That increases the computational complexity of the algorithm and prevents achieving the maximal target code performance, and it is the main drawback in comparison with reordering transformations.

Each technique is manually created for the code that should be optimized. Adapting such a technique even to a slightly different problem can require additional work that can be time-consuming and not always possible.

After an extensive analysis of many nonreordering techniques mentioned above, we did not find any one that exactly implements the problem presented with the code in Listing 1. Without extensive research, it is not clear how any of those techniques can be adapted to implement exactly the same problem that implements the code in Listing 1.

In the class of reordering transformations, we examined the PLUTO [17] and TRACO compilers [15]. PLUTO is based on affine transformations and automatically generates tiled and/or parallel code. TRACO uses the transitive closure of dependence graphs to tile and/or parallelize input code. For the code in Listing 1, both PLUTO and TRACO are unable to generate any tiled and/or parallel code. In Section 3, we presented the reason why affine transformations fail to generate any parallel and/or tiled code for the serial code in Listing 1.

Below, we discuss some nonreordering transformations, which allow for generation of code implementing algorithms similar to that realizing with the code in Listing 1 but not exactly the same. Without extensive research, it is not clear how to adapt any of those techniques to generate target code fulfilling the same calculations as those performed with the code in Listing 1.

Karp et al. [20] discussed parallelism in recurrence equations. They proposed a decomposition algorithm that decides if a system of uniform recurrence equations (SURE) is computable or not. If so, multidimensional schedules can be derived and applied to extract parallelism.

Papers [3,4] introduced a recursive doubling strategy to compute recurrence equations in parallel. Recursive doubling envisages the splitting of the computation of a function into two subfunctions whose evaluation can be performed simultaneously in two separate processors. Successive splitting of each of these subfunctions allows for the computation over more processors.

Maleki and Burtscher [21] introduced two phase approach to compute recurrence equations. The first phase iteratively merges pairs of adjacent chunks by correcting the values in the second chunk of each pair. The second phase produces the resulting chunks in a pipelined mode to compute the final solution.

Sung et al. [22,23] proposed the idea to divide the input into blocks and decompose the computation over each block. Interblock parallelism is exploited to enhance code performance.

Nehab et al. [24] also suggested splitting the input data into blocks that are processed in parallel by modern GPU architectures and overlapped the causal, anticausal, row, and column filter processing.

Marongiu and Palazzari [25] addressed the parallelization of a class of iterative algorithms described as the system of affine recurrence equations (SARE). It introduces an affine timing function and an affine allocation function that perform a space-time transformation of the loop nest iteration space. It considers algorithms dealing with only uniform dependence vectors, while the approach presented in this paper deals with nonuniform vectors.

Ben-Asher and Haber [26] defined recurrence equations called "simple indexed recurrences" (SIR). In this type of equation, for extending capabilities, ordinary recurrences are generalized to *X*[*g*(*i*)] = *opi*(*X*[ *f*(*i*)], *X*[*g*(*i*)]), where *f* and *g* are affine functions *opi*(*x*, *y*) is a binary associative operator. In that paper, the authors proposed a parallel solution to the SIR problem. This case of recurrences is simpler than that considered in our paper and any tiled code is not considered.

Summing up, we may conclude that in the class of reordering transformations, there does not exist any technique allowing for parallelizing and/or tiling the code in Listing 1. In the class of nonordering transformations, to our best knowledge, no technique has been published to generate parallel tiled code implementing the problem addressed in this paper.

The main contribution of our paper is presenting a novel technique, which for the first time allows us to parallelize and tile the examined loop nest implementing computing general linear recurrence equations by means of reordering transformations. The novelty consists in adding an additional phase to classical reordering transformations: to source code, we first apply a reordering schedule that respects all data dependencies; then, we apply classical affine transformations to the serial code obtained in the first phase. This increases target code generation time but does not introduce any additional computations to the source code. Generated target code is still within the class of reordering transformations.

#### **5. Results**

The primary reason for writing a parallel program is speed. We strive that the parallel program execution should be completed at a shorter time in comparison with that of the serial one. We need to know what is the benefit from tiling and parallelism. For this purpose, we need to compute the parallel program speedup.

The speedup of a parallel program over a corresponding sequential program is the ratio of the compute time for the sequential program to the time for the parallel program. The value of speedup shows how efficient is a parallel program.

Perfect linear speedup occurs when the value of speedup is the same as the number of threads used for running a parallel program. In practice, perfect linear speedup seldom occurs because of parallel program overhead and the fact that all computations of an original program cannot be parallelized.

According to Amdahl's law, the parallel code speedup, *S*, is limited to *S* <= 1/*s*, where *s* is the serial fraction of code, i.e., the fraction of code that cannot be parallelized. For example, if *s* = 0.2, the maximal speedup is 5 regardless of the number of threads used for running a parallel program.

To evaluate the performance of the parallel tiled code presented in Listing 4, we carried out experiments aimed at measuring the execution time of the original program and parallel one (for the different number of threads) and next calculated the speedup of the parallel program.

Below, we present the results of experiments carried out with the codes shown in Listing 1 (serial code) and Listing 4 (parallel code). As we mentioned in the previous section, we cannot find any related parallel code that fulfills exactly the same computations as those executed with the code in Listing 1. Thus, we limited our experiments to the codes mentioned above.

To carry out experiments, we used a processor Intel Xeon X5570, 2.93 GHz, 2 physical units, 8 (2 × 4) cores, 16 hyper-threads, and an 8 MB cache. Executable parallel tiled code was generated by means of the g++ compiler with the -O3 flag of optimization.

Experiments were carried out for ten different lengths of the problem defined with parameter *N* from 1000 to 5000 for the codes presented in Listing 1 (serial code) and Listing 4 (parallel code).

All of the source code to perform the experiments and the program to run the tested codes can be found at https://github.com/piotrbla/livc (accessed on 24 August 2021).

We carried our experiments to choose the optimal size of a tile. The size of a tile is optimal if (i) all data associated with that tile can be held in cache, (ii) those data occupy almost the entire capacity of cache, and (iii) tiled code execution time is minimal provided that the conditions (i) and (ii) above are satisfied.

For this purpose, we fulfilled three trials whose results are shown in Figure 1. The curve "trial1" represents how the time of tiled program execution depends on the block size along axis *h*0 when the block size along axis *w*0 is fixed equal to 16. After first trial 1, we chose the best size along the *h*0 axis equal to 32.

**Figure 1.** Time for different tile sizes. All phases.

The curve "trial2" demonstrates how the time of tiled program execution depends on the block size along axis *w*0 when the block size along axis *h*0 is equal to 32 (the result of trial 1). After trial 2, we chose the best size along axis *w*0 equal to 24.

The curve "trial3" shows how the time of tiled program execution depends on the block size along axis *h*0 when the block size along axis *w*0 is equal to 24 (the result of trial 2). After trial 3, we chose the best tile size along axis *h*0 equal to 54. Finaly, as the best size of a 2D tile in the parallel tiled code, we chose 24 × 54.

For the best tile size, Table 1 presents execution times and speedup of the serial program in Listing 1 and the parallel tiled one presented in Listing 4 for 32 OpenMP threads used. Figure 2 depicts the data presented in Table 1 in a graphical way. As presented, the execution time of parallel tiled program grows practically in a linear manner exposing considerable speedup (the ratio of the serial program execution time to that of the corresponding parallel one) presented in Figure 3.

**Table 1.** Time in seconds and speedup for Intel Xeon X5570 and 32 OpenMP threads.

**Figure 2.** Time for Intel Xeon X5570 and different problem sizes.

**Figure 3.** Speedup for Intel Xeon X5570 v3 and 32 OpenMP threads.

Figure 4 presents how parallel tiled code speedup depends on the thread number for the maximal problem size used for experiments, i.e., for *N* = 5000. The parallel tiled code speedup grows practically linear for the number of threads 1 to 12. Linear speedup for the number of threads ≥ 12 is prevented with the serial fraction of code (Amdahl's law) parallel loop initialization fulfilled with a single thread and serial input–output operations. Speedup is also limited with parallel program overhead—there is thread synchronization in the examined parallel code, after each wavefront, barrier synchronization is inserted because the following wavefront can be executed after completing the calculations of the previous wavefront.

**Figure 4.** Speedup for Intel Xeon X5570 for different threads number.

We may conclude that the generated parallel tiled code implementing computingintensive general linear recurrence equations and presented in Listing 4 can be successfully run on modern multicore machines with a large number of cores.

#### **6. Conclusions**

We presented an approach to generate parallel tiled code for computing general linear recurrence equations (GLREs) presented in Listing 1. That code is computing-intensive and must be run on modern multicore computers to reduce execution time. We demonstrated how to transform that code to obtain the modified code shown in Listing 3, which exposes dependencies such that there exist two linear independent solutions to the time partition constraints formed on the basis of those dependencies. This allows us to apply the affine transformation framework and generate parallel 2D tiled code computing GLREs presented in Listing 4. The parallelism is achieved using the wavefront technique and presented with the code that conforms to the OpenMP standard. To our best knowledge, the target parallel tiled code generated by us and presented in Listing 4 is the first to allow for enumerating 2D tiles and the first that does not require any additional computations in comparison with those of the original serial code. This code is derived by means of tiling the loop nest iteration space. Our experiments with the resulting parallel tiled code show that the code significantly outperforms the original GLREs computing serial code. The code performance improvement is achieved due to the parallelism and better locality of the target code.

**Author Contributions:** Conceptualization and methodology, W.B. and P.B.; software, P.B.; validation, W.B. and P.B.; data curation, P.B.; original draft preparation, W.B.; writing—review and editing, W.B. and P.B.; visualization, P.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Source code to reproduce all the results described in this paper can be found at: https://github.com/piotrbla/livc (accessed on 24 August 2021).

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Abbreviations**

The following abbreviations are used in this manuscript:

GLRE General Linear Recurrence Equations

TRACO Compiler based on the TRAnsitive ClOsure of dependence graphs

#### **References**

## *Article* **Design of Automatic Correction System for UAV's Smoke Trajectory Angle Based on KNN Algorithm**

**Pao-Yuan Chao \*, Wei-Chih Hsu and Wei-You Chen**

Department of Computer and Communication Engineering, National Kaohsiung University of Science and Technology (NKUST), Kaohsiung 807618, Taiwan

**\*** Correspondence<EMAIL_ADDRESS>
**Abstract:** Unmanned aerial vehicles (UAVs) have evolved with the progress of science and technology in recent years. They combine high-tech, such as information and communications technology, mechanical power, remote control, and electric power storage. In the past, drones could be flown only via remote control, and the mounted cameras captured images from the air. Now, UAVs integrate new technologies such as 5G, AI, and IoT in Taiwan. They have a great application value in a high-altitude data acquisition, entertainment performances (such as night light shows and UAV shows with smoke), agriculture, and 3D modeling. UAVs are susceptible to the natural wind when spraying smoke into the air, which leads to a smoke track offset. This study developed an autocorrect system for UAV smoke tracing. An AI model was used to calculate smoke tube angle corrections so that smoke tube angles could be immediately corrected when smoke is sprayed. This led to smoke tracks being consistent with flight tracks.

**Keywords:** unmanned aerial vehicle; machine learning; UAV smoke show; mobile networks; artificial intelligence

#### **1. Introduction**

Flexible, safe, stable, high-speed, and low-cost UAVs have been developed in the past few years. This was achieved due to the continuous development of the modules, including the process materials, electric power storage, and sensors [1–12]. So far, UAVs have been widely used in civilian, commercial, and government units. Industrial UAVs can be applied in environmental monitoring, infrastructure inspection, disaster or accident rescue, agriculture, forestry, fishery, animal husbandry management, spatial information measurement, land and guard patrol, media communication, telecommunications services, home delivery logistics, and the military. Most application fields can be further subdivided. For example, environmental monitoring can be divided into the monitoring and investigating of air pollution, oil pollution, nuclear pollution, marine pollution, and river pollution, and even includes the study of weather changes. Many types of infrastructure are subjects for inspection, including roads, railways, transmission towers, and oil fields. Regarding the rescue, drones can be used for video recording, a real-time image transmission, and material delivery. Regarding the environmental conditions, there are waters, mountainous areas, or buildings. Agriculture, forestry, fishery, and animal husbandry management includes pesticide or fertilizer spraying and the observation of crops, trees, pastures, and fish farms. The work in spatial information includes aerial mapping, a terrain attribute classification and survey, a national land survey, urban planning, a land survey and development, water control and flood control planning, and 3D real scene modeling. Guard patrol includes a coastal patrol, criminal chasing, and general security work. Regarding media communication, in addition to providing real-time news about disaster areas and war zones, they can be applied to business and tourism marketing. There are diverse applications in the military, and they can be used as reconnaissance aircraft, target aircrafts, and bombers. Therefore, industrial UAVs have unlimited business opportunities. In Ghana,

**Citation:** Chao, P.-Y.; Hsu, W.-C.; Chen, W.-Y. Design of Automatic Correction System for UAV's Smoke Trajectory Angle Based on KNN Algorithm. *Electronics* **2022**, *11*, 3587. https://doi.org/10.3390/ electronics11213587

Academic Editor: Juan M. Corchado

Received: 23 September 2022 Accepted: 31 October 2022 Published: 3 November 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

there have been about 275,000 UAVs flying and delivering medical kits containing vaccines [13]. The edge computing technology is used in UAVs for the power transmission line inspection [14].

On many major holidays and celebrations in Taiwan, people can see the colorful smoke from fighters flying by the military in the sky. In addition to the purchase of fighters, such activities often cost a lot of money (military aircraft maintenance, personnel training, and aviation gasoline), which consumes gasoline and imposes a burden on the environment. Therefore, the use of UAVs carrying smoke tubes for aerial smoke spraying has grown in recent years. Today, most UAVs are mainly powered by electricity and do not emit exhaust gas like gasoline engines. Moreover, UAVs cost less than traditional fighters in smoke spraying. For example, in Taiwan, according to the Regulations of Drone of the Civil Aeronautics Administration [15], if pilots hold G2 licenses (flying more than 400 feet above the ground or water, operating beyond the range of visibility, dropping, or spraying objects) and UAVs are registered in the UAV system of the Civil Aeronautics Administration and apply for airspace in advance, regulatory restrictions can be exempted and UAV shows with smoke can be performed.

UAVs cannot carry large smoke tubes and smoke tubes installed at different positions are susceptible to the natural wind and lead to a smoke track offset. In this study, detectors and smoke tube correctors were installed above the UAV. An AI model was used for training to correct the offset tracks. In this way, small low-altitude UAVs could achieve the same visual effects as traditional high-altitude fighters in UAV shows with smoke, in a way which is cheap and environmentally friendly.

#### **2. Materials and Methods**

#### *2.1. Unmanned Aerial Vehicle*

The UAV used in this study has the following basic flying parts: a flight controller, motor, electronic transmission, frame, and GPS positioning module. In addition, a servo motor, smoke tube, Raspberry Pi, 4G communication module, and lithium battery are installed. The overall weight is about 2.5 kg. Considering the motor load in the flight process and flight time, a UAV with a wheelbase of 450 mm and an EDU450 carbon fiber frame was selected [16], as shown in Figure 1.

**Figure 1.** EDU450 carbon fiber frame.

A Pixhawk2 CUBE consisting of two processors was used as the flight controller of the UAV. In the Pixhawk2 CUBE, the main processor was STM32F427 V3, and the coprocessor was STM32F1. The built-in sensor included a tri-axis accelerometer (L3GD20), an accelerometer and magnetometer (LS303D), a gyroscope (MPU9250), and a barometer (MS5611). With a weight of 73 g, this light and efficient flight controller supported opensource flight control software PX4 and Ardupilot. Three–eight axis multi-rotor models and multiple interfaces, including the Mavlink interface, I2C interface, and PWM signal output system, are used in this study. Raspberry Pi 3B+ has a USB interface which can be used for Internet access apart from the built-in WiFi card. If the 4G network card fails, it can switch to the built-in WiFi [17], a backup network, as shown in Figure 2.

**Figure 2.** System configuration.

#### *2.2. Smoke Tube Position*

In this study, a four-axis multi-rotor UAV was used for the experimentation. The smoke from the smoke tube is vulnerable to the natural wind and tube position, leading to a smoke track offset. In that case, the audience is unable to enjoy the performance. As shown in Figure 3, when the smoke tube was placed directly below or above the UAV, the smoke from the smoke tube was affected by the downdraft generated by the propeller. Regardless of the rotation of the servo motor and adjustment of the smoke spraying direction, the smoke would be affected by the airflow and sprayed downward. In this study, self-made 3D material parts were attached to the UAV. An additional extension area was built to install the servo motor and smoke tube for adjusting the spraying angle. This is to prevent the influence of downdraft generated by the propeller as much as possible and to make the smoke from the smoke tube be sprayed backward, as shown in Figure 4a. To avoid the pendulum effect and the excessive output energy of the rear motor, the length of the extension area is adjusted to 29 cm after multiple outdoor flight tests. The aerial testing shows that the smoke track is significantly improved, as shown in Figure 4b. After improving the smoke emission direction, the track angle correction was discussed below, so that the audience could enjoy the best effect of the smoke track.

**Figure 3.** Spraying effects at different smoke tube positions, should be listed as: (**a**) below the UAV; (**b**) above the UAV.

**Figure 4.** Spraying effects at the redesigned smoke tube position, should be listed as: (**a**) the smoke tube extension frame is at the rear; (**b**) effect after modification.

#### *2.3. AI Model Selection and Design*

K-nearest neighbors (KNN) are one of the most popular machine learning algorithms. It has been widely used in HPC applications, such as image/video retrieval, big data analysis, machine learning, and computer vision [18,19]. It is a nonparametric statistical method for regression and classification. The K-nearest training samples in the feature space ware input [20,21] and the k-value were used to determine which classification group the data were nearest to. The classification criteria were decided by majority voting, and the Euclidean distance was used to calculate the distance, as shown in Equation (1).

$$\mathbf{P} = \sqrt{\left(\mathbf{x}\_2 - \mathbf{x}\_1\right)^2 + \left(\mathbf{y}\_2 - \mathbf{y}\_1\right)^2} \tag{1}$$

In the KNN classifier, the output was the group classification, and its neighbors determined the category corresponding to the input object by a majority voting. KNN adopted the vector space modal for a classification, and objects in the same category were highly similar. The similarity could be calculated by the known category cases to evaluate the possible categories of the input objects. The training samples were multi-dimensional eigenspace vectors, in which each training sample had a classification label. The algorithm included eigenvector access and training sample labels at the training stage.

The KNN classifier assigns a weight of 1/k to the k-nearest neighbors and zero to all other neighbors. This can be applied to the weighted nearest neighbor classifier. The weight wni is given to the nearest neighbor i. In (2), a similar result holds for the strong consistency of the weighted nearest neighbor classifier [21].

$$\sum\_{i=i}^{n} \mathbf{w}\_{\text{ni}} = 1 \tag{2}$$

Let Cwnn n denote the weighted nearest classifier with the weight {wni}n i=1. According to the regularity condition of the category distribution, the excess risk has (3) an asymptotic expansion.

$$\mathbf{R}\_{\rm R}(\mathbf{C}\_{\rm n}^{\rm wmn}) - \mathbf{R}\_{\rm R}\left(\mathbf{C}^{\rm Bayes}\right) = \left(\mathbf{B}\_{1}\mathbf{s}\_{\rm n}^{2} + \mathbf{B}\_{2}\mathbf{t}\_{\rm n}^{2}\right)\left\{1 + o(1)\right\}\tag{3}$$

$$\mathbf{s}\_{\mathbf{n}}^{2} = \sum\_{i=1}^{n} \mathbf{w}\_{\mathbf{n}i}^{2} \tag{4}$$

$$\mathbf{t}\_{\mathbf{n}} = \mathbf{n}^{-\frac{2}{d}} \sum\_{i=1}^{n} \mathbf{w}\_{\mathbf{n}i} \left\{ \mathbf{i}^{1 + \frac{2}{d}} - (\mathbf{i} - \mathbf{1})^{1 + \frac{2}{d}} \right\} \tag{5}$$

The optimal weighting method { w∗ ni }n i=1 is used to balance the two items above. Let k∗ = Bn 4 d+4 , (6) correspond to i = 1, 2, ... , k∗, and w∗ ni = 0 correspond to i=k∗ + 1, . . . , n. After using the optimum weight, the dominant term in the asymptotic expansion of excess risk is O n− 4 d+4 .

$$\mathbf{w}\_{\rm ri}^{\*} = \frac{1}{\mathbf{k}^{\*}} \left[ 1 + \frac{\mathbf{d}}{2} - \frac{\mathbf{d}}{2\mathbf{k}^{\*2/\mathbf{d}}} \left\{ \mathbf{i}^{1+2/\mathbf{d}} - (\mathbf{i}-1)^{1+2/\mathbf{d}} \right\} \right] \tag{6}$$

At the classification stage, k is a user-defined constant. A vector without a category label (query or test point) will be classified into the most frequently used category among the k sample points nearest to the point.

When flying in the air, the UAV is affected by the crosswind wwind and deviates from its route. At this point, to correct the route, the flight controller will give a Roll value to adjust the pitch angle θ1 of the UAV, as shown in Figure 5. When the airframe is corrected, an angle adjustment θ2 is given to the smoke tube to make the smoke tube turn to the windward face to face the direction the wind comes from, as shown in Figure 6. According to the angle adjustment θ1 of the flight controller and the angle θ2 of the smoke tube, the direction and magnitude of the wind the UAV is exposed to in the air can be known, as shown in Equation (7).

$$\mathbf{w}\_{\text{wind}} \rightarrow \; \theta\_1 + \; \theta\_2 \tag{7}$$

**Figure 5.** Roll angle correction θ1.

**Figure 6.** Smoke tube angle correction *θ2*.

Based on the above conclusion, the direction and magnitude of the wind are related to the value of θ1. The operator can adjust θ2 according to the value of θ1 when flying the UAV. θ1 and θ2 will be trained by the machine learning-based KNN classification method, and their relationship is shown in Equation (8).

$$\mathbf{w} \propto \theta\_1 \to \theta\_2 \tag{8}$$

#### **3. Experimental Results and System Validation**

During its flight, the UAV is affected by the natural wind and deviated from its course. At this point, the flight controller corrects the pitch angle in real-time to make the UAV return to its course. Its flight direction was mainly changed by correcting the pitch, yaw, and roll parameters. In addition to the airspeed meter sensor installed above the UAV, the three sensing values of the pitch, yaw, and roll parameters of the UAV can be used to learn the changes in the wind fields in the air.

Five angles, namely −60◦, −30◦, 0◦, 30◦, and 60◦ were designed in this study. These angles are the output y to be estimated by the KNN model. The input features include the pitch, yaw, roll, and airspeed meter values read from the flight controller. A KNN model was built for training.

Finally, the trained model was stored in the Joblib package and then put into the Raspberry PI in the UAV. Later, the designed system was used to read the model so that the real-time wind speed data read could be directly put into the AI model in the Raspberry PI for calculation. The results could be transmitted to the flight controller via the system to adjust the servo motor that controlled the smoke tube direction. It helped adjust the smoke tube to the optimal angle. The correction flow chart is shown in Figure 7.

**Figure 7.** Autocorrection flow chart for smoke trailing.

As the UAV which was selected could not be fitted with a larger smoke tube, each spraying took about 30 to 40 s. The operator could collect about five pieces of data on each flight, which is not much. During the training, the pitch, yaw, roll, and airspeed meter data were put into the KNN model, and the accuracy was 50%. To improve the accuracy, more data is required. Therefore, the pitch and yaw were discarded, and the roll value was kept. The UAV was designed to fly back and forth in a straight line automatically. In the case of the deviation caused by a crosswind, the UAV could return to its route mainly by correcting the roll value. A roll value of 0 indicates that the UAV flew horizontally without any roll. A positive roll value indicates that the wind blew from the left side of the UAV. A higher value represents a higher wind speed. On the contrary, a negative roll indicates that the wind blew from the right side of the UAV. A smaller value reflects a higher wind speed. In this study, the roll value was collected to determine the speed and direction of the wind. This is to make up for the limitation of the airspeed meter under the breeze. The roll value was sensitive and thus could detect detailed data. In this study, 64 data items from the database were put into the AI training model, and the accuracy was 71%, as shown in Figure 8.

**Figure 8.** KNN training accuracy.

The roll data read by the UAV was recorded and put into the KNN model for testing. The test data and predicted angles are shown in Table 1. The table shows that the rolls and angle corrections were as expected. The smoke tube should have been shifted to the right when the roll was positive and left when the roll was negative.

**Table 1.** KNN model test data and results.

Images taken behind the smoke tube show that the smoke tube could automatically change the direction and angle according to the direction and speed of the wind. Figure 9 shows that the UAV deviated to the left due to the wind from the right side. The roll of the angle correction to the right given by the flight controller was 0.061363. The smoke tube should be adjusted 30◦ to the right according to the calculation by the AI model. Figure 10 shows that the UAV deviated to the right due to the wind from the left side. The smoke tube should be adjusted 30◦ to the left, according to the calculation by the AI model.

**Figure 9.** The smoke tube shifted to the right.

**Figure 10.** The smoke tube shifted to the left.

The smoke from the smoke tube would be adjusted to the direction of the windward face. Based on the wind speed, the smoke tube would be adjusted to 30◦ or 60◦. In this case, when the windward face of the UAV was in the front and rear directions, the smoke tube angle would not be adjusted.

Based on the observation of the actual flight, an accuracy of 71% indicated a significant improvement. Figure 11 illustrates the case without correction by the AI model. From the audience's angle, it could be clearly seen that the track (yellow arrow) of the smoke from the smoke tube was offset due to the natural wind. Figure 12 shows the case with a correction by the AI model. From the audience's angle, the smoke from the smoke tube could be observed (yellow arrow). Despite the influence of the wind field in the air, the smoke track could be manipulated to be almost the same as the flight route.

**Figure 11.** Smoke track correction without the AI model.

**Figure 12.** Smoke track correction with the AI model.

#### **4. Conclusions and Future Work**

In this study, Raspberry PI was used as the microcomputer for transmissions with the server. The flight data were received to control the smoke tube and sensor of a quadaxis UAV. To avoid the influence of the airflow, 3D-printed parts were used to refit the UAV. Its frame was extended to install a smoke tube and an electronic igniter so that the smoke tube could be lit for shows with smoke. As the smoke from the smoke tube was susceptible to the natural wind and became offset, a servo motor was installed to adjust the direction of the smoke from the smoke tube. A manned aircraft flew to record the angle adjustment, wind direction, and wind speed. The KNN was used to train a modified AI model. After applying the AI model to the Raspberry PI, the UAV emitted smoke in the air. The Raspberry PI and the flight controller could directly read the wind field data, and the angle could be immediately calculated and then sent back to the flight controller if it needed to be corrected. In this way, the spraying angle of the smoke tube could be adjusted immediately to make the smoke track the same as the flight route as much as possible. According to the results, the correction accuracy was 71%, which can demonstrate the difference between before and after the correction.

In the past, fighters sprayed smoke in the air to celebrate major festivals, which was expensive and polluted the environment. Our design is expected to make UAV shows with smoke possible on small occasions so that such events can be enjoyed on many occasions other than major festivals. The UAV designed by us is powered by electricity. Compared with a fuel-powered aircraft, it causes less environmental pollution and is cheaper.

The architecture in this study was designed for single UAVs. If the information of multiple UAVs can be displayed simultaneously on the web page and the crowd control can be carried out through function buttons on the web page, multiple shows with smoke can be performed simultaneously to spray. As for the collection of the wind speed data, this study collected various data, such as the roll, pitch, yaw, and anemometer values. Due to insufficient sample data, only the roll data were used for the machine learning. If more data can be collected in the future and all data collected can be imported for machine learning, the AI model will be more effective, and the overall correction effect will be perfect.

**Author Contributions:** Conceptualization, P.-Y.C. and W.-Y.C.; methodology, W.-C.H. and W.-Y.C.; resources, P.-Y.C. and W.-C.H.; writing—original draft preparation, P.-Y.C. and W.-Y.C.; writing review and editing, W.-C.H. and W.-Y.C.; visualization, W.-Y.C.; supervision, P.-Y.C. and W.-C.H. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

## *Article* **Prediction of Offshore Wave at East Coast of Malaysia—A Comparative Study**

**Mohammad Azad 1,\* and Md. Alhaz Uddin 2**

**Abstract:** Exploration of oil and gas in the offshore regions is increasing due to global energy demand. The weather in offshore areas is truly unpredictable due to the sparsity and unreliability of metocean data. Offshore structures may be affected by critical marine environments (severe storms, cyclones, etc.) during oil and gas exploration. In the interest of public safety, fast decisions must be made about whether to proceed or cancel oil and gas exploration, based on offshore wave estimates and anticipated wind speed provided by the Meteorological Department. In this paper, using the metocean data, the offshore wave height and period are predicted from the wind speed by three state-of-the-art machine learning algorithms (Artificial Neural Network, Support Vector Machine, and Random Forest). Such data has been acquired from satellite altimetry and calibrated and corrected by Fugro OCEANOR. The performance of the considered algorithms is compared by various metrics such as mean squared error, root mean squared error, mean absolute error, and coefficient of determination. The experimental results show that the Random Forest algorithm performs best for the prediction of wave period and the Artificial Neural Network algorithm performs best for the prediction of wave height.

**Keywords:** prediction; artificial neural network; support vector machine; random forest; regression; offshore wave; wind speed

#### **1. Introduction**

Human activities in the offshore region are continually increasing due to oil and gas exploration. Adverse conditions can often occur due to environmental disasters during offshore operations. Predicting wave characteristics is a crucial prerequisite for offshore oil and gas development (Figure 1), considering the safety of lives and the avoidance of economic damage. Wave height and period are typically significantly increased by the wind associated with storms passing across the ocean's surface. Weather forecasting departments usually predict wind forces rather than wave periods and heights. There are several empirical approaches for estimating wave height and wave period from wind force [1–5]. Estimating wave height and period is inherently inaccurate and random, making it difficult to simulate using deterministic equations [6]. The numerical approach for calculating wave height and period is a difficult and complex procedure that, despite substantial breakthroughs in computational tools, produces solutions that are neither dependable nor consistently applicable. Machine learning methods are perfect for modeling inputs with corresponding outputs since they do not necessitate an understanding of the underlying physical mechanism [7].

Several studies have been performed using artificial neural networks (ANN) to measure important wave heights and mean-zero-up-crossing wave period history for different locations in the seas. These parameters were predicted three, six, twelve, and twenty-four hours in advance using two different neural network methods [8–10]. The time series of

**Citation:** Azad, M.; Uddin, M.A. Prediction of Offshore Wave at East Coast of Malaysia—A Comparative Study. *Electronics* **2022**, *11*, 2527. https://doi.org/10.3390/electronics 11162527

Academic Editor: Giuseppe Ciaburro

Received: 25 April 2022 Accepted: 23 July 2022 Published: 12 August 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

these wave parameters has been investigated using simulations in Portugal's western coast area. Time series with wave height have been disintegrated into multi-resolution time series using wavelet transformation hybridized with ANN and wavelet transformation [11,12]. As the input of the ANN, the multi-resolution time series has been used to predict the important wave height at an unlike multi-step lead time near Mangalore, India's west coast. Mandal et al. [13] expected wave heights from observed ocean waves off the west coast of India, Marmugao. To predict wave height, recurrent neural networks with a resilient propagation (Rprop) update algorithm have been implemented.

**Figure 1.** Offshore oil and gas development

The time series of wave height and mean-zero-up crossing wave period were simulated using data from the local wind. The application of various types of machine learning models [14–17] was performed to improve the accuracy of the prediction. The authors [18] used ANN to estimate wave height from wind force at 10 chosen locations in the Baltic Sea. The WAM4 wave model was used to figure out the time series for the waves that had been forecasted in the past. There were two different techniques, Feed-forward Back-propagation (FFBP) and Radial Basis Function (RBF), used to develop a machine learning system for predicting wave heights at a particular coastal point in deeper offshore areas [19–21]. Data on wave height, average wave period, and wind speed were obtained from remotely sensed satellite data on India's west coast. Tsai et al. [22] employed the ANN with a back-propagation algorithm to estimate wave height and cycle from wind input based on the wind-wave relationship [23]. Time records of waves at either station may be predicted based on data from the neighboring station. Several deterministic neural network models have been developed by Deo et al. [24] to predict wave height and wave periods from generated wind speed. However, the model can offer adequate results in deep water and open areas, and the prediction periods are extensive.

This research uses wind force/metocean data to correctly predict wave height and wave period. Critical comparisons are performed for the current investigation to provide more accurate findings in predicting the wave parameters using three machine learning algorithms: ANN (Artificial Neural Network), SVM (Support Vector Machine), and RF (Random Forest).

The main contributions of the paper are (i) performing the experiments on the data sets obtained on the east coast of Malaysia and (ii) understanding the best predictive models for future prediction.

We arrange the remaining parts of the manuscript as follows. First, we explain the data collection procedures, the three regressor models, and performance metrics used for the comparison among the models in Section 2. Then, we show the experimental results and discuss the findings in Section 3. Finally, Section 4 contains a concise conclusion.

#### **2. Materials and Methods**

In this section, first, the data collection procedure is described, followed by the regression methods, and finally the evaluation metrics are described.

#### *2.1. Data Collection*

On a 2° × 2° grid in South China, environmental data are collected along Malaysia's east coast, including the basins of Sabah (longitude 114.39° E, latitude 5.83° N) and Sarawak (longitude 111.82° E, latitude 5.15° N) (Figure 2). The Malaysian Meteorological Service and Fugro OCEANOR provided the data that was acquired from satellite altimetry [25,26] using oceanographic SEAWATCH meteorological (metocean) buoys and sensors, calibrated and corrected by Fugro OCEANOR [27].

**Figure 2.** Data collection location of east coast of Malaysia [26]. Reprinted from Renewable Energy, 88, Omar Yaakob,Farah Ellyza Hashim,Kamaludin Mohd Omar,Ami Hassan Md Din,Kho King Koh, Satellite-based wave data and wave energy resource assessment for South China Sea, 359-371, 2016, with permission from Elsevier.

The summary of common statistical values of the collected data is given in Table 1.

**Table 1.** Summary of basic statistics of the collected data.

There are 1460 data samples and three features: wind speed (m/s), wave height (m), and wave period (s).

#### *2.2. Methods*

There are a lot of algorithms [28–33] that can be used for regression analysis. In this study, three well-known and commonly used regression algorithms (ANN, SVM, and RF) were chosen to predict wave height and period from wind force. The wind force is employed as an input for training the network, while the wave height and period are used as outputs. Before discussing the details of each method, Table 2 shows the advantages and disadvantages of each method [34]:

The ANN algorithm is the most frequently utilized algorithm in such problems. The difference between the state-of-the-art ANN algorithm and two commonly used algorithms, SVM and RF, is then determined. Note that ANN is a parametric regression algorithm while SVM and RF are nonparametric regression algorithms. However, because the data set is small enough, it is not possible to employ advanced deep learning techniques (such as convolutional neural networks or recurrent neural networks).

#### 2.2.1. Artificial Neural Network (ANN)

Inspiring by biological neurons, the researchers created artificial neurons [35]. It is possible to create a network of artificial neurons to predict the desired functions (i.e., the target variables).

The basic idea behind a single artificial neuron is that it takes an input function *I*, which is a product of weight vectors with the sample input vector, and feeds it into an activation function *g* to produce an output (Figure 3). Generally, the usual practice is to create a network of neurons based on many layers: input layers, hidden layers, and output layers. The input layers basically consist of one or more neurons based on the supplied inputs (for the present study, only one input is used), then there can be one or more hidden layers (only one hidden layer is used), and finally, there can be one output layer consisting of the expected outputs (Figure 4). The general trend is to use a back-propagation algorithm to update the weights involved in each connection between neurons so that it is possible to minimize the errors of regression between true output and the predicted output.

**Figure 3.** The architecture of a single artificial neuron.

**Figure 4.** The architecture of ANN.

The hyperbolic tangent function, tanh(*x*) = *ex*−*e*−*x ex*+*e*−*x* is used as the activation function. The steps of the Python implementation are described in Algorithm 1.

#### **Algorithm 1** The Python implementation of the ANN algorithm

**Input:** Data set containing features and target

**Output:** Prediction of ANN for the given data set

7. Predict using the obtained regressor and calculate the performance metrics on the test data set;

#### 2.2.2. Support Vector Machine (SVM)

The basic idea behind regression using SVM (i.e., also popularly known as Support Vector Regression) is that given input examples {(*i*1, *o*1), ... ,(*in*, *on*)} ⊂ *X* × *R*, where *X* denotes the example space (e.g., *X* = *R* for our problem), a function *f*(*x*) is obtained that has at most deviation from the actual outputs oi for all the input examples (Figure 5) [36].

Generally, it is easy to describe the case of linear functions [37] *f* , *f*(*x*) =< *w*, *x* > +*b* , with *w* ∈ *X*, *b* ∈ *R*

where < , > denotes the dot product in *X*. Furthermore, it is possible to rewrite this problem as a convex optimization problem:

minimize 1 2 ||*w*2|| subject to

$$y\_i - < w\_\prime \\ x > -b \\ \le \epsilon$$

$$ +b \\ -y\_i \le \epsilon$$

The steps of the Python implementation are described in Algorithm 2.

**Input:** Data set containing features and target

**Output:** Prediction of SVR for the given data set

1. Read the data set using the Python read\_excel function;

2. Extract features and target variables;

3. Scale the features and target using the StandardScaler function;

6. Derive the SVR using the above results;

7. Predict using the obtained regressor and calculate the performance metrics on the test data set;

#### 2.2.3. Random Forest (RF)

A decision tree has been widely used from the very beginning of machine learning research to find the best hypothesis as the classifier and regressor. A decision tree is a tree-like structure where a sequence of actions is taken from the root to the leaf nodes based on the values of each decision node in the concerned path. A sample decision tree is depicted in Figure 6, where the course of action of playing outside is taken based on the value of whether it is training outside or not.

**Figure 6.** A simple decision tree model.

There are a number of different variants of decision tree algorithms, from ID3 [38], C4.5 [39], and CART [40] to ensembles of decision trees (e.g., Random Forest [28])) that have been proposed to improve the accuracy of the classifiers and regressors. For a single decision tree construction, the most widely used splitting criteria, "Gini index" or "Entropy", can be mathematically stated as follows [32,38–41]:

$$p\_t = \frac{N\_t(T)}{N(T)}$$

$$\text{• } \quad \text{Entropy } \text{ent}(T) = -\sum\_{t \in D(\underline{T})} p\_t \log\_2(p\_t);$$

$$\text{• } \quad \text{Gini index } \operatorname{gini}(T) = 1 - \sum\_{t \in D(T)} p\_t^2.$$

where *T* is the data set and *D*(*T*) is the set of labels in the data set *T*. In addition, *N*(*T*) represents the number of samples and *Nt*(*T*) represents the number of samples with the label *t*.

A Random Forest (RF) is an ensemble of decision tree models such that each tree in the ensemble is built based on the bootstrap samples of the training data set. Furthermore, during the construction of the decision trees, the random forest selects only a subset of the features at each split point. In this way, the constructed decision trees are more different from each other to reduce the correlation among them and to have a better prediction. In contrast to the single decision tree (CART), random forests do not use any pruning of the tree.

Like any ensemble algorithm, the random forest also takes the average among all decision tree regressors to produce the final predictions (Figure 7).

**Figure 7.** The average prediction of *n* decision trees in a random forest prediction model

The steps of the Python implementation are described in Algorithm 3.

#### **Algorithm 3** The Python implementation of the RF algorithm

#### **Input:** Data set containing features and target

**Output:** Prediction of RF for the given data set

7. Predict using the obtained regressor and calculate the performance metrics on the test data set;

#### *2.3. Performance Metrics*

Four performance metrics [42] have been used for the comparison of the results.

For each sample, the error (*ei*) is the difference between actual output (*oi*) and predicted output (*o*ˆ*i*). The average of the squares of the error (*ei*) is the *MSE*. The average is taken by dividing the summation of the square of all the errors by the total number of samples (*n*).

$$MSE = \frac{\sum\_{i=1}^{n} (o\_i - \delta\_i)^2}{n} \tag{1}$$

Root mean squared error (*RMSE*) is the square root of the average difference in squared error between actual output (*oi*) and predicted output (*o*ˆ*i*). It is the most commonly used metric.

$$RMSE = \sqrt{\frac{\sum\_{i=1}^{n} (o\_i - \delta\_i)^2}{n}} \tag{2}$$

Mean absolute error (*MAE*) is the average absolute difference of error between actual output (*oi*) and predicted output (*o*ˆ*i*). It is the most commonly used metric. The advantage of *MAE* is that it is softer to outliers. The reason is that it does not have any square term associated with the equation, therefore the penalty for the outlier points is not that much compared to *MSE* and *RMSE* where there are heavy penalties imposed by the square term.

$$MAE = \frac{\sum\_{i=1}^{n} |(o\_i - \delta\_i)|}{n} \tag{3}$$

The coefficient of determination or *R*2 explains the degree to which the independent variables explain the variation of the output variable. It is a measure of how well new samples can be predicted by the model through the proportion of explained variance. If *oi* is the true value of the *i*-th sample and *o*ˆ*i* is the corresponding output or predicted value, the estimated *R*2 can be calculated as:

$$\mathcal{R}^2 = 1 - \frac{\sum\_{i=1}^n (o\_i - \delta\_i)^2}{\sum\_{i=1}^n (o\_i - \delta)^2} \tag{4}$$

where *o*¯ = ∑*n i*=1 *oi n* .

#### **3. Results and Discussion**

The experiments have been executed based on the available field data from the east coast of Malaysia, which consists of 1460 samples, one input variable, wind speed (m/s), and one target variable, either wave height (m) or wave period (s). The regression results are compared among the three regressors (ANN, SVM, and RF). The Python programming environment (version 3.6) with scikit-learn (version 0.24.1) is used for the implementation.

Initially, the ANN regressor is trained using the training data and then the accuracy is validated using the testing data. The network is chosen with one hidden layer of 4 neurons, and each neuron uses the hyperbolic tangent function. The weight of the ANN is optimized using "lbfgs" in the family of quasi-Newton methods. Furthermore, the SVM regressor is trained using the same training data and validated using the same testing data as mentioned above. For SVM, the standard radial basis function (RBF) is used as the kernel function. Finally, the RF regressor is trained and validated in a similar fashion. For RF, 400 trees are used for the number of trees, 3 is used for the maximum depth of the tree, and 1460 × 0.05 = 73 samples are used for the training of each decision tree in the ensembles. Each regressor is compared by the state-of-the-art performance metrics.

#### *3.1. Cross Validation Results*

It is a common practice to keep a part of the available data for testing and the remaining part for training. However, the model evaluation is particularly dependent on specific pairs of (training and testing) fractions, and the result can be overfitting. To overcome this problem, cross validation [43] is used to test the model's fitness to predict new data that was not used to train the model and how it can generalize to unknown data.

There are many variants of cross-validation. The standard method is to use *k*-fold cross-validation. In general, the procedure is to use *k* rounds in cross-validation. It splits the data into *k* subsets. In a single round, it completes the analysis on one subset (the training set), and then validates the performance on the remaining *k* − 1 subsets (the validation set). In the next round, another subset is chosen for training and the remaining is chosen for validation. To reduce variability, these steps are repeated *k* times for *k* subsets so that each subset is selected for training. Finally, an average among the *k* validation results is reported as the fitness of the model's predictive performance.

Nevertheless, it is common practice to repeat the k-fold cross-validation process multiple times and report the average performance among all folds and all repeats. This approach is called repeated *k*-fold cross-validation.

The average and standard deviation of the regression performance of wave height were determined using 10-fold cross validation that was done for three methods as presented in Table 3. The four-performance metrics are reported in the same way as in the preceding section.

**Table 3.** Regression results for 10-fold cross validation with three times repeated for wave height (m).

It is clearly evident that the method of ANN has a minimum average value of mean absolute error, mean square error, and root mean square error. Besides, it has the maximum average value of the coefficient of determination. As a result, this approach performs best for predicting wave height.

Similar findings for wave period have been presented in Table 4 using 10-fold cross validation that has been done for the consecutive three methods.

**Table 4.** Regression results for 10-fold cross validation with three times repeated for wave period (s).

It is apparent that the RF approach has the lowest mean absolute error, mean square error, and root mean square error. Furthermore, it has the maximum value of the coefficient of correlation. As a result, this approach is the most accurate for predicting wave periods.

*R*2 is close to 1, which indicates that regression predictions perfectly fit the data. However, in our case, the experimental results show the lower value of *R*2 due to the inconsistencies in the collected data set.

#### *3.2. Graphical Representation of a Sample Training and Testing Results*

To illustrate the model's performance graphically, one sample is chosen randomly as the pair of (training, testing), where the training is 80% and testing is 20% of the data. In Table 5, the results for both training and testing data are shown in the case of the wave height regression problem. It is evident that the mean squared error, mean absolute error, and root mean squared errors are the smallest for RF compared to others in the training data. However, these metrics are the smallest for ANN compared to others in the testing data. Furthermore, the coefficient of correlation (*R*2) is largest for RF compared to others in the training data, but it is largest for ANN compared to others in the testing data. As a result, the ANN approach is superior at predicting wave height in the future.

**Table 5.** A sample regression results for wave height (m).

The regression results of ANN, SVM, and RF for wave height are depicted in Figure 8. Figure 8a shows the ANN regression findings for wave height. The black dots are actual testing data, and the green dots are the predicted values. It is clear that the predicted values follow the pattern of the actual testing data. Nevertheless, the regression results for the SVM regressor are shown in Figure 8b. The black dots are actual testing data, and the blue

dots are the predicted values. It is clear that the predicted values follow the pattern of the actual testing data.

**Figure 8.** *Cont.*

**Figure 8.** The regression result for wave height (m); (**a**) shows results using ANN, (**b**) shows results using SVR, and (**c**) shows results using RF. The predicted values follow the pattern of the actual testing data in all cases. For RF, the pattern is not as smooth as in ANN or SVM, because RF is not a single decision tree method; rather it is an ensemble method that works by taking the average votes. ANN prediction is smoother than SVR.

The regression results of wave height based on the RF regressor are presented in Figure 8c. The black dots are actual testing data, and the orange dots are the predicted values. It is clear that the predicted values follow the pattern of the actual testing data. The pattern is not as smooth as SVM or ANN because RF is not a single decision tree method; rather it is an ensemble method that works by taking the average of different decision trees' predictions.

In addition, in the case of wave period, Table 6 shows the findings for both training and testing data. It is evident that the mean squared error, mean absolute error, and root mean squared error is the smallest for RF compared to others in both training and testing data. Nevertheless, the coefficient of determination (*R*2) is the largest for RF compared to others in both training and testing data. Therefore, for future prediction of wave period, the RF method is best.

**Table 6.** A sample regression results for wave period (s).

The regression results of ANN, SVM, and RF for wave period are depicted in Figure 9. For ANN, in Figure 9a, the actual testing data is shown in black dots and the predicted values are shown in green dots. It is clear that the predicted values follow the pattern of the actual testing data for ANN. For SVM, in Figure 9b, the actual testing data are represented by black dots, while the actual predicted values are represented by blue dots. It is clear that the predicted values follow the pattern of the actual testing data for SVM.

**Figure 9.** The regression result for wave period (s); (**a**) shows results using ANN, (**b**) shows results using SVR, and (**c**) shows results using RF. The pattern is not directed positively as like wave height in Figure 8 and rather going horizontally which indicates lower values of *R*2 compared to wave height.

Finally, for RF, the real testing data is represented in black dots, whereas the actual predicted values are shown in orange dots in Figure 9c. It is clear that the predicted values follow the pattern of the actual testing data for RF. The pattern is not as smooth as SVM or ANN because RF is not a single decision tree method; rather it is an ensemble method that works by taking the average of different decision trees' predictions.

#### *3.3. Comparison with Standard Non-Parametric Kernel Regression*

The non-parametric method of kernel regression (KR) in statistics is used to calculate the conditional expectation of a random variable. The goal is to find a non-linear relationship between two random variables, *I* and *O* [44]. The problem under consideration can be modeled using a standard non-parametric regression problem. As an example, *I* can be wind speed and *O* can be wave height. In Table 7, the previous results are compared with the results of KR for the same sample.

It is clear that KR results are not far from those of the SVM, ANN, or RF. In fact, ANN produces the best results in the case of wave height prediction.

#### *3.4. Overall Discussion*

The goal of this study is to understand the best predictive model among ANN, SVM, and RF for the prediction of wave height and period from the wind speed. A detailed experiments are performed in terms of cross-validation and sample training and testing results. The summary of this result is given in the Table 8.

It is evident from Table 8 that the Random Forest (RF) method is truly performing well across two prediction problems. It is the second best for the prediction of wave height and the best for predicting wave period. Nevertheless, it has the advantage of being a non-parametric method. Moreover, the underlying tree structure has the advantages of interpretability and usage. However, to be specific, RF performs best for the prediction of wave period while ANN performs best for the prediction of wave height.

#### **4. Conclusions**

This study carried out three different and well-known machine learning algorithms for the prediction of offshore waves. These approaches are used to predict the wave height and wave period from the given wind forces. Multiple accuracy ranges are obtained in terms of the mean absolute error, mean square error, root mean square error, and coefficient of determination. Overall, these performance measures show average behavior. However, it is possible to compare the three employed methods and analyze the results. The regression analysis in the random forest performs best for the prediction of wave period, and the artificial neural network performs best for the prediction of wave height. Furthermore, it was compared with the standard non-parametric kernel regression and found to have a similar result.

In situations when a traditional analysis would be challenging, these machine learning techniques can produce very quick and reasonable predictions. These studies can benefit the community as a measurement of safety and precaution in the critical marine environment.

In this regard, there are numerous potential future research directions. One disadvantage of using such metocean data is the existence of discrepancies. Future work should incorporate advanced algorithms in addition to the aforementioned models to address such inconsistencies. Additionally, more machine learning techniques should be used to find the ideal answer for this specific prediction problem. The time dependence of the metocean data, which can be examined using time series analysis tools, is another interesting subject of research.

**Author Contributions:** Conceptualization, all authors; methodology, all authors; software, M.A.; validation, M.A.; formal analysis, all authors; investigation, all authors; resources, all authors; data curation, M.A.U.; writing, all authors; visualization, all authors; supervision, all authors.; project administration, M.A.U.; funding acquisition, M.A. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** It is available with a justifiable and reasonable request.

**Acknowledgments:** The authors are grateful to Jouf University for their assistance with this research. Moreover, the authors would like to express their gratitude to all of the volunteers and anonymous reviewers for their suggestions.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

### *Article* **Classification of Task Types in Software Development Projects**

**Włodzimierz Wysocki 1,\*, Ireneusz Miciuła 2 and Marcin Mastalerz 3**

**Abstract:** Managing software development processes is still a serious challenge and offers the possibility of introducing improvements that will reduce the resources needed to successfully complete projects. The article presents the original concept of classification of types of project tasks, which will allow for more beneficial use of the collected data in management support systems in the IT industry. The currently used agile management methods—described in the article—and the fact that changes during the course of projects are inevitable, were the inspiration for creating sets of tasks that occur in software development. Thanks to statistics for generating tasks and aggregating results in an iterative and incremental way, the analysis is more accurate and allows planning the further course of work in the project, selecting the optimal number of employees in task teams, and identifying bottlenecks that may decide on faster completion of the project with success. The use of data from actual software projects in the IT industry made it possible to classify the types of tasks and the necessary values for further work planning, depending on the nature of the planned software development project.

**Keywords:** software; knowledge management; reasoning; information extraction; rule mining; knowledge acquisition; engineering

#### **1. Introduction**

The contemporary intensity of changes in the surrounding reality due to technological progress makes management by economic units extremely demanding and difficult. Therefore, it is necessary to react quickly and effectively to changes that create new conditions for business activity. In the era of knowledge-based economy, information, information systems co-created by it, and related information technologies are extremely valuable and inextricably linked with knowledge [1]. The usefulness of IT systems is enormous and directly influences the increase of the possibilities of management units by reflecting their innovativeness and technological potential, which is of great importance when making strategic business decisions.

Globalization, the era of e-economy, and computerization are the topicality of modern economic activity [2]. That is why the production of an IT product (IT product) is so fundamental. Managing the production of IT products is a new scientific and technological discipline that has emerged at the interface between computer science and management engineering. All economic entities that base their activities on IT solutions are interested in the practical results of research in this discipline. In addition, it should be stated that, indirectly from all kinds of improvements in software development and management, any activity that uses any IT software will benefit. The digitization of the economy results in the need for reliability and security of emerging applications, IT systems, or transaction services. A lack of reliability of IT software can cause enormous losses. That is why it is so important to develop and create IT software that is as reliable as possible. Therefore,

**Citation:** Wysocki, W.; Miciuła, I.; Mastalerz, M. Classification of Task Types in Software Development Projects. *Electronics* **2022**, *11*, 3827. https://doi.org/10.3390/ electronics11223827

Academic Editor: Maria Liz Crespo

Received: 3 November 2022 Accepted: 19 November 2022 Published: 21 November 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

the work efficiency of teams producing IT software—characterized by a cyclical nature—is so important. Moreover, the literature highlights the dependence of the implementation effectiveness of subsequent software development stages on current project stages and the maturity level of the IT product.

The software development process is constantly evolving; there are technological changes and ways of organizing the work of project teams. In the knowledge economy era, most economic activities require systems and all kinds of IT products [3]. This forces the appropriate acceleration of software development, while at the same time maintaining the required quality. That is why it is so important to have the right team to implement often complex and innovative IT solutions. On the other hand, numerous project teams cause problems in effective communication in the group and cause the need to optimize work organization. To ensure this optimization, tools and IT products designed for managing software development processes are often used.

Traditional approaches to IT project management have become insufficient due to the high variability of requirements and the need to change the work organization of project teams [4]. It is particularly noticeable in innovative software development projects for the e-economy, where flexibility in the production process and originality of solutions are required. Therefore, despite the support of software development processes with increasingly better tools and systems, the pursuit of optimal use of project resources including human teams—is still a considerable challenge. Since the beginning of the 21st century, the agile approach has been increasingly used in the practice of software development. The literature on the subject shows a significant improvement in the number of successful projects and the optimization of the use of resources necessary in software development through the use of continuously developed agile methods. This is particularly important due to the constant changes taking place in the requirements for the software being developed [5]. Therefore, traditional methodologies cannot be successfully applied to the numerous changes that are natural in innovation projects because they are extremely static and, apart from the first phase, consist of a fixed number of tasks. The next stage of the methodology development is based on iterations, which assumes that we learn about new requirements in the development process; requirements change and the existing ones are detailed. However, in the current agile methodologies, the requirements are incomplete and we do not know the details; this is also due to the nature of the innovative projects [6]. We have a set of requirements in the form of epics and user stories. Then, in the course of the software development process, creating new requirements becomes the cause of new tasks. Another source of new implementation tasks is the tests of current versions that detect defects or reveal new implementation possibilities. Automation of this process by the task generator in agile methodologies allows for the creation of appropriate increments of tasks; the number and the nature of which depends on the size and type of the project. On the other hand, the iterative measurement of tasks allows for subsequent iteration planning, with the assumption of the invariability of certain aspects for the best estimation. Therefore, planning methods based on summary work are nowadays of great practical importance for the management of manufacturing processes.

Project management systems for software development collect a lot of data about tasks and information related to them necessary to perform the work and the necessary exchange of information between the project team [7]. However, there is still room for improvement when it comes to providing information that will optimize IT project management. Because the data and information collected in the systems are semantically poor, it is necessary to analyze the works that are performed during the implementation of tasks to correctly classify them, which will also allow the estimation of their size, thus increasing the quality of planning. The isolation and classification of the characteristic types of tasks, and the work performed within them, will allow for the automatic execution of analyses and estimations to create insightful statistics and reports necessary for optimal forecasting decisions. In addition, taking into account the knowledge about the current technological and business components will give an appropriate insight into the current state of work in the project and will allow for more optimal planning of the software development process.

This work aims to identify the impact of the type of tasks on efficiency and to indicate those factors that positively affect the effectiveness of software development teams. The model of task types was built based on the data analysis of real software projects from the financial sector, managed by the Jira system. Connecting the model with the Jira system enables easy data acquisition for analysis and increases its commercial potential. A separate abstract layer of the model, in combination with a dedicated database, supports the possibility of creating interfaces to other IT project management systems. The article attempts to classify project tasks, thus determining the necessary expenditure on tasks of various types. In conjunction with the results of research on tasks created during the production process and cyclical works, it allows for planning projects. Linking task types with project team roles allows the simulation of project work and supports project management by identifying bottlenecks in the manufacturing process and avoiding over-employment.

The article discusses the basic concepts of the context in which software development projects are implemented and the development of methodology for the optimal management of these processes. At the same time, the principles and practices, as well as the lifecycle of software development projects, are characterized, which is important for the possibility of introducing improvements in this process. The concept of the programming process model is inspired by the agile approach to managing software as it is developed. Extending the original approach to the roles of team members and task types allows more precise work planning in the project and management of the project team, including detecting bottlenecks or unused resources. The article discusses the programming process model that is oriented toward the manual recognition of task types, thus enabling effective support for planning tasks carried out in innovative IT projects.

#### **2. Literature Review**

Software development in IT companies is mainly carried out through appropriate organization of the work of project teams. Unfortunately, in projects involving the search and creation of new solutions where the variability and unusual nature of ideas and requirements are natural, there is a constant need to improve the management of this process. Due to the innovative nature of IT projects, project teams are burdened with high risk. This forces the search for optimal project management methods and techniques that will allow for better control and use of the company's resources. Nowadays, analytical methods and tools embedded in IT systems are a great help in supporting project management.

Traditional software development project management methodologies are currently being replaced by, or supplemented with, agile methodologies. The development of project management methods in the IT industry has resulted from dissatisfaction with the small number of successful projects. A big problem was discovered in the management of innovative projects, where rigid and strongly formalized traditional software development methodologies did not work and even made it difficult to introduce changes and the proper functioning of projects [8]. Agile methods are used primarily due to the desire to reduce the number of errors, to shorten the time needed to create finished products, or to reduce the total production costs [9]. However, as noted by J. Shore and S. Warden, when making decisions regarding the implementation of agile practices, it is difficult to unequivocally state greater successes in the implementation of IT projects [10]. This is due to the variety of interpretations of the concept of "IT project success" and the fact that the use of the same methodologies in different companies from the same industry results in different results [11]. There is a common view in the literature that the key to success in project management is learning about appropriate techniques and tools [12]. However, the instrumental layer of an appropriate project management methodology alone will not ensure its success if not properly applied by its members, because the most important part of software development projects are the people [13]. Therefore, an appropriate balance should be found between the hard (budget, schedule, and implementation time) and soft (communication, changes, motivation, and competencies) elements of the project. Therefore, regardless of the methodology used, project management is closely related to people management, and practice shows that most of the problems affecting project success result from the omission of the "purely human aspects" of team management [14]. Accordingly, this article aims to identify the main problems related to people management in the success of an IT project by enabling effective task planning carried out in innovative, and thus variable, software development processes.

Since the beginning of the 21st century, we have witnessed a dynamic development in agile methodologies [15]. This is due to the adaptation of the management process to projects with a high degree of innovation. In these cases, hard methodologies such as PRINCE2 and PMI PMBoK do not work due to the detailed and long-term planning stage that is not able to take into account future changes [16,17]. With these characteristics of the projects, the too-high level of standardization of activities does not work, because there are often unsuccessful attempts at establishing the project lifecycle and detailed requirements already at the project initiation stage. Even though this approach gives a lot of comfort to the implementers, in the case of innovative projects, it does not meet the real needs that will be known only at further stages of implementation. The progress of work on the innovative project reveals the necessity to repeatedly verify many assumptions, actions, and plans, because, only as the implementation progresses, knowledge and actual visualization of the developed solutions in the shape of the final product are acquired. Additionally, certain paths for reaching specific results turn out to be ineffective only after the verification and testing phase. External changes should also be mentioned, i.e., changes that are beyond the control of project teams, e.g., changes in regulations and legal norms or the current market situation. Often, consistent plan implementation leads to the creation of functionalities that will not be adapted to reality or that will be useful only after costly modifications. Therefore, this article attempts to classify the characteristic types of tasks in software development projects that will allow for more effective planning and adaptation of the required work to dynamically changing reality.

All factors that make it difficult, or even impossible, to precisely define and describe the results of the project at the stage of its definition result from the following elements that occur when creating software (especially innovative software) [18,19]:

In addition, the implementation of innovative IT projects is associated with the risk of other irregularities. Kruchten (2004) lists several main problems that can affect any project, which are [20]:

The above factors and elements of project implementation management, and the inability to solve them in accordance with the traditional approach, contributed directly to the development of methods for making traditional methodologies more flexible and, as a result, prepared the ground for the Manifesto for Agile Software Development that announced in 2001 the principles for agile software development methodologies [21,22]:

The principles of an adaptive approach to software development project management indicate the direction for project teams, while specific practice is necessary for the actual implementation of works [23]. The process structure and specific practices create a minimal flexible framework for self-organizing teams. IT tools are essential to accelerate software development and to reduce costs. Contracts are crucial for the development of the customer– supplier relationship. Documentation supports communication [24]. However, the key issue is to provide the project team with feedback to answer the question of where the team currently is in the software development process [25]. This is possible thanks to iteration and incremental software lifecycles. The essence of the iterative process is the frequent delivery of working pieces of software (successive increments) that implement selected sets of functions that together make up the usefulness of the final product. The iterative software development cycle leads to a management style where long-term plans are fluid, while a stable plan can be created for a short period of time. Iterative and incremental software development leads to completely new relationships with the business client and different principles of the project team's functioning.

The concept of the iterative and incremental programming development cycle as a remedy for dilemmas of the e-economy era has resulted in the creation of new methodologies for IT project management [26]. These methodologies, called agile, do not cut off completely from document-oriented formalized traditional methodologies, but have specific features (adapted to the requirements of modern software development projects) [27,28]:

Organizations are increasingly implementing digital transformation plans, due to the threat of disruptions, in order to keep pace with the growing pace of business. Agile software development plays a huge role in this process. Many of the digital workflows in use today are based on agile principles. Thanks to a flexible scalable IT infrastructure, cloud computing is evolving in line with the needs of agile software development. The DevOps concept removes the traditional distinction between software development and operations. The software is used as a tool in SRE–DevOps implementation and systems management and automation of operational duties [29]. CI/CD methodologies confirm that software will change frequently and provide tools that help developers deliver new code faster [30]. Agile methodologies come in a variety of forms to meet the needs of any project [31]. Even though agile approaches are different, all of them are based on the key ideas contained in the agile approach. For this reason, any framework or behavior that complies with these principles is termed agile. Regardless of the specific agile approaches that a team decides to apply, the benefits of an agile methodology can only be fully realized through the collaboration of all parties involved [32]. In recent years, a large number of agile software development project management methodologies have emerged. The most popular are [33–36]:

time delivery can never be compromised and that design modifications are always to be expected. This belief is based on eight principles and a business-based methodology.

• Lean development, a "lean" software development. The basic idea of the approach is the elimination of losses understood as elements that do not add any value to the product. The aim of such action is to deliver the finished product to the customer as soon as possible. Lean software development is based on the values and principles of adaptive project management. The term "lean software development" is considered by Mary Poppendieck and Tom Poppendieck, who, in their book *Lean Software Development: An Agile Toolkit,* presented, among others, the seven main principles of lean management and a set of 22 techniques supporting the approach.

Agile management methodologies are a group of methodologies that are characterized by an adaptive and variable approach to managing software development projects. Additionally, agile methodologies were developed much earlier than the agile manifesto itself. The first works on adaptive project management methods date back to the 1980s (an example is the rapid application development methodology) and the concept of agile methodologies was introduced in the mid-1990s [37]. The idea of a time frame is a well-defined process dedicated to software development control at the lowest level in an iterative cycle with several review points. Reviews help ensure the quality and efficiency of software development. By delivering the software on time at the lowest level, the timely production at the highest level (i.e., the project level) is ensured [38]. The basic principle of the project plan is to prepare a schedule of planned increments and, within them, the planned time frames, which will create a complete project schedule capable of changes with emerging new requirements. The use of the time frame technique, together with the MoSCoW prioritization technique, ensures no delays in project implementation and the delivery of ready-made software that will meet business goals within a given time [39].

Projects with a high degree of innovation are very difficult to include in a complete schedule and scope of work. Therefore, adaptive methodologies describe functionalities (i.e., independent elements of the subsystem), which in subsequent releases can be quickly changed and handed over for implementation. Agile methodologies, as opposed to traditional methodologies, rule out the validity of long-term planning [40,41]. Therefore, the plans are speculative and not deterministic. This allows you to adapt to all types of changes that appear during the implementation of the project. Additionally, the distinguishing factor of agile methodology is a strong emphasis on the cooperation and integration of the project team, because only in this case is the smooth flow of information and effective communication ensured. The general scheme of the project lifecycle in the case of agile methodologies is based on five phases, indicated by J. Highsmith [18,42]:

Agile project management methodologies require an appropriate level of project maturity for organizations and project teams [43]. In addition, agile methods continue to be improved in order to most effectively respond to the needs of managing still-innovative software development projects. Hence, there are so many methods that are still looking for new, more perfect, solutions; although, they undoubtedly already seem to be better adapted than the classic methods to dynamically changing project environments. Implementation based on iteration allows for effectively adapting signals that come from both inside the project and from the external environment.

Classifying specific types of tasks in software development projects will allow for determination of the requirements, the time necessary, and the expenditure that will be needed to implement them. This will allow for more effective work on further adaptive project planning methods. Linking task types with the roles of the project team will allow you to simulate project work, which will support project management by identifying bottlenecks in the software development process. In addition, it will allow for avoiding over employment and will allow support for quick "what if" simulations. At the same time, the introduction of the task classification in terms of business and technological components, in conjunction with the employee competency model, will allow for automatically selecting the composition of the project team and optimally managing the team; this has the highest priority in agile methodologies. This is of great importance for the optimization of software development process management and leads to the development of perfect methods in response to the dynamism of real-world changes. It will have replicative and predictive capabilities to plan the project, to simulate it during changes, and to detect bottlenecks during the entire process. In this article, experiments were carried out on many real IT projects to classify task types and to attempt to find the relationship between the nature of the planned project and the specificity and number of these tasks.

#### **3. Materials and Methods**

The model of task types was built on the basis of the data analysis of real software projects from the financial sector, managed by the Jira system. Table 1 contains summary data of these projects. In four projects, the team worked in accordance with traditional methodologies, to which more and more elements of the agile approach were introduced. In two projects, the teams worked in accordance with the scrum approach. It is quite a large data set that allows for the analysis of phenomena related to modern manufacturing processes. The total number of project issues exceeds 30,000, which makes it possible to use artificial intelligence methods. The total duration of the projects is approximately 22 years, which does not mean that historical records are dealt with. Project teams worked on most projects in parallel.

**Table 1.** Projects, the data of which were used for research on task types.

In Jira, teams use issues to describe and track specific tasks to be done. Issues are the basic building blocks of projects. The issue can be of a specific type. The most common type of issue is task, which in practice is used for many purposes. The second most frequent type is bug, which describes the defects found in the developed software and the work needed for fixing them. The full list of types used in the projects with their frequency is shown in Figure 1.

**Figure 1.** The occurrence frequencies of issue types in the researched projects.

The most commonly used issue types are task, bug, and subtask. The summary determines the idea of the actions to perform. The description field contains a detailed description of the work to do. A significant feature of the issue is its state. Tracking issue state changes allows for the recognition of the performed work. The workflow describes the values of the state field and the allowed transitions between them. Figure 2 shows a typical simplified workflow used in the analyzed projects.

**Figure 2.** A typical simplified issue workflow used in projects.

As shown in Figure 2, transitions between issue states correspond to performing basic actions in implementing new system functions or fixing bugs. For example, the DevOps (development and operations) methodology helps to establish cooperation between developers and operators to automate the continuous delivery of new software, which is expected to contribute to shortening the development cycle and to creating high-quality software [44,45]. Another development of DevOps is the concept of development, security, and operations (DevSecOps), which at the same time is designed to integrate security methods with the software development process, where security measures are built in to ensure the integrity and availability of the application [46].

A valuable feature of the Jira system is storing the history of changes in the value of issue elements, including the state. Storing history supports the tracking of issue execution. Not all issues use the state to track work. The Jira system allows employees to register working hours devoted to work on an issue. For some types of work, there is no need to keep track of their state, as they are repetitive works performed as needed. In this case, the possibility of registering working hours is sufficient.

The presented approach is based on a precise division of the development process into activities and roles in line with the RUP methodology, adapted to hybrid and agile processes [47,48]. The previous series of articles described the agent–object model of the manufacturing process (AGOMO), which was first used to assess the maturity of RUP processes, then to plan hybrid water–scrum–fall processes [49,50]; it became an inspiration for the research presented here.

As we showed in the previous section, the types of Jira issues are insufficient to clearly define the work's purpose and type. The model was created in order to fill the gap that prevents linking the work carried out with the composition and competencies of the project team. The combination of these two areas will allow for a more precise quantitative analysis of the projects' works and the detection of the causes of the observed phenomena. It is a way to optimize the efficiency of development processes and to increase the level of maturity of project teams and organizations [51].

Software projects consist of tasks that a project team performs to build an IT system [52,53]. The task represents the Jira issue. Each of the tasks performed has a clearly defined goal. This goal can be, for example, implementing a requirement, testing a component or the entire module, administering environments, or managing a project; it determines its type. Members of the project team perform the tasks. Roles define the competencies and responsibilities of those carrying out the tasks. One person can have many roles; one role can have many people.

The task of implementing system functions requires performing some essential subtasks. They include implementation, i.e., the creation of component source code, code review, creation of unit tests, testing and verification of functions, and acceptance tests. Such tasks correspond to the issue, the state of which changes according to the workflow shown in Figure 2. The implementation tasks that consist of subtasks are named stateful tasks in the model.

Each task and subtask contains the number of project team members' working hours. Task and subtask types enable association efforts with the roles of project team members. They help to recognize bottlenecks or over-employment in completed IT projects and to avoid them during project planning.

Issues whose states do not change during the project are represented in the model by recurring tasks. They have been called recurring because a single task of this type can be performed as needed for the project's duration. For example, this might be a project development management task performed by an administrator when defects arise or when new software components need to be configured. For recurring tasks, you can calculate the average labor intensity in a period and, on this basis, assume what the fixed costs of the project are [54]. Genuinely recurring tasks are, for example, meetings such as daily stand-up meetings, workshops with clients, or steering group meetings.

The authors analyzed the tasks of the studied projects in terms of their type. The results, in the form of graphs showing the percentage number and the sizes of stateful and recurring tasks, are shown in Figure 3. The chart on the left shows the ratio of the number of stateful tasks to the recurring tasks in the projects. The chart on the right shows the ratio of the effort of stateful tasks to recurring tasks.

**Figure 3.** The number and effort of stateful and recurring tasks in the researched projects.

The case of the P2 project is significant, where the number of stateful tasks is greater than the recurring tasks (amounting to over 90%), while the effort of these tasks is slightly over 40%. It is very interesting, because stateful tasks are responsible for implementing the functions of the built system and for fixing the detected defects in it. From the developer's point of view, this workload should probably be the greatest. From the point of view of a researcher of software development processes, this phenomenon arouses great curiosity. In order to satisfy it, it is first necessary to precisely define the types of tasks in the projects, which this article implements on the basis of actual and implemented projects from practical activity.

#### **4. Results**

The classification of task types is based on the division of tasks into stateful and recurring tasks. By definition, stateful tasks represent work related to the implementation of system components and the fixing of defects. Stateful tasks are processed by a subtasking algorithm. Recurring tasks are responsible for the remaining works. Recurring tasks can be classified on the basis of summary and extended-text descriptions included in the issue. Initially, these tasks were classified manually. Subsequently, simple classification algorithms were created based on keyword searches. Recurring tasks are also broken down into subtasks that correspond to the registered work.

Figure 4 shows a diagram of the classification algorithm activity. The algorithm first checks how many times the issue has changed its state; if at least three times (more than open and close), the workflow is the basis for dividing the task into subtasks. The split algorithm tries to create subtasks (e.g., development, testing, code review). For this, it uses the value of the state before and after, the time of the change, the role of the person, and the registered works and then assigns completed work to the created subtask. The result is a stateful task.

If the task did not change state or if the algorithm did not detect any subtasks, the task type is determined by searching for keywords in the text description. When it fails, the tasks go to a spreadsheet, where they are manually classified. Then, the algorithm checks who was working on the task. If many people worked on a task on one day, it is a recurring group task (e.g., stand-up and other meetings). If one person worked on a task on one day, it is a recurring individual. The way of dividing the recurring task into subtasks depends on the individual/group classification. If the programmer and the tester alternately execute a task, it is stateful; the algorithm divides them into subtasks according to the roles of the team members.

**Figure 4.** Algorithm of task classification.

The algorithm for breaking stateful tasks into subtasks is very complex due to the many ways that users use Jira to work on a project issue. Searching for keywords is not very accurate and supporting it by manual classification makes it impossible to use the classification in practice. Therefore, it is planned to replace these solutions with the NLP (natural language processing) classification. However, there are some good points to manual ranking. During the work on the algorithm, the process of the manual analysis of recurring tasks detected more than 50 types of tasks that were aggregated into three groups containing 14 main types of tasks.

Stateful tasks consist of tasks for implementing new features and for fixing bugs reported by customers and testers. Recurring tasks are divided into three main groups: implementation, meetings, and organizational tasks. A complete list of the main task types is provided in Table 2. The task types listed in the table form a hierarchical structure divided into types and groups.

The task classification algorithm—classifying tasks performed alternately by programmers and testers as stateful tasks—increased the number and effort of state tasks. An updated version of the graphs in Figure 3 is included in Figure 5. The changes are significant. For example, for project P3, the percentage of stateful tasks increased from 58% to 66% and the percentage of stateful task efforts increased from 24% to 43%.

The development of task types and an algorithm enabling automatic classification of the researched projects allowed for a detailed analysis of the tasks of the researched projects. Below, we present effort charts (Figure 6) for six groups of task types in the researched projects.

The P1 project was implemented in a hybrid methodology. It is the longest project among the researched, with a duration of more than 8 years. The project went through many phases of implementation, delivery, and maintenance. Therefore, the values and the ratio of effort in the project were averaged. They can serve as a benchmark when compared with other hybrid and traditional designs.

The P2 project was carried out in the hybrid methodology. It is characterized by a large number of hours used for various types of meetings. On the other hand, the very small scope of work in the management field suggests that the project may have been managed collectively. The very little work involved in fixing defects found by clients indicates that time spent in meetings was well spent.

**Table 2.** Types of stateful and recurring tasks used by the model.

**Figure 5.** Updated numbers and effort of stateful and recurring tasks.

The P3 project was carried out in a hybrid methodology. It has high administrative and management costs. The ratio of repairing defects found by customers to repairing defects found by testers is interesting. There is no such tendency in the other projects, except perhaps for the P4 project. This may indicate inaccurately defined requirements or difficult contact with the customer.

**Figure 6.** Effort of the six main groups of task types in projects.