1 Introduction

Optimization over a simplex appears naturally in many applications: machine learning, signal processing, control, and game theory [1]. In such problems, the feasible set is determined by linear constraints that form a simplex.

Branch and Bound (B &B) algorithms aim at finding the set of global minimum points up to a certain accuracy. These types of algorithms have the advantage of providing an exhaustive search and the requested precision of the solution. Their computational burden increases with the dimension of the problem and the accuracy. Therefore, they are usually applied to low-dimensional problems. The overview book [11] discusses why partitioning a box-constrained area into simplices may have an advantage over using box-shaped partition sets. Their focus is mainly on Lipschitz-based bounds to guarantee the required accuracy. In this paper, we consider the feasible set as a simplicial set, that is a way to avoid dealing with linear constraints.

The use of longest edge bisection in simplicial branch and bound has a long tradition and has been studied thoroughly, see [3, 11, 13]. Recently, new ways of bounding the objective over a simplex have been studied [8,9,10, 15], where bounds on the first and second derivatives are considered. Bounding rules using gradient information were found to outperform other bounding rules, allowing monotonicity to be taken into account.

Monotonicity considerations to remove subsets and to reduce dimension have a long tradition in Interval Arithmetic (IA) based B &B, see [2, 12]. The basic property is to be able to remove interior boxes where the function is monotone over the box and to reduce the dimension with respect to monotone components, when a box facet is at the boundary of the search space.

The ideas of monotonicity have been extended for simplicial partition sets recently [5, 8, 9, 14]. One of our questions is how to handle the mathematical properties in an efficient way in a simplicial branch and bound. The main question is how to extend the ideas of a simplicial dimension reduction. This requires identifying whether a corresponding facet is at the boundary of the feasible set. The latter will be defined more precisely with the terminology of border faces.

The main contributions of this study are:

  • A mathematical theory to support new monotonicity tests over a simplex that can be used independently on how reliable bounds of the objective and its derivatives were obtained. Here, we use Interval Arithmetic because it is a simple way to calculate the bounds.

  • A method to determine when a facet of a simplex can be considered border of the search space.

  • A computational algorithm.

  • Computational illustration of the new monotonicity tests on artificially generated problems.

To investigate the research questions, Sect. 2 introduces the notation and mathematical properties that can be used in a simplicial branch and bound. Section 3 describes several ways to discover directions in which the objective function is monotonic. Section 4 presents a method for determining whether a facet of a simplex belongs to the border of the feasible set. Section 5 focuses on the algorithm under investigation and various ways to implement it. Section 6 reports on the behaviour of the algorithm for various instances. Our findings are reported in Sect. 7.

2 Notation and Mathematical Properties

2.1 Notation

We consider minimization of a continuously differentiable function \(f:\mathbb {R}^n\rightarrow \mathbb {R}\), over a feasible set \(\Delta \), which is an \(n-\)simplex. An \(n-\)simplex can be defined by its vertices, \(v_0,\ldots ,v_n \in \mathbb {R}^n\), a set of \(n+1\) affine independent vectors. Namely, \(\Delta :={{\,\textrm{conv}\,}}(\mathcal {W})\) with set of vertex \(\mathcal {W}:=\{v_0,\ldots ,v_{n}\}\). The considered problem is

$$\begin{aligned} \min f(x),x\in \Delta . \end{aligned}$$
(1)

Actually, the problem can be posed more generically, where the number of vertices of \(\Delta \) is smaller, and thus the feasible set is not full-dimensional compared to the objective function [16]. The latter is relevant in mixture models, where quantities, probabilities or prices are normalized in such a way that the decision variables are on the standard simplex. We study a branch and bound algorithm to enclose all minimum points of f on \(\Delta \). In contrast to the algorithms described in [11], the used partition sets are \(m-\)simplices S, where \(m\le n\) can be not full-dimensional \(S:={{\,\textrm{conv}\,}}(\mathcal {V})\) with \(\mathcal {V}\) a set of \(m+1=|\mathcal {V}|\) vertices. The branch and bound algorithm works with a set \(\Lambda \) of partition sets, which as a whole includes all global minimum points.

The studied algorithm uses the longest edge bisection (LEB), where the longest edge (vw) of a partition set S is bisected using mid-point \(x:=\frac{v+w}{2}\) leading to two new simplices with vertex sets \(\mathcal {V}\setminus \{v\} \cup \{x\}\) and \(\mathcal {V}\setminus \{w\} \cup \{x\}\). This is sketched in Fig. 1. Two Longest Edge bisections (LEB) are performed iteratively over the feasible area \(\Delta \) leading to new points \(x_3\) and \(x_4\) and a partition of 3 simplicial partition sets.

Fig. 1
figure 1

Three partition sets generated by two longest edge bisections, using bisection points \(x_3\) and \(x_4\). We focus on monotonic directions in S

The algorithm also considers dimension reduction due to monotonicity considerations where the set \(\mathcal {V}\) of vertices of \(m-\)simplex S is reduced and S is replaced by one (or more) of its facets \(F:={{\,\textrm{conv}\,}}(\mathcal {V}\setminus \{v\})\) for some \(v \in \mathcal {V}\). Notice that F is an \((m-1)-\)simplex. It may be clear that for \(m=0\), the \(0-\)simplex \(S={{\,\textrm{conv}\,}}(\{v\})\) is an individual point and has no faces. Its dimension cannot be reduced.

The centroid of \(m-\)simplex \(S= {{\,\textrm{conv}\,}}(\{v_0,v_1,\ldots ,v_m\})\) is given by \(c:=\frac{1}{m+1}\sum _{j=0}^m v_j\) and the relative interior is defined by

$$\begin{aligned} {{\,\textrm{rint}\,}}(S)=\{x=\sum _j\lambda _jv_j, \lambda _j>0, j=0,\ldots ,m,\sum _{j=0}^m \lambda _j=1\}. \end{aligned}$$
(2)

The relative boundary of a simplex S is defined by removing the relative interior from it. Considering a simplicial partition set S, we are interested in studying whether its (simplicial) facets F are border with respect to the feasible set \(\Delta \).

Definition 1

Given the feasible area as an \(n-\)simplex \(\Delta \). An \(m-\)simplex S with \(m<n\) is called border with respect to \(\Delta \) if there exists an \(m-\)simplex face \(\varphi \) of \(\Delta \), such that \(S\subseteq \varphi \).

One of the main questions is how to determine whether a facet of a simplex is border in a numerically efficient way, see Sect. 4. Border facets and the concept of monotonicity are used to reject a simplex or to reduce its simplicial dimension.

To find monotonic directions, we will use the information of a so-called enclosure G of the gradient \(\nabla f(x)\subseteq G:=[\underline{G},\overline{G}], \forall x \in S\) where \(\underline{G}\) and \(\overline{G}\) denote the lower and upper bound of interval box G. The interval vector G can be calculated by automatic interval differentiation over the so-called interval hull of a simplex S. The interval hull, represented as \(\Box S\), is the smallest box which encloses S [4, 17]. Now consider a directional vector d as the difference of two points \(x,y \in S\), \(d=x-y\). The corresponding directional derivative \(d^T\nabla f(x),\ x\in S\), is included in the inner product

$$\begin{aligned} d^TG=\left[ \underline{d^TG},\overline{d^TG}\right] = \left[ \sum _{i=1}^{n} \min \{d_i\underline{G}_i, d_i \overline{G}_i\}, \sum _{i=1}^{n} \max \{d_i\underline{G}_i, d_i \overline{G}_i\}\right] . \end{aligned}$$
(3)

2.2 Mathematical Properties on Monotonicity

First, we discuss properties related to full-dimensional simplices, where \(m=n\) and then generalize the properties to the general case, where \(m\le n\). For the interior of a full-dimensional \(n-\)simplex, the monotonicity properties are the same as for any full-dimensional convex set:

Proposition 1

Let \(S\subset {{\,\textrm{rint}\,}}(\Delta )\) be an \(n-\)simplex with gradient enclosure \(G {\ \supseteq }\)\({\nabla f(x), \forall x \in S}\). If \(\exists \, i\in \{1,\ldots , n\}\) with \(0\notin G_i\), then S does not contain a global minimum point of (1).

For the proof, see Proposition 2 in [8].

Proposition 2

Let \(S\subseteq \Delta \) be an \(n-\)simplex with gradient enclosure \(G {\ \supseteq \nabla f(x),}\)\({ \forall x \in S}\). If \(\exists \, i\in \{1,\ldots , n\}\) with \(0\notin G_i\) then \({{\,\textrm{rint}\,}}(S)\) does not contain a global minimum point of (1).

For the proof, see Proposition 2 in [8].

Corollary 1

Let \(S\subseteq \Delta \) be an \(n-\)simplex with corresponding gradient enclosure \(G {\, (\supseteq \nabla f(x), \forall x \in S)}\) as partition set in a branch and bound algorithm. If \(\exists \,\,i\,\in \{1,\ldots , n\}\) with \(0\notin G_i\) and S has no border facets, then S can be rejected.

The argument is that the relative boundary of S may contain a global minimum point, but the same point is enclosed in the relative boundary of another partition set.

If an \(n-\)simplex S is not interior to \(\Delta \), it can be replaced by its border facets. However, not all border facets are interesting partition sets. It is possible to show that the relative interior of a facet F cannot contain a minimum point. For this property, we could use the so-called normal vector p with respect to facet F. Let \(S= {{\,\textrm{conv}\,}}(\mathcal {V})\), \(v\in \mathcal {V}\), \(\mathcal {U}=\mathcal {V}\setminus \{v\}\) and facet \(F= {{\,\textrm{conv}\,}}(\mathcal {U})\). Consider \(u\in \mathcal {U}\) and construct a set of vectors \(\mathcal {Y}:=\{w-u, w\in \mathcal {U}\setminus \{u\} \}\) with corresponding matrix representation Y. The normal vector of F in the direction of v can be expressed by defining (see e.g. [3])

$$\begin{aligned} p:=\left[ I- Y(Y^TY)^{-1}Y^T\right] (v-u). \end{aligned}$$
(4)

The projection procedure (4) is sketched in Fig. 2. Vector p points into the same half-space where S is located. From each point in \({{\,\textrm{rint}\,}}(F)\) a direction d with \(p^Td>0\) points into \({{\,\textrm{rint}\,}}(S)\). Considering the gradient enclosure \(G{\ (\supseteq \nabla f(x), \forall x \in S)}\) component-wise, we can come to several conclusions.

Fig. 2
figure 2

Normal vector calculation via projection of an edge \(v-u\)

Proposition 3

Consider an \(n-\)simplex \(S={{\,\textrm{conv}\,}}(\mathcal {V})\) with a facet F generated by removing vertex v from \(\mathcal {V}\) with normal vector p as in (4). If for one component i we have \(G_ip_i<0\) (either \(\overline{G}_i<0\) and \(p_i>0\), or \(\underline{G}_i>0\) and \(p_i<0\)), then \({{\,\textrm{rint}\,}}(F)\) cannot contain a minimum point of (1).

Proof

See Proposition 4.3 in [9]. \(\square \)

In the B &B procedure, if the search space is a full-dimensional simplex, the partition sets are full-dimensional \(n-\)simplices, up to a dimension reduction occurs to its facets. Using the property of \(0\notin G\), if \({{\,\textrm{rint}\,}}(S)\) can be removed, we can reduce it to its border facets. However, for lower dimension \(m-\)simplices, \(m<n\), the condition \(0\notin G\) is a necessary condition to have monotonicity, but not sufficient. Therefore, we have to focus on bounds of directional derivatives. The most general result for an \(m-\)simplex is the following.

Proposition 4

Let \(S\subseteq \Delta \) be an \(m-\)simplex with \(m\le n\) and gradient enclosure \({G\ \supseteq \nabla f(x), \forall x \in S}\). If \(\exists \, x,y\in S\), such that direction \(d=x-y\) has corresponding directional derivative bounds (3) with \(0\notin [\underline{d^TG},\overline{d^TG}]\) then \({{\,\textrm{rint}\,}}(S)\) does not contain a global minimum point of (1).

Proof

Consider \(z\in {{\,\textrm{rint}\,}}(S)\). As z is in the relative interior, there exists a feasible direction d in which lower function values can be found, i.e. \(\exists \, \epsilon \in \mathbb {R}\) small enough, such that \(z+\epsilon d \in S\) and \(f(z+\epsilon d)<f(z)\). So z cannot be a minimum point of f. \(\square \)

The elaboration for the algorithm depends on the choice of the direction d and the way to compute it.

Corollary 2

Let \(S\subseteq \Delta \) be an \(m-\)simplex as partition set in a branch and bound algorithm with corresponding gradient enclosure \({G\ \supseteq \nabla f(x), \forall x \in S}\). If the conditions of Proposition 4 apply and S has no border facets, then S can be rejected.

The argument is that the relative border of S may contain a global minimum point, but the same point is enclosed in the relative boundary of another partition set.

Given that the minimum is not in \({{\,\textrm{rint}\,}}(S)\), we have to decide which of the facets to focus on.

Proposition 5

Given \(m-\)simplex \(S={{\,\textrm{conv}\,}}(\mathcal {V})\), \(1\le m\le n\) and a facet F generated by removing vertex v from \(\mathcal {V}\), consider direction \(d=v-y\), where \(y\in S\). Consider the enclosure of the gradient \(G\ \supseteq \nabla f(x), \forall x \in S\). If \(\underline{d^TG}>0\), i.e. the directional derivative in the direction \(d=v-y\) is positive, then the facet F contains all minimum points in S, i.e. \({{\,\textrm{argmin}\,}}_{x\in S}f(x)\subseteq F\).

Proof

Consider an ordered vertex set \(\mathcal{V} =\{v_0,\ldots ,v_m\}\) and a facet \(F={{\,\textrm{conv}\,}}(\{v_1,\ldots ,v_m\})\), so the removed vertex is \(v_0\). First, we show that \(\forall x\in S\setminus F,\ \exists \, \rho > 0\) such that \(z=x-\rho d\in F\). We have \(x,y\in S\), such that there exist vectors \(\lambda ,\gamma \) with \(x=\sum _{j=0}^m \lambda _j v_j\) with \(\lambda _0>0\) and \(y=\sum _{j=0}^m \gamma _j v_j\). Looking for the value \(\rho \), such that \(z=x-\rho (v_0-y)\in F\), we have that \(\lambda _0-\rho (1-\gamma _0)=0\). Hence, taking \(\rho =\frac{\lambda _0}{1-\gamma _0}\), we obtain \(z\in F\). According to the mean value theorem, there exists a point \(\xi \) in between x and z such that

$$ f(x) = f(z+\rho d) = f(z) +\rho d^T\nabla f(\xi ) \ge f(z) + \rho \underline{d^TG} > f(z). $$

This means that each point \(x\in S\setminus F\) corresponds to a point \(z\in F\) with a lower function value.\(\square \)

Remark 1

Notice that if there are several facets fulfilling the condition of Proposition 5, then the consequence is also valid for the intersection of the facets, which is a face of S. So if Proposition 5 applies to F and H, then it also applies to \(F\cap H\).

Proposition 6

Given \(m-\)simplex S fulfilling the conditions in Proposition 5. If F is a non-border facet, then global minimum points of (1) cannot be in \({{\,\textrm{rint}\,}}(F)\), i.e. global minimum points can only be found at the relative boundary of F.

Proof

According to Proposition 5, the minimum points of f over S can be found on facet F. We first show that \(x\in {{\,\textrm{rint}\,}}(F)\) does not contain global minimum points. Consider a point \(x\in {{\,\textrm{rint}\,}}(F)\). Because F is non-border, there exists an \(\epsilon >0\), such that \(z(\epsilon ):=x-\epsilon (v-y)\in \Delta \) and \(z(\epsilon )\notin S\). As \(\underline{d^T G}>0\) and f is continuously differentiable, there exists an \(\epsilon >0\) (sufficiently small) such that for all \(\delta \in [0,\epsilon ]\), we have \(d^T \nabla f(z(\delta ))>0\). Therefore, following the mean value theorem, there exists a point \(\delta \in [0,\epsilon ]\) such that \(f(z(\epsilon ))=f(x)-\epsilon (v-y)^T\nabla f(z(\delta ))<f(x)\). This means that in a neighbourhood of x outside S, we can find a feasible point with a lower function value than x. The relative boundary of F could contain global minimum points of (1).\(\square \)

Remark 2

First, let us note that the boundary of F is shared by other partition sets, thus we can disregard simplex S in a branch and bound process.

Second, notice that there can be several facets fulfilling the condition of Proposition 5. However, we need only one non-border facet fulfilling the condition to remove simplex S from consideration.

Corollary 3

Let \(m-\)simplex S with centroid c fulfil the conditions of Proposition 5. Consider direction \(d=v-c\). Consider the enclosure of the gradient \(G\ \supseteq \nabla f(x), \forall x \in S\). If \(\underline{d^TG}>0\), then facet F contains all minimum points in S, i.e. \({{\,\textrm{argmin}\,}}_{x\in S}f(x)\subseteq F\). Moreover, if F is a non-border facet, then simplex S cannot contain the global minimum of (1).

This follows from Propositions 5 and 6 by using the centroid \(y=c\).

Example 1

For the illustration of the concept, consider the simplices in Fig. 1. It shows three partition sets generated by bisection, using bisection points \(x_3\) and \(x_4\). Consider \(S={{\,\textrm{conv}\,}}(\mathcal {V})\) with \(\mathcal {V}=\{v_1,x_3,x_4\}\), where we assume that the orange direction provides a monotonously increasing direction. According to Proposition 4, the interior of S does not contain a global minimum point. The blue direction provides an increasing monotonicity for facet \(F={{\,\textrm{conv}\,}}(\mathcal {V}\setminus \{v\})\) with \(v=x_4\). According to Proposition 5, facet \(F={{\,\textrm{conv}\,}}(\{v_1,x_3\})\) contains all minimum points of S. Specifically, they are at lower dimensional faces of F (the vertices), as \({{\,\textrm{rint}\,}}(F)\) cannot contain minima either. Now, Proposition 6 states, that we even can remove S, as there is another partition set sharing F, that encloses all minimum points.

An alternative is to use the projection point \(y=(v-p)\in S\) of the normal vector p of (4) as illustrated in Fig. 2 for the choice of y. However, this is only correct if the simplex has angles of less than \(\frac{\pi }{2}\). This is typically the case when we follow the longest edge bisection process.

Corollary 4

Given an \(m-\)simplex partition set \(S={{\,\textrm{conv}\,}}(\mathcal {V})\) in the bisection and dimension reduction process and a facet F generated by removing vertex v from \(\mathcal {V}\) with normal vector p according to (4). Consider the enclosure of the gradient \(G\ \supseteq \nabla f(x), \forall x \in S\). If \(\underline{p^TG}>0\), then F contains all minimum points over S. Moreover, if F is a non-border facet, then simplex S cannot contain the global minimum of (1).

This follows from Propositions 5 and 6 by using the normal vector p in place of \(d=v-y\).

Proposition 7

Given \(m-\)simplex \(S={{\,\textrm{conv}\,}}(\mathcal {V})\), \(1\le m\le n\), consider a facet F generated by removing vertex v from \(\mathcal {V}\) and gradient enclosure \(G \ \supseteq \nabla f(x), \forall x \in S\). If for direction \(d=v-y\) with \(y\in S\) we have \(\overline{d^TG}<0\), then \({{\,\textrm{rint}\,}}(F)\subset S\) does not contain a minimum point of (1).

Proof

Following the proof of Proposition 5, considering that for any \(\xi \in S\), \(\nabla f(\xi ) \le \overline{d^TG}\), we have

$$ f(x) = f(z+\rho d) = f(z) +\rho d^T\nabla f(\xi ) \le f(z) + \rho \overline{d^TG} < f(z) $$

and the result follows. \(\square \)

3 Checking on Monotonicity

To prove that there exists a monotone direction in an \(n-\)simplex, it is sufficient to have \(0\notin G\). This is a necessary, but not sufficient condition for a lower dimensional \(m-\)simplex, \(m<n\). To prove that such a direction exists, according to Proposition 4, we need to find a direction \(d=x-y\), with \(x,y \in S\) corresponding to a positive lower bound of the directional derivative, \(\underline{d^TG}>0\), where

$$\begin{aligned} \underline{d^TG} = \sum _{i=1}^{n} \min \{d_i\underline{G}_i, d_i \overline{G}_i\}. \end{aligned}$$
(5)

Note that it is enough to find a direction d with positive directional derivative bound, as \(-d\) will give a negative directional derivative. Finding such a direction can be done by searching for the steepest monotone direction \(\max _{d} \underline{d^TG}\) by an LP and check if it is positive. We reformulate the \(\min \) part by introducing variable \(z_i=\min \{d_i\underline{G}_i, d_i \overline{G}_i\}\). In order to decrease the number of variables, we can fix x to the centroid in the directional vector, i.e. \(d=c-y\), and solve

$$\begin{aligned} {\begin{matrix} &{}\max \sum _{i=1}^n z_{i} \\ {{\,\mathrm{\quad s.t. \quad }\,}}&{} z_{i}\le (c_i-y_i)\underline{G}_i,\ i=1,\ldots ,n\\ &{} z_{i}\le (c_i-y_i)\overline{G}_i,\ i=1,\ldots ,n\\ &{} \sum _{j=1}^m \lambda _j=1 \\ &{} y=\sum _{j=0}^m\lambda _jv_j\\ &{}\lambda _j\ge 0, j=0,\ldots ,m \end{matrix}} \end{aligned}$$
(6)

If there is no monotone direction, then \(y=c\) and \(\sum _{i=1}^m z_{i}=0\). If a monotone direction exists, y will lay on the relative boundary of S to maximize the directional derivative.

After proving a positive direction exists, one can validate Proposition 5 for a specific facet F of S. Consider the vertex set of F given as \(\{v_1,\ldots ,v_m\}:=\mathcal {V}\setminus \{v_0\}\). Focusing on the direction \(v_0-y\) with \(y=\sum _{j=1}^m \lambda _jv_j\), we can demonstrate there is a (maximum) positive directional derivative solving the LP

$$\begin{aligned} {\begin{matrix} &{}\max \sum _{i=1}^n z_{i} \\ {{\,\mathrm{\quad s.t. \quad }\,}}&{} z_{i}\le \underline{G}_i (v_{0i}-y_i),\ i=1,\ldots ,n\\ &{} z_{i}\le \overline{G}_i (v_{0i}-y_i),\ i=1,\ldots ,n\\ &{} \sum _{j=1}^m \lambda _j=1 \\ &{} y=\sum _{j=1}^m \lambda _jv_{j} \\ &{}\lambda _j\ge 0, j=1,\ldots ,m \end{matrix}} \end{aligned}$$
(7)

and checking the result is positive. It is the question whether it is worthwhile to solve (7) for each facet of a partition set. The focus is on having a direction \(d=v_0-y\) where the function is monotonic and to consider the corresponding facet.

We can also search for a facet according to (7) over all facets simultaneously, solving an Integer Programming problem where a binary variable \(\delta _j\) selects the facet corresponding to the most positive directional derivative.

$$\begin{aligned}&\max \sum _{i=1}^n z_{i}\nonumber \\ {{\,\mathrm{\quad s.t. \quad }\,}}&z_{i}\le \underline{G}_i d_i,\ i=1,\ldots ,n\nonumber \\&z_{i}\le \overline{G}_i d_i,\ i=1,\ldots ,n\nonumber \\&\sum _{j=0}^m \lambda _j=1 \\&d = \sum _{j=0}^m \delta _jv_j -\sum _{j=0}^m\lambda _j v_{j}\nonumber \\&\sum _{j=0}^m \delta _j = 1\nonumber \\&\lambda _j\ge 0, j=0,\ldots ,m\nonumber \\&\delta _j\in \{0, 1\}, j=0,\ldots ,m\nonumber \end{aligned}$$
(8)

Problem (8) can be used to replace LP problem (7) if such a monotonously increasing direction towards one of the vertices exists. If this is not the case, there may still exist an increasing direction according to (6).

The two steps, looking whether a monotone direction exists and identifying which facet contains all the minima (if any), can also be done in one step. Thus, instead of solving LP (6) and (7) \(m+1\) times (or MILP (8)), we can solve directly one MILP as follows. Let \(d=x-y\), where \(x,y\in S\), so \(x=\sum _{j=0}^m \lambda _j v_j\) and \(y=\sum _{j=0}^m \mu _j v_j\). Consider direction \(d=x-y=\sum _{j=0}^m (\lambda _j-\mu _j)v_j\)

$$\begin{aligned}&\max \sum _{j=0}^m \delta _j \nonumber \\ {{\,\mathrm{\quad s.t. \quad }\,}}&z_{i}\le \underline{G}_i d_i,\ i=1,\ldots ,n\nonumber \\&z_{i}\le \overline{G}_i d_i,\ i=1,\ldots ,n\nonumber \\&\sum _{j=0}^m \lambda _j=1 \nonumber \\&\sum _{j=0}^m \mu _j=1\nonumber \\&d = \sum _{j=0}^m (\lambda _j-\mu _j)v_j \\&\sum _{i=1}^n z_{i} \ge \epsilon \nonumber \\&\mu _j \ge \delta _j \quad \forall j=0,\ldots , m\nonumber \\&\sum _{j=0}^m \delta _j \le 1.\nonumber \\&\mu _j, \lambda _j\in \left[ 0, 1\right] , j=0,\ldots ,m\nonumber \\&\delta _j\in \{0, 1\}, j=0,\ldots ,m\nonumber \end{aligned}$$
(9)

A solution with \(\mu _k=1\) and \(\mu _j= 0, j\ne k\) represents a direction d pointing to vertex \(v_k\). A solution with \(\lambda _k=1\) and \(\lambda _j=0, j\ne k\) represents a direction pointing from vertex \(v_k\). We connect \(\mu _j\) with binary variables \(\delta _j\) such that \(\delta _j=1\) implies \(\mu _j=1\). The inequality \(\sum _{i=1}^n z_{i} \ge \epsilon \) with \(\epsilon >0\) assures d is a monotone direction. Moreover, it forces \(\mu \ne \lambda \), because otherwise \(z=0\) as well. If no monotone direction exists, (9) has no feasible solution.

The objective is to maximize the sum of \(\delta _j\), meaning that we aim at finding a monotone direction with \(\mu _k=1\) corresponding to \(\delta _k=1\). In this case, we know that facet \(F={{\,\textrm{conv}\,}}(\mathcal{V}\setminus \{v_k\})\) contains all minima according to Proposition 5. If the objective is zero, there is no facet containing all the minima, i.e. there is no \(\delta _k=1\), but there is a monotone direction d.

4 Keeping Track of Border Facets

Fig. 3
figure 3

Faces of the feasible set with bisection points \(x_3\) and \(x_4\)

To know the border status of a given facet, we use a labelling system to find out which minimum dimensional face of \(\Delta \) the facet or simplex is included in. This is done by assigning to each face \(\varphi \) of feasible set \(\Delta ={{\,\textrm{conv}\,}}(\mathcal {W})\) a label \(\mathcal {B}(\varphi )\), starting with the vertex faces labeled \(\mathcal {B}(v_j)=0..010..0\) where the only 1 corresponds to bit j, for \(j=0,\ldots ,n\). For each face \(\varphi \), which is a convex combination of vertices \(\mathcal{V}\subseteq \mathcal{W}\), the corresponding bit-string \(\mathcal {B}(\varphi )\) has a value 1 for each vertex \(v\in \mathcal{V}\) in the same position as in \(\mathcal {B}(v)\). For instance, in Fig. 3, the edge \((v_1,v_2)\) has label 011, and the simplex \(\Delta \) has label 111. In fact, for an \(m-\)simplex face \(\varphi \) of the feasible set, its label is given by the bitor of the label of all its vertices. The complete face graph is given in Fig. 4 for a 4-vertex simplex.

Fig. 4
figure 4

Face graph for 4-vertex simplex. Binary labels in boxes for each face

After bisecting the original set of vertices \(\mathcal {W}\), we store the bisection points in set X. For instance in Fig. 3, we have \(X=\{v_0, v_1, v_2, x_3, x_4\}\). We label all generated points \(x\in X\), which serve as vertices for the partition sets. The label of point x is the same as the label of the minimum dimensional face \(\varphi \) of the feasible set x is in. For instance, the label of \(x_3\) is 101, the same as the label of face \({{\,\textrm{conv}\,}}(v_0,v_2)\). In the bisection, a new vertex \(x=\frac{v+w}{2}\) gets label \(\mathcal {B}(x)={\textbf {BitOr}} (\mathcal {B}(v),\mathcal {B}(w)).\)

Given an \(m-\)simplex \(S={{\,\textrm{conv}\,}}(\mathcal {V})\), the question is what is the label of the (smallest dimensional) face \(\varphi \) it is included in. This is determined by the label \(\mathcal {B}(\varphi )={\textbf {BitOr}} (\mathcal {B}(\mathcal {V}))\), to be interpreted as a bitor on all its vertex labels. The number of ones of a bitstring \(\mathcal {B}(\varphi )\) is denoted by \(|\mathcal {B}(\varphi )|\), giving the number of vertices of \(\varphi \). According to Definition 1, \(m-\)simplex S is border if there exists an \(m-\)simplex face \(\varphi \) of the feasible set including S (\(m<n\)).

Proposition 8

Consider an \(m-\)simplex \(S={{\,\textrm{conv}\,}}(\mathcal {V})\) with \(m<n\).

If \(|{\textbf {BitOr}} (\mathcal {B}(\mathcal {V}))|=m+1\), then simplex S is border.

Proof

Consider the face \(\varphi \), which is the minimal dimensional face containing S, i.e. label \(\mathcal {B}(\varphi )={\textbf {BitOr}} (\mathcal {B}(\mathcal {V}))\). As \(|{\textbf {BitOr}} (\mathcal {B}(\mathcal {V}))|=m+1\), we have \(|\mathcal {B}(\varphi )|=m+1\), thus \(\varphi \) is an \(m-\)simplex. Therefore, S is enclosed by an \(m-\)simplex face of the feasible set and thus is border. \(\square \)

For example in Fig. 3, the edge \((x_3,v_2)\) is border, as \(|{\textbf {BitOr}} (\mathcal {B}(\{x_3,v_2\}))|=|{\textbf {BitOr}} (101,001)|=|101|=2\) corresponding to face \({{\,\textrm{conv}\,}}(\{v_0,v_2\})\). In contrast, edge \((x_3,x_4)\) is not border, because \(|{\textbf {BitOr}} (\mathcal {B}(\{x_3,x_4\}))|=|111|=3\ne 2\). In fact, the minimum dimensional face it is included in, is \(\Delta \) itself.

As we have seen in the theoretical results, for monotonicity testing the border certification of facets is of importance. One way to do so is to follow Proposition  8. One of our research questions is whether it might be more efficient in an algorithm to store the border status of the facets of a simplex in a boolean vector \(\beta \) indicating for each facet whether it is border or not. Then, we can transfer the status \(\beta \) to facets of newly generated simplicial subsets instead of evaluating the border status \(|{\textbf {BitOr}} (\mathcal {B}(\mathcal {V}))|=m+1\) for all facets. Therefore, we introduce a border status for each facet of \(m-\)simplicial subset S. Each facet F corresponds to a vertex \(v\in \mathcal {V}\) which is dropped from the vertex set S to obtain F. Let \(F_v\) denote \({{\,\textrm{conv}\,}}(\mathcal {V}\setminus \{v\})\) with its marker \(\beta (F_v)\in \{true, false\}\). Storing this information with subset S requires storing a binary vector \(\beta =(\beta _0, \ldots , \beta _m)=(\beta (F_{v_0}), \ldots , \beta (F_{v_m}))\) ordered in the same way as the vertex set and the corresponding facets.

4.1 Setting the Border Status for Facets in Longest Edge Bisection Refinement

Fig. 5
figure 5

Part of the search tree where LEB splits simplex S into two new simplices. Blue edges are border and red ones are non-border. The border status of vertex \(v_j\) and corresponding facet \(F_j\) is given by j(t) if \(\beta (F_{v_{j}}\!) = true\) and j(f) if \(\beta (F_{v_{j}}\!) = false\). The dotted facet is the new facet generated by the bisection

Consider \(m-\)simplex S with ordered vertex set \(\{v_0, \ldots , v_i, \ldots , v_j, \ldots , v_m\}\) and the corresponding facet border status vector \((\beta _0,\ldots ,\beta _i,\ldots ,\beta _j,\ldots ,\beta _m)\). Longest Edge Bisection (LEB) along edge \((v_i,v_j)\) generates the new point \(x=\frac{v_i+v_j}{2}\) and two new partition sets \(S_1={{\,\textrm{conv}\,}}(\{v_0,\ldots ,x,\ldots ,v_j,\ldots ,v_m\})\) and \(S_2={{\,\textrm{conv}\,}}(\{v_0,\ldots ,v_i,\ldots ,x\ldots ,v_m\})\). Both new simplices inherit the border status \(\beta \) from parent S, but we set the border status of the new facet to false: for \(S_1\), we set \(\beta _j:=false\) and for \(S_2\) we put \(\beta _i:=false\). So, apart from the newly generated shared facet, all other facets inherit the border status from the bisected simplex. Figure 5 shows the first bisection of the feasible set being 3-simplex \(\Delta \) with the border status of its facets.

4.2 Setting Border Status for Facets in Simplicial Dimension Reduction

Consider \(m-\)simplicial subset S with ordered vertex set \(\{v_0,v_1,\ldots ,v_m\}\) and corresponding border status \((\beta _0,\beta _1,\ldots ,\beta _m)\). Let us focus on dropping one of the vertices and reducing the simplex to the corresponding border facet. Without loss of generality consider dropping \(v_1\). The remaining facet \(F_{v_1}= {{\,\textrm{conv}\,}}(\{v_0,v_2\ldots ,v_m\})\) is border when \(\beta _1=true\). Otherwise \(F_{v_1}\) is not border. What will be the new border status \((\beta _0,\beta _2\ldots ,\beta _m)\)?

Proposition 9

Given \(m-\)simplex \(S={{\,\textrm{conv}\,}}(\mathcal {V})\) and border facets \(F={{\,\textrm{conv}\,}}(\mathcal {V}\setminus \{v\})\) and \(E={{\,\textrm{conv}\,}}(\mathcal {V}\setminus \{w\}), v\ne w\). The intersection \(F\cap E={{\,\textrm{conv}\,}}(\mathcal {V}\setminus \{v,w\})\) is border.

Proof

According to Proposition 8, the including faces \(\varphi _F,\varphi _E\) corresponding to F and E have labels such that \(|\mathcal {B}(\varphi _F)|=|\mathcal {B}(\varphi _E)|= m\), because F and E are border \((m-1)\)-simplices. The intersection \(F\cap E\) is an \((m-2)-\)simplex, as we removed two vertices from \(\mathcal {V}\).

Assume that \(F\cap E\) is not border, then \(|\mathcal {B}(\varphi _{F\cap E})|>m-1\). As \(\varphi _{F\cap E}\subset \varphi _F\), we have that \(|\mathcal {B}(\varphi _{F\cap E})|\le m\) bits. Consequently, \(|\mathcal {B}(\varphi _{F\cap E})|= m\).

As \((F\cap E)\cup \{w\}=F\), also \({\textbf {BitOr}} (\mathcal {B}(\varphi _{F\cap E}),\mathcal {B}(w)) = \mathcal {B}(\varphi _F)\). From \(|\mathcal {B}(\varphi _{F\cap E})|= |\mathcal {B}(\varphi _{F})|= m\) follows that \(\mathcal {B}(\varphi _{F\cap E}) = \mathcal {B}(\varphi _F)\). Similarly, we obtain \(\mathcal {B}(\varphi _{F\cap E}) = \mathcal {B}(\varphi _E)\). Therefore, \(\mathcal {B}(\varphi _F)=\mathcal {B}(\varphi _E)\) and \(\varphi _E=\varphi _F\). Two different border facets F and E cannot be included in the same m-face of the feasible set. Thus, \(F\cap E\) must be border.\(\square \)

The consequence of Proposition 9 is that a border status \(\beta (F_v)=true\) remains true after reduction and \(\beta (F_v)=false\) requires recalculation using Proposition 8. Figure 6 shows an illustrative example.

Fig. 6
figure 6

Reduction of \(S_{15}\) by removing vertex \(v_1\) (denoted by 1(t) on the left). Before reduction, \(S_{15}\) is included in face \(\varphi \) with label \({\textbf {BitOr}} (\mathcal {B}(\mathcal {V}))=1111\). So \(S_{15}\) contains interior points in \(\Delta \). After reduction, \({\textbf {BitOr}} (\mathcal {B}(\mathcal {V}\setminus \{v_1\}))=1011\). The border status of facet \(F_{v_0}\) is true and stays the same. Vertex \(v_2\) (\(v_3\)) becomes vertex \(v_1\) (\(v_2\)). Because the border status of their facets is false, the new status requires recalculation

Algorithm 1
figure a

SimplicialBranch-and-Bound\((f,\ \Delta ,\ \alpha )\)

5 A Simplicial Branch and Bound Algorithm Exploiting Monotonicity

The former section described several ways to deal with the monotonicity and border certification. Basically, one can evaluate the border status of facets, or store their status within the simplex object. This section describes the branch and bound algorithm that exploits the theoretical monotonicity results. It is based on a traditional scheme where a set \(\Lambda \) of subsets includes all global minimum points. The stopping criterion is taken as the classical concept of [6] on the gap between the incumbent upper bound \(\overline{f}\) of the best point found so far and a lower bound \(\underline{f}\) of the lowest bound among all subsets with an accuracy \(\alpha \). This concept corresponds to guarantee to find a point which is not worse than an accuracy \(\alpha \) from the global minimum function value. Of course one can consider other criteria based on the distance from the minimum points or the relative function value. However, the focus is on investigating the impact of the different monotonicity rules.

There are several ways to define the lower bound \(\underline{f}\) over simplicial subsets [9]. We evaluate several ways to certify the monotonicity of the objective over a subset. Algorithm 2 describes the monotonicity test using center to vertex (C) and vertex to vertex (V) directions. In line 3, we generate the directions to check for all vertices: i) the center to vertex directions and ii) all vertex to vertex directions without repetition. For each direction, we check its monotonicity, and we either discard the simplex (line 7, return false), reduce to some of its facets (line 13 or 17, return false) or return true to divide simplex S.

Algorithm 2
figure b

MonotonicityTest\((f,\ S,\ G,\ \Lambda )\)

Algorithm 3 describes the monotonicity test using LP. This can be used instead of Algorithm 2, but it is better to use it when the CV based directions did not provide a monotone direction. Lines 6-14 seek for a monotone direction towards a vertex, such that we can reduce to as few facets as possible. These lines can be changed to solving MIP(8), and we will refer to exactly that by monotonicity test using MIP(8).

Algorithm 3
figure c

MonotonicityTest-LP\((f,\ S,\ G,\ \Lambda )\)

In Algorithm 4, the most compact version of the monotonicity test is described, where MIP(9) is solved, and from its solution (or lack of it), we know if the simplex can be removed, reduced, or should be divided.

Algorithm 4
figure d

MonotonicityTest-MIP\((f,\ S,\ G,\ \Lambda )\)

6 Numerical Illustration

Algorithms were run on an Intel(R) Core(TM) i7-4770K CPU and 16GB of RAM running Fedora 37 Linux distribution. The algorithm was coded with g++ (gcc version 11.2.1-4) and it uses Kv\(-\)0.4.55 for Interval Arithmetic and Automatic Differentiation. Kv uses Boost libraries. Algorithms were compiled with-O3 -DNDEBUG -DKV_FASTROUND options. For the Linear Programming, we use PNL 1.11.0 compiled with-DCMAKE_BUILD_TYPE = Release and-DWITH_MPI = OFF, as a C++ wrapper to LPsolve 5.5.2.0-29. The default Lp_solve parameters are used in PNL. We also use the IBM ILOG CPLEX Studio 22.1.0 for LP and MIP models. All results were obtained with a termination accuracy \(\alpha =10^{-6}\).

In Sect. 6.1, we first use a small instance to illustrate why the use of a simplicial B &B (sBB) is sometimes much more efficient than using standard Interval based B &B (iBB) with linear constraints. Then in Sect. 6.2, we compare variants of the simplicial algorithm on test instances introduced in [8].

In all experiments, the directions from the centre to the vertices (C) are checked first. In case a border facet with a minimum has not been found, vertex to vertex directions (V) are checked. As soon as a direction with a positive directional derivative is found, the search is stopped.

6.1 Small Illustrative Instance sBB versus iBB

Interval branch and bound has a long tradition and knows many implementations with acceleration devices [2, 12]. Several variants exist that can handle linear constraints. In this experiment, we take a small instance from [10] where a quadratic convex function is considered over a simplicial feasible set. The vertex structure of the feasible set is simply translated to a set of linear inequalities for iBB. In instance Ex6-B (see Appendix A.2.1), we define the constraints such that the minimum point lies on the boundary of the feasible set. In instance Ex6-In (see Appendix A.2.1), the minimum point is interior.

The idea of the design of experiments is to have an instance as small as possible to focus on the difference of having an interior minimum in a full dimensional case and of having the minimum on the boundary. For the interior optimum, we expect an established iBB to be far more efficient. However, we would like to observe what happens if the optimum is at the boundary, like in the instances we will use in Sect. 6.2. We compare the behavior of a standard iBB using linear constraints with the simplicial sBB which has been presented.

Both iBB and sBB, use the same inclusion information in the experiments. Interval Arithmetic uses boxes (interval vectors). Therefore, the interval hull of a simplex S (\(\Box S\)) implies an overestimation of the area of S. The iBB has to deal with linear constraints. Linear constraints derived from a set of vertices may have errors in a computer due to rounding. So, an interval is considered feasible if the enclosure of all constraints is negative or smaller than a given precision. For the precision, we use \(10^{-14}\).

Fig. 7
figure 7

Graphical output of iBB on instance Ex6-B. Right graph is a zoom over the minimum

Fig. 8
figure 8

Graphical output of sBB on instance Ex6-B. Right graph is a zoom over the minimum

In this comparison, we use the centred form, which is the interval extension of the first-order Taylor expansion, using the centre of the set (+CFc), to try to improve objective function bounds that follow from direct use of Interval Arithmetic, see [8]. Algorithm sBB applies a monotonicity test based on C and V directional derivatives. Bounds of the derivatives obtained by Interval Arithmetic are component-wise. Therefore, they can be used directly on boxes in iBB.

In Figs. 7 and 8, we can see the generated sets for Ex6-B using iBB and sBB. It is clear that the number of sets generated by sBB is much smaller. Table 1 describes the meaning of the used colors in the figures. For a more comprehensive comparison Table 2 provides the numerical results of iBB and sBB for this small illustrative instance.

Table 1 Colors used in the Figures
Table 2 Comparing Interval versus Simplicial B &B for Ex6

Numerical results for Ex6-B, the test instance with the minimum at the boundary of the feasible set, show that simplex reduction in sBB provides a large advantage compared to the iBB results. At the final stages of the running algorithms, iBB must divide and check the feasibility of a lot of boxes near the global minimum while sBB is just dividing the reduced simplices (edges).

Numerical results for Ex6-In, the test instance with an interior minimum, show that somewhat more simplex than box evaluations are required. The main reason is that to reach the desired accuracy, a box in iBB and a simplex S in sBB (with bounds determined by enclosure \(\Box S\)), sBB requires more longest edge divisions than iBB and that number increases with the dimension. Another division method and devising specific and sharper bounding for simplices could balance the scales. See Figs. 9 and 10 for the graphical outputs.

Fig. 9
figure 9

Graphical output of iBB on instance Ex6-In. Right graph is a zoom over the minimum

Fig. 10
figure 10

Graphical output of sBB on instance Ex6-In. Right graph is a zoom over the minimum

6.2 Comparing Variants of the Monotonicity Test in the Simplicial B &B

The implemented algorithm uses the centred form based on the highest vertex as the center for the computation of the lower bound, denoted by +CFvs, that was introduced in [8]. Several variants of Algorithm 2 are described and tested. Several directions \((v-y)\) can be used. Additionally, the best direction can be returned by solving LP (7) in Algorithm 3 or a MIP in Algorithm 4. The symbols used for the different evaluated directional derivatives and models for sBB are given in Table 3.

Table 3 Notation for generated and evaluated directions

For CVLP and CVMIP, LP and MIP models are evaluated only if CV did not find a positive directional derivative. This is done to reduce the number of LP and MIP models to be solved due to their high computational time. To reduce the computational time, we sort the facets by having \(\min \{\underline{d_v^TG},\overline{d_v^TG}\}\) first, during the CV directions search. In case of a successful LP(6) search, LP(7) visits the facets in the predefined order and stops as soon as a positive directional derivative is found. So, some facets might not be checked. If LP(6) is infeasible, the simplex is divided and LP(7) or MIP(8) are not evaluated. If LP() or MIP(8) do not find a monotone directional derivative, then \(d_v^TG=0\).

LP(6) and LP(7) were solved using the PNL library and CPLEX in CVLP. Data in the LP model was reused as much as possible in PNL, but models are built each time and data is not reused in CPLEX for all models. The time limit for each model was set to 1 s and verbose was off, for both PNL (lp_solve) and CPLEX.

Table 4 Test problems. Traditional test functions optimized over a simplicial feasible set. The portfolio problem described as a quadratic program (UPQP) can be found in Appendix A.1. An asterisk at n indicates the selected dimension for a varying dimension test instance. A description of the test instances with corresponding minima can be found in www.sfu.ca/~ssurjano/optimization.html

We use the same set of test problems as in [8], given in Table 4. Traditional test functions are solved over a simplicial feasible set. The feasible simplex is determined such that the global optimizer is at the centre of its most left facet, and its size is half of the enclosing bounding box.

For instance, problem GP2, with enclosing box \(([-2,2],[-2,2])\) and minimum at \((0,-1)\), would have a corresponding feasible set \(\Delta ={{\,\textrm{conv}\,}}(\{(-0.5, -2), (1.5, -2),(0.5, 0.0)\})\). However, in order to avoid an early hit of a global minimum point due to LEB and the clustering effect [7], the instances were adapted moving the simplicial feasible area a bit, such that the minimum point does not coincide with the centroid of the facet F it is included in. The simplex was shifted such that the minimum gets closer to a vertex v of F with 0.1 times the distance between the centroid of F and vertex v. For instance, problem GP2, with simplicial feasible area \(\Delta ={{\,\textrm{conv}\,}}((-0.5, -2), (1.5, -2), (0.5, 0.0))\) and minimum at \((0,-1)\), is solved for the feasible set \(\Delta ={{\,\textrm{conv}\,}}((-0.55, -2.1), (1.45, -2.1), (0.45, 0.1))\). On one hand, the global minimum is not hit by one of the newly generated points and delays successful RangeUp and Cutoff tests. On the other hand, avoiding the clustering effect reduces the number of evaluated simplices. Results for instances with a shifted centroid are given in Tables 6, 7, 8, 9 and 10. The meaning of the symbols used in the tables to present the indicator values are given in Table 5.

Table 5 Meaning of headers in the following tables
Table 6 Results for the CV variant of the sBB algorithm

The numerical results illustrate the usual explosion of the number of sub-simplices as the dimension increases in spatial branch and bound with few accelerating devices. Specifically, we are using just monotonicity (Alg. 1, line 11), RangeUp (Alg. 1, line 21) and Cutoff (Alg. 1, line 18) tests. Another characteristic is that the number of simplices where \(0\not \in G\), the ratio (N/N2) is small compared to the number of evaluated simplices (ES) for most cases. This is intrinsic to the problem at hand, but it is also due to the overestimation in using interval hull \(\Box S\) to get bounds of f(S) and G using Interval Arithmetic and Automatic Differentiation. The number of simplices where a monotonic directional derivative was found (OkCV) was neither specifically large in comparison to the number N/N2. The reason is that only few easy directions were checked using CV and \(0\notin G\) for \(\Box S\) does not imply the existence of a monotonic directional derivative in S. Extreme cases are G7 (OkCV=0) and L8 where OkCV=343 out of 135,590. Still, the use of the new described monotonicity test reduces the search space for most instances. In general, the shifted minimum instances are more challenging to be solved. Again, the number of N/N2 in comparison to the number of ES can be considered low in general as well as between N/N2 and OkCV for the most challenging instances. Anyway, the use of the described monotonicity test reduces the search space in all instances. However, just one simplex is rejected for instance G7.

Table 7 shows that adding directions generated by LP(6) + LP(7) results in a slightly smaller number of evaluated simplices, but a higher running time. Instances ST5 and H6 could not be solved within 5 min.

Table 7 Results using monotonicity tests evaluating CV and LP generated directions using solver PNL (CVLP)

Using CPLEX instead of PNL (Table 8) results in a larger running time for most of the instances apart from instance ST5 that could now be solved in less than 5 min.

Table 8 Results using monotonicity tests evaluating CV and LP generated directions using solver CPLEX (CVLP)

Using CV + LP(6) + MIP(8) (see Table 9) does not improve the running time. Instances ST5, H6 and DP7 could not be solved in less than 5 min.

Table 9 Results using monotonicity tests with CVMIP8 directions

Table 10 illustrates an increase in running time when using MIP models instead of LP models or not using models at all. For some instances, the running time is better than that reported in Table 9, but for other instances, it is worse.

Table 10 Results using monotonicity tests with CVMIP9 directions

7 Conclusions

We found that the concept of a directional derivative given bounds on the gradient is of importance to decide if a simplicial feasible set can be reduced to some of its facets or eliminated from consideration. Moreover, the concept of a border simplex serving as a facet of a subset is relevant. The paper summarizes the theoretical results and their proofs. Furthermore, it has been described how solving a Linear Program or Mixed Integer Linear Program can be used to detect with certainty the monotonicity in a subset.

The methods using LPs and MIPs have a theoretical character and due to numerical aspects can be hard to handle when subsets get small. Solving LPs and MIPs may require more computing time than easier tests, which are more practical to handle. Several numerical experiments are reported on (adapted) benchmark instances from literature.

A potential improvement of the described method can be found in the overestimation of the range of the gradient over a simplex, because the interval hull of the simplex is used to apply Interval Arithmetic and Automatic differentiation. Better and probably more expensive methods to get sharper bounds can be studied.

Regarding the use of LP and MIP models to search for a positive directional derivative bound in a simplex, the numerical results show that the number of positive results can be improved at an increasing running time. Their effectiveness could be increased by having sharper bounds on the gradient. Due to the large amount of simplices generated in a spatial branch and bound with a small number of accelerating devices, the time spent on each simplex should be short. We plan to develop a simple local search for a positive directional derivative. This could be less effective, but may result in a faster algorithm.