# Learn Continuously, Act Discretely:
Hybrid Action-Space Reinforcement Learning For Optimal Execution

Feiyang Pan∗1 Tongzhe Zhang∗2 Ling Luo1 Jia He1 &Shuoling Liu2
1Huawei EI Innovation Lab
2E Fund Innovation Lab
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Optimal execution is a sequential decision-making problem for cost-saving in
algorithmic trading. Studies have found that reinforcement learning (RL) can
help decide the order-splitting sizes. However, a problem remains unsolved:
how to place limit orders at appropriate limit prices?

The key challenge lies in the “continuous-discrete duality” of the action
space. On the one hand, the continuous action space using percentage changes
in prices is preferred for generalization. On the other hand, the trader
eventually needs to choose limit prices discretely due to the existence of the
tick size, which requires specialization for every single stock with different
characteristics (e.g., the liquidity and the price range). So we need
continuous control for generalization and discrete control for specialization.
To this end, we propose a hybrid RL method to combine the advantages of both
of them. We first use a continuous control agent to scope an action subset,
then deploy a fine-grained agent to choose a specific limit price. Extensive
experiments show that our method has higher sample efficiency and better
training stability than existing RL algorithms and significantly outperforms
previous learning-based methods for order execution.

##  1 Introduction

Since computer algorithms encompass the whole trading process in modern
financial markets, optimal execution has attracted significant research
attention for decades. For the buy-side, optimal execution aims to save cost
and reduce market impact when executing large orders. It is often viewed as a
stochastic sequential decision-making process because large orders need to be
divided into a series of small sub-orders within a period of time to be filled
without influencing the market too much. Traders or algorithmic trading
systems should decide the prices and volumes of the split sub-orders during
this process.

There are two fundamental types of transactions: market orders and limit
orders. Market orders are transactions meant to be executed as quickly as
possible at the market price. Therefore, when using market orders, the agent
only decides the size of each sub-order. Many studies showed that modern RL is
powerful for such a problem Nevmyvaka et al. (2006); Hendricks and Wilcox
(2014); Ning et al. (2018); Fang et al. (2021).

However, there are two disadvantages for such methods: 1) market orders have
no control over the price and might be filled at disappointing prices if the
bid-ask spreads are wide, and 2) the volume of sub-orders given by RL agents
might fluctuate too much because they work as a sort of market timing,
potentially causing large market impact. So in this paper, in contrast to
previous work, we focus on optimal execution with limit orders rather than
market orders.

A limit order is an order to buy or sell a stock at a specific price or
better, which can only be filled if the stock’s market price reaches the limit
price. Therefore, an agent for optimal limit order placement should control
both limit prices and volumes of sub-orders.

A challenge emerges when using a learning-based agent to decide the limit
price: the action space can be viewed as both discrete and continuous. From
the perspective of the original problem, the choices of limit prices are
discrete due to the existence of the tick size (the minimum price increment).
For example, for a $10.00 stock, one may set bid prices at $10.00 or $9.99,
but not at $9.995. So an ideal end-to-end agent should use discrete control to
select the actual prices precisely, but it is usually impossible to generalize
between stocks with different price-level and liquidity. For example, an
action of $-1$ tick means $-1\%$ for a $1.00 stock but $-0.01\%$ for a $100.00
stock.

On the other hand, people usually use percentages to understand the rise and
fall of stock prices, which is a continuous representation. As training deep
RL agents with high-dimensional and noisy state dynamics requires a lot of
trial-and-error, continuous control that determines the percentage change is
preferred for its generalization ability across stocks with different prices.
For example, a bid at -0.1% results in $9.99 for a $10.00 stock, and results
in $99.90 for a $100.00 stock. However, such a method is more inclined to
learn the strategy of stocks with high prices and liquidity, and fails to
specialize in determining the ticks for illiquid low-priced stocks.

In light of these observations, we propose Hybrid Action-space Limit Order
Placement (HALOP), which trains an agent to learn continuously and act
discretely. First, for generalization, HALOP employs a continuous-control
agent to scope an action subset with high-dimensional market dynamics as
input. Next, a fine-grained discrete-control agent chooses a specific discrete
limit price from the action subset, which specializes well in stocks with
different characteristics. Both agents are trained end-to-end to maximize the
long-term return, i.e., the excess return in execution.

We conduct extensive experiments and show that our new optimal execution
method with limit orders can beat the market, and significantly improves the
excess return upon previous learning-based order execution methods.

##  2 Related Work

###  2.1 Optimal Execution

####  2.1.1 Non-machine-learning Order Execution

Traditional optimal execution research tends to assume that the market price
movements follow some stochastic process such as the Brownian motion, and then
use stochastic optimal control method to derive the volume trajectory
analytically Almgren and Chriss (2001); Huberman and Stanzl (2000); Bertsimas
and Lo (1998). However, practitioners rarely use these methods to place orders
in stock markets, because the assumptions may not hold in the real-world. The
most widely used trade execution strategies are based on pure rules or
statistical rules. For example, the time-weighted average price (TWAP)
Bertsimas and Lo (1998) strategy divides a large order into equal-sized sub-
orders and executes each within equally divided time intervals.

The volume-weighted average price (VWAP) Kakade et al. (2004) strategy first
estimates the average volume traded for each time interval from historical
data, then splits order based on the estimates.

Although simple enough, TWAP and VWAP are still popular these days because
their execution costs are always close to the market. Specifically, we let the
TWAP strategy with market order be the benchmark in this paper.

####  2.1.2 Reinforcement Learning for Order Execution

As order execution is fundamentally a problem of making decisions under
uncertainty and the actual trading data is full of noise, RL becomes the best
choice to solve such a problem. Nevmyvaka Nevmyvaka et al. (2006) is the
pioneer for RL in optimal execution where the agent is trained by Q-Learning
Watkins and Dayan (1992) to chooses the limit price. However, they only
considers a small limited set of discrete actions (a few bid and ask prices)
without considering the more important percent changes in actions.

With the development of deep learning and deep RL in the past few years, some
studies use deep RL to learn to execute orders with high-dimensional market
data as inputs Hendricks and Wilcox (2014); Ning et al. (2018); Lin and Beling
(2019, 2020); Fang et al. (2021). While Hendricks and Wilcox (2014) applies RL
to modify a given volume trajectory suggested by the traditional order
execution model Almgren-Chriss Hendricks and Wilcox (2014). Without any market
assumptions, several variations of Deep Q-Network (DQN) Mnih et al. (2013) are
proposed to choose discrete volumes, which could address the high dimensions
and the complexity of the finance market and trading signals with the deep
neural network Ning et al. (2018); Lin and Beling (2019). Instead of manually
designed attributes, a PPO-based optimal execution framework is designed to
make decisions based on raw level-2 market data Lin and Beling (2020). Policy
distillation paradigm is also deployed in order execution Fang et al. (2021),
which uses a distilled PPO agent to choose optimal order-splitting volumes.
However, most of the previous work studies the problem of volume-splitting,
which only addresses the order execution problem partially.

##  3 Hybrid Action-Space Optimal Execution

###  3.1 Problem Setting

In this part, we first describe the general setting of the optimal execution
problem with limit orders, then formulate it as a sequential Markov Decision
Process (MDP). In particular, we focus on setting the limit prices of sub-
orders rather than the volumes for clarity, although the proposed method can
work seamlessly together with previous volume-oriented methods.

####  3.1.1 General Problem Description

We simplify the problem as a discrete-time decision process following Cartea
et al. (2015); Ning et al. (2018); Fang et al. (2021), where the agent
interacts with the environment multiple times within a predefined time period
based on discrete time space. So consider a horizon of $T$ timesteps
$t=1,\dots,T$, each represents a $1/T$ of a predefined total time periods. For
example, if there is a total time period of three hours and $T=90$, it has
$90$ small time periods of two minutes. At the $t$-th step (the beginning of
the $t$-th small time period), the environment reveals the recent market
dynamics to the agent, consisting of a series of history quotes and top 5
bid/ask prices and volumes. Then, the agent sends a limit order with a limit
price rounded to the tick size and a volume. Next, at the end of the $1/T$
period, the environment tells the agent whether (or how much) the order has
been filled, and withdraw the rest of it. Then it moves to the $(t+1)$-th
timestep. Specifically, at the end of the $T$-th period, if there is still
inventory remained, the system automatically send a market order to execute it
at once.

Figure 1: Illustration of different types of policies for different types of
action-spaces.

In addition, as we focus on the price, we assume there is a predefined
schedule of trading volume for the agent, denoted by
$v^{*}_{1},\dots,v^{*}_{T}$, which sums to one $\sum_{t=1}^{T}v^{*}_{t}=1$.
For example, it can be the volumes of classic financial strategies like TWAP
or VWAP, or the outputs of any volume-oriented algorithmic trading agent Fang
et al. (2021). As the limit orders are not guaranteed to be executed, let
$\tilde{v}_{1},\dots,\tilde{v}_{T}$ be the executed volume within each period,
which has
$\Delta_{t}=\sum_{\tau=1}^{t}v^{*}_{\tau}-\sum_{\tau=1}^{t}\tilde{v}_{\tau}\geq
0$. So at step $t$, we assume that the agent tries to catch up with the
schedule, by placing a new order with a volume

$v_{t}=v^{*}_{t}+\Delta_{t-1}.$ (1)

####  3.1.2 MDP Formulations for Limit Order Placement

Now we formulate the problem as a finite-time MDP, and put forward the
continuous-discrete duality of the problem.

The state consists of two types of information, the public state and the
private state $s_{t}=(s^{\textrm{pub}}_{t},s^{\textrm{priv}}_{t})$. Similar to
Nevmyvaka et al. (2006); Lin and Beling (2020); Fang et al. (2021), the public
state includes past market dynamics obtained from the environment. Private
observations are real-time status of the running agent, including the
remaining inventory $1-\sum_{\tau=1}^{t-1}\tilde{v}_{\tau}$; the amount of
unfulfilled orders $\Delta_{t-1}$; and the remaining time fraction $1-t/T$.

The action, by default, should be an actual price rounded to the tick size.
For example, for a $10.00 stock, the action of a bid should be like $9.97.
However, such an action changing over time cannot directly be used to
formulate an MDP. So we put forward two types of action representations:

  * •

Action measured by ticks
$a^{\textrm{ticks}}_{t}\in\mathcal{A}^{\textrm{ticks}}\subset\mathbb{Z}$,
i.e.,

$a^{\textrm{ticks}}_{t}=(\textrm{LimitPrice}_{t}-\textrm{CurrentPrice}_{t})/\textrm{TickSize}$

It is natural to represent price changes in ticks based on the current
price111Here “current price” refers to the latest executed price.. So it is a
discrete number, e.g., a bid at $9.97 for a $10.00 can be represented as
$a^{\textrm{ticks}}=-3$, if the tick size is one cent.

  * •

Action measured by percentage
$a^{\textrm{pct}}_{t}\in\mathcal{A}^{\textrm{pct}}\subset\mathbb{R}$, i.e.,

$a^{\textrm{pct}}_{t}=(\textrm{LimitPrice}_{t}/\textrm{CurrentPrice}_{t}-1)\times
100\%$

It is more general and easier to understand when measuring price movements by
percentages or basis points. So any action can has a continuous representation
by basis points, which is a continuous variable. E.g., the same bid can be
represented as $a^{\textrm{pct}}=-0.3\%$.

So apparently one can calculate $a^{\textrm{pct}}$ from any
$a^{\textrm{ticks}}$, but not from the opposite direction if without rounding.
Note that although here we only discuss the limit price, one can always append
the volume $v_{t}$ as the second term of the action.

Finally, we design the reward feedback. Like in real trading systems, we give
the agent a reward only at the end of the last period. With the order schedule
and the execution status of $t=1,\dots,T$, the reward is

$R=\mathrm{D}\cdot\bigg{[}\sum_{t=1}^{T}v^{*}_{t}p^{*}_{t}-\bigg{(}\sum_{t=1}^{T}\tilde{v}_{t}\tilde{p}_{t}+\Delta_{T}\tilde{p}_{-1}\bigg{)}\bigg{]},$
(2)

with $p_{t}^{*}$ the market price of buying $v^{*}_{t}$ in the $t$-th period,
$\tilde{p}_{t}$ the agent’s average execution price, $\Delta_{T}$ the final
remained inventory, $\tilde{p}_{-1}$ the execution price of the final market
order. The multiplier $\mathrm{D}$ denotes the trading direction,
$\mathrm{D}=1$ for buying and $\mathrm{D}=-1$ for selling. So using market
orders yields a reward of 0.

###  3.2 Policy Optimization with Hybrid Action-spaces

In this part, we put forward our policy optimization method that balances
generalization with continuous control and specialization with discrete
control. We use a general stochastic policy optimization method, PPO Schulman
et al. (2017), as our base learner.

There are two basic types of action-spaces, continuous and discrete, and the
corresponding stochastic policies, as shown in Figure 1(a)(b). For the
continuous action space $\mathcal{A}^{\textrm{pct}}$ to choose percentage
changes in price, one can use a common parametric Gaussian policy with its
mean $\mu$ and scale $\sigma$ learned with a policy neural network, i.e.,

$a^{\textrm{pct}}_{t}|s_{t}\sim\mathcal{N}(\mu(s_{t}),\sigma(s_{t})).$

For the discrete action space $\mathcal{A}^{\textrm{ticks}}$, a common choice
is the Softmax policy, which randomly choose actions from a multinoulli
categorical distribution

$a^{\textrm{ticks}}_{t}|s_{t}\sim\textrm{Categorical}\big{\\{}\textrm{Softmax}(\bm{f}_{1:m}(s_{t}))\big{\\}},$

where $\bm{f}_{1:m}(s_{t})$ denotes the policy network’s output logits for $m$
candidate actions, and $\textrm{Categorical}\\{\cdot\\}$ is the categorical
(multinoulli) distribution.

Although straight-forward, these policies are hard to optimize because the
environment’s state dynamics are too complex, the rewards are noisy and
sparse, and the action space is too large especially for the discrete policy.
To further restrict the search space, we propose a new family of policies: the
_discretized Gaussian_ policy, and its simplified version, the _Gaussian-and-
Softmax_ policy.

####  3.2.1 Discretized Gaussian Policies

Consider that we have a policy network that outputs the mean $\mu$ and scale
$\sigma$ of some Gaussian distribution $\mathcal{N}$, and a increasing series
of predefined locations (candidate actions) $a^{*}_{1},\dots,a^{*}_{m}$ for
discretization.

Let $a^{*}_{0}=-\infty$ and $a^{*}_{m+1}=+\infty$, for $1\leq k\leq m$

$\displaystyle d_{k}$
$\displaystyle=\textrm{P}(\frac{a^{*}_{k-1}+a^{*}_{k}}{2}<a<\frac{a^{*}_{k}+a^{*}_{k+1}}{2})$
(3)
$\displaystyle=\Phi(\frac{a^{*}_{k}+a^{*}_{k+1}-2\mu}{2\sigma})-\Phi(\frac{a^{*}_{k-1}+a^{*}_{k}-2\mu}{2\sigma})$
(4)

where
$\Phi(x)=\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{x}\exp\\{-\frac{u^{2}}{2}\\}du$
is the CDF of the standard normal distribution.

So now we can build a discretized Gaussian policy by using a categorical
distribution to choose the percentage action, i.e.,

$a^{\textrm{pct}}_{t}|s_{t}\sim\textrm{Categorical}\big{\\{}d_{1},\cdots,d_{m}\big{\\}}$
(5)

where the $a^{*}$s for calculating $d$s can be determined by the actual tick
actions. In this case, such a policy can take the place of the original
continuous Gaussian polices.

However, Eq.(4) involves an integral which is intractable for a feed-forward
neural network. Therefore, we use uniform samples from
$U(\frac{a_{k-1}+a_{k}}{2},\frac{a_{k}+a_{k+1}}{2})$ and calculate the
averaged Gaussian density so that the output action probability can be
differentiable. I.e.,

$\displaystyle\hat{d}_{k}=\frac{a_{k+1}-a_{k-1}}{2}$
$\displaystyle\hat{\mathbb{E}}\bigg{[}\frac{1}{\sigma\sqrt{2\pi}}\exp\left\\{-\frac{(a-\mu)^{2}}{2\sigma^{2}}\right\\}$
$\displaystyle\big{|}\,a\sim
U(\frac{a_{k-1}+a_{k}}{2},\frac{a_{k}+a_{k+1}}{2})\bigg{]}$ (6)

In this way, any Gaussian policy can be transformed into a trainable discrete
policy, as demonstrated in Figure 1 (c). For readability, for a continuous
policy $\pi(\cdot|s)=\mathcal{N}(\mu(s),\sigma(s))$, we write the discretized
policy as $\textrm{Discretize}(\mathcal{N}(\mu(s),\sigma(s)))$.

####  3.2.2 Gaussian-and-Softmax Policy

A straight-forward simplification of the discretized Gaussian policy is to
replace the sampling of Eq.(6) with a deterministic density function. So we
can directly use a series of unnormalized logits $\bm{l}=(l_{1},\dots,l_{m})$
to get a Softmax policy,

$\displaystyle
a^{\textrm{pct}}_{t}|s_{t}\sim\textrm{Categorical}\big{\\{}\textrm{Softmax}(\bm{l})\big{\\}}.$
(7)

where $l_{k}=-\frac{(a_{k}-\mu)^{2}}{2\sigma^{2}}$. We write
$\textrm{GSoftmax}(\mathcal{N}(\mu(s),\sigma(s)))$ to denote such a Gaussian-
and-Softmax policy.

####  3.2.3 HALOP with Two Stages

From the derivation, now we can get a discretized policy from any continuous
Gaussian policy. It motivates us to feel free to design a two-stage
generalization-first-specialization-second method under the same policy family
and network structure.

Stage 1: learn the most generalized knowledge without distinguishing different
types of stocks. Therefore, we have the following configuration. Firstly, the
input states in this stage only contains the standardized public state
$s^{\textrm{pub}}_{t}$ (by preprocessing all the prices and volumes to the
same scale for different stocks), without the private execution status that
reflects the characteristics of stocks. Second, the agent outputs the
discretized percentage changes of price movements according to the realizable
limit prices in order that the policy learned can be smoothly connected with
the second stage. Formally, let $p^{\textrm{c}}_{t}$ be the current price of a
stock, and $\theta_{1}$ be the parameters of the policy network in stage-1,

$\displaystyle\small\mathcal{A}_{\textrm{S1}}=$
$\displaystyle\bigg{\\{}a^{\textrm{pct}}\bigg{|}a^{\textrm{pct}}=\frac{a^{\textrm{ticks}}\cdot\textrm{TickSize}}{p^{\textrm{c}}_{t}},a^{\textrm{ticks}}\in\mathcal{A}^{\textrm{ticks}}\bigg{\\}},$
(8) $\displaystyle a^{\textrm{pct}}_{\textrm{S1},t}\sim$
$\displaystyle\pi_{\textrm{S1}}(\cdot|s^{\textrm{pub}}_{t})=\textrm{Discretize}(\mathcal{N}(\mu_{\theta_{1}}(s^{\textrm{pub}}_{t}),\sigma_{\theta_{1}}(s^{\textrm{pub}}_{t}))),$
(9) $\displaystyle
a^{\textrm{ticks}}_{\textrm{S1},t}=a^{\textrm{pct}}_{\textrm{S1},t}\cdot
p^{\textrm{c}}_{t}/\,\textrm{TickSize}.$ (10)

So the resulted output of the policy is a discrete action
$a^{\textrm{ticks}}_{\textrm{S1},t}$ in the tick action space, which
represents a general guess of the limit price. We would like to use this
action to scope a subset of actions around it, as illustrated in Figure 1(d).

Stage 2: After generalization, now for the specialization. We use the action
offered by the policy of stage 1 as a medium point and search a fine-grained
action around it in stage 2. We predefine a fixed-length window to scope the
action subset of stage 2, i.e., consider a window-size of $2K+1$, the action
subset of stage 2 is
$\mathcal{A}_{\textrm{S2}}=\big{\\{}-K,\dots,0,\dots,K\big{\\}}$. Therefore,
as illustrated in Figure 1(e), for any stock with arbitrary price, the action-
space in stage-2 is always $2K+1$ integers, encouraging the agent to focus on
localized optimization. The state in stage 2 $s_{t}$ includes the raw public
state without standardization and the private states, which contains the
absolute value of prices, volumes, and the agent’s execution status,
characterizing each single stock. For stabilized learning, we scale the prices
and volumes by taking the logarithms of the absolute values. Formally, let
$\theta_{2}$ denote the parameters of the policy network in stage-2, we have

$\displaystyle\small a_{\textrm{S2},t}\sim\pi_{\textrm{S2}}(\cdot|s_{t})$
$\displaystyle=\textrm{GSoftmax}(\mathcal{N}(\mu_{\theta_{2}}(s_{t}),\sigma_{\theta_{2}}(s_{t})))$
(11) $\displaystyle
a^{\textrm{ticks}}_{\textrm{S2},t}=a^{\textrm{ticks}}_{\textrm{S1},t}+a_{\textrm{S2},t}.$
(12)

Figure 2: The architecture of HALOP.

###  3.3 Architectures and Implementation Details

In this part, we detail the choice of neural network architecture for learning
from high-dimensional sequence data, as well as the learning details of PPO.

The overall network structure is shown in Figure 2. The input of the
representation learning network is a sequence of high-dimensional tick
snapshot data, i.e., the public states, with tens of real-valued features and
hundreds of timesteps. First, we use two stacked neural blocks for sequence
encoding, each block consists of a residual network with an 1-D convolutional
neural network (1DCNN) to capture temporal patterns and shorten the sequence
length, and a multi-head self-attention layer Vaswani et al. (2017) to learn
the correlation among time steps. Finally, we get a dense representation
$h(s)$ for the state by attentive pooling the output of the self-attention
matrix.

Second, we use an actor-critic architecture Schulman et al. (2017); Pan et al.
(2019) for the policy optimization. For both stage-1 and stage-2, we use two
fully-connected layers as the policy and value networks. In stage-1, for
generalization, the networks directly take the public
$h(s^{\textrm{pub}}_{t})$ as input. In stage-2, we concatenate the private
state to $h(s^{\textrm{pub}}_{t})$ to get a new dense vector as input. The two
stages share the same representation learning network for better sample
efficiency.

For policy optimization, we regard the order execution process of each
(trading day, stock) pair as a single _episode_ , which has $T$ steps. Then,
we group the episodes of all the stocks in each trading day as an $epoch$. For
example, if there are 250 trade days and 300 stocks, then it groups into 250
epochs and $250\times 300$ episodes. During training, we iterates between
rollouts and updates in a batch RL manner. At each round, we randomly select
one trading day from the training set, and run parallel rollout for all the
episodes in that trading day. After finished, we perform the on-policy policy
update with the collected data like the standard PPO algorithm Schulman et al.
(2017).

Method | TWAP | VWAP | ODP
---|---|---|---
Return | PnL | Std |  $t$-value | Return | PnL | Std |  $t$-value | Return | PnL | Std |  $t$-value
Market Order | 0.00 | 0.00 | 0.00 | 0.0 | 0.14 | 0.14 | 8.72 | 3.1 | 1.64 | 0.92 | 65.34 | 4.8
PPO-Gaussian | 1.41 | 1.33 | 4.61 | 58.7 | 1.38 | 1.31 | 8.94 | 29.7 | 2.20 | 1.12 | 64.55 | 6.5
PPO-Softmax | 2.12 | -1.46 | 35.88 | 11.3 | 2.36 | -1.17 | 35.47 | 12.8 | 3.34 | -0.53 | 66.73 | 9.6
HALOP Stage-1 | 3.69 | 1.46 | 17.05 | 41.5 | 3.87 | 1.68 | 17.51 | 42.4 | 4.51 | 1.45 | 63.67 | 13.6
HALOP | 4.41 | 2.98 | 6.87 | 123.2 | 4.52 | 3.12 | 9.24 | 93.8 | 5.04 | 2.46 | 63.42 | 15.3
Table 1: Main results.

##  4 Experiments

Our experiments are designed to investigate the following questions:

Q1. Can HALOP beat the market?

Q2. Can HALOP improve upon other volume-oriented order splitting strategies?
Is HALOP stable and robust?

Q3. What benefits can we get from HALOP compared with vanilla continuous or
discrete control?

Q4. Can HALOP address the trade-off between generalization and specialization
better than other methods?

###  4.1 Experiment Setup

####  4.1.1 The Simulation Environment

We conduct experiments in simulation environments which is built over
historical high-frequency transaction data of stocks of the CSI 300 index in
the China A-shares. Specifically, we test all the methods in the buying
direction and directly use the raw market information as the agent’s state,
which consists of the top 5 bid/ask prices and volumes, and the current (last)
price. As described in the previous section, we use the data from January 2010
to December 2019 to build the training episodes, and use the data from January
2020 to June 2020 to build the testing episodes.

At each episode (i.e., for one stock on one trading day), our agent’s
execution mission starts at 10:00:00 (half an hour after opening) and ends at
14:30:00 (half an hour before close) which is an 180-minute trading hour. The
time period is then split into $T=90$ intervals evenly. At each step, the
agents observes the market information and send an order to the market. We
simulate a three second communication delay so as to make the experiment close
to reality. For example, an order sent at 10:02:00 can be executed after
10:02:03. Moreover, to avoid potential large market impact, our environment
does not execute any large order whose size is larger than $1/10$ of the
agent’s total inventory of the day. So on receiving an order with size
$v_{t}$, the environment executes no more than $\min(v_{t},1/10)$ and cancels
the rest of it, if there is.

####  4.1.2 Evaluation Metrics

Our evaluation metrics are designed around two core concepts in quantitative
finance: the excess return and the risk.

We use the average excess return over TWAP-with-market-order as the Return
metric. That is, the TWAP strategy using market orders always get a return of
0, and other methods get positive returns if it beats TWAP. The returns are
displayed in basis points (bps, 1 bp equals to $0.01\%$).

Next, we compute the standard deviation of excess returns as a risk indicator,
denoted as Std. Smaller standard deviation of returns means that the algorithm
is more stable throughout different time and stocks.

With the computed Return and Std metrics, we calculate the $t$-value for
paired student’s $t$-test, to show the significance of the excess return. For
any method, given the metric Return and Std, the $t$-value can be computed by

$t\textrm{-value}=\frac{\textrm{Return}}{\textrm{Std}/\sqrt{n-1}},$ (13)

where $n$ is the number of episodes. A larger $t$-value indicates more
significant results. When $n$ is large enough, $t>3.29$ means that there is a
improvement over the benchmark in confidence level 99.9%.

Finally, we need a more realistic metric that takes both return and safety
constaints into consideration, which is our PnL (or Profit and Loss) for order
execution. The profit is the excess return. For the loss, we take a regulation
constraint into consideration: the cancellation rate for institutional
investors cannot be greater than 50%. So when the cancellation rate exceeds
50%, we penalize the agent by -5 bps. So

$\textrm{PnL}=\textrm{Return}-\sum_{i}^{n}5\cdot\mathbb{I}\big{[}\Delta_{T}+\sum_{t=1}^{T}v_{t}>2\big{]},$
(14)

with $v_{t}$ defined in Eq.(1) and $\Delta_{T}$ the size of market order at
the end of an episode.

###  4.2 Compared Methods

There are two categories of methods we compared: the price- and volume-
oriented optimization methods. For the first part, we compare the following
methods:

  * •

Market Order which place market orders at all steps, meant to be executed as
quickly as possible.

  * •

PPO-Gaussian which is PPO Schulman et al. (2017) with a Gaussian policy
(Figure 1(a)).

  * •

PPO-Softmax which is PPO Schulman et al. (2017) with a Softmax policy (Figure
1(b)).

  * •

HALOP is our proposed method with two stages.

  * •

HALOP Stage-1 is HALOP with only one stage (using the discretized Gaussian
policy). It takes both public and private states as inputs. We test it to see
whether it can improve upon the ordinary Gaussian policy.

For a fair comparison, all the RL algorithms are tested with the same network
architecture (shown in Figure 2). The number of actions in HALOP stage-2 is
set to $7=2\times 3+1$.

To see how the agents perform under different volume schedule, we test the
following volume-oriented methods:

  * •

TWAP which evenly splits the order to $T$ pieces and execute the same amount
of shares at each step.

  * •

VWAP which sends orders in proportion to an estimated volume ratio of each
time period of the day, which is estimated from the previous 21 trading days.

  * •

OPD Fang et al. (2021) which leverages RL with policy distillation to
determine the size of each sub-order.

###  4.3 Main Experimental Results

Table 1 reports the overall performance.

First, to answer question Q1, it is clear that our HALOP method can indeed
outperform the market orders with a large margin. When using the basic TWAP,
HALOP beats Market Order with 4.41 bps and a PnL of 2.98, with a $t$-value of
123.2 which proves the significance of improvements.

Second, to answer Q2, we would like to see how HALOP performs under other
volume allocation methods. We firstly see that the superiority of HALOP
againest Market Order still holds with either TWAP, VWAP or ODP. Also,
surprisingly, we find that when using HALOP, the Std of ODP decreased
comparing to the original ODP with market order, indicating that HALOP can
even stabilize the order execution process in certain circumstances.

Third, to answer Q3, we need to compare HALOP with ordinary PPOs as well as
the lite version of HALOP with only stage-1. We find that though PPO-Gaussian
and PPO-Softmax all yields positive excess returns, the improvements are not
as large as that of HALOP. Moreover, while PPO-Softmax gets a higher excess-
return than PPO-Gaussian, it has negative PnL and very large Std. It indicates
that although the Softmax policy is capable of finding the best action, it is
unstable and unreliable. We also see that HALOP Stage-1 yields significantly
better result than the ordinary PPOs, showing its ability of taking the
advantages of both Gaussian policy and discrete control. Finally, the
comparison between HALOP and HALOP Stage-1 shows that the specialization in
stage-2 indeed helps improve both excess return and stability.

###  4.4 Grouping Study

We further design a grouping study to see whether HALOP address the trade-off
between generalization and specialization well. We group the test episodes by
stock close prices: “low-priced” for price less than 10.00 CNY, “medium-
priced” for price between 10.00-50.00 CNY, and “high-priced” for larger than
50.00 CNY.

Figure 3: Results of grouping study

As shown in Figure 3, we find PPO-Gaussian more in favor of high-priced stocks
than low-priced stocks. Next, we observe that HALOP Stage-1 has a strong
ability of generalization in all the groups: the excess returns are all over
3.5. Further, the two-stage HALOP further improves the results especially for
low-priced stocks, which indicates that the fine-grained discrete control in
stage-2 indeed achieves better specialization for different types of stocks.
Specifically, for penny stocks which are known to have less liquidity, the
discretization in HALOP Stage-2 is shown to be more beneficial.

##  5 Conclusion

This paper focuses on optimal execution with limit orders by setting better
limit prices. We propose a novel reinforcement solution, Hybrid Action-Space
Order Placement, to address the generalization-specialization trade-off by
combining both advantages of continuous and discrete control. The proposed
method has two stages: in stage-1, the agent selects a target percentage
change in price with a discretized Gaussian policy and scopes a subset of
discrete actions thereby; in stage-2, the agent uses a fine-grained Gaussian-
and-Softmax policy to select a tick-based action. Extensive simulation-based
experiments show that our method can beat the market and improve upon other
volume-oriented order splitting strategies. It shows stable and robust
improvements comparing with vanilla continuous or discrete control methods.
Specifically, the discretized Gaussian policy of stage-1 helps improve the
generalization ability, and meanwhile in stage-2 the agent specializes in
policy for different types of stocks especially for low-priced stocks. In a
broader sense, we believe that HALOP offers a reference of applying RL in
real-world applications where the trade-off between generalization and
specialization should be taken into consideration.

## References

  * Almgren and Chriss [2001] Robert Almgren and Neil Chriss.  Optimal execution of portfolio transactions.  Journal of Risk, 3:5–40, 2001.
  * Bertsimas and Lo [1998] Dimitris Bertsimas and Andrew W Lo.  Optimal control of execution costs.  Journal of Financial Markets, 1(1):1–50, 1998.
  * Cartea et al. [2015] Álvaro Cartea, Sebastian Jaimungal, and José Penalva.  Algorithmic and high-frequency trading.  Cambridge University Press, 2015.
  * Fang et al. [2021] Yuchen Fang, Kan Ren, Weiqing Liu, Dong Zhou, Weinan Zhang, Jiang Bian, Yong Yu, and Tie-Yan Liu.  Universal trading for order execution with oracle policy distillation.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 35, pages 107–115, 2021.
  * Hendricks and Wilcox [2014] Dieter Hendricks and Diane Wilcox.  A reinforcement learning extension to the almgren-chriss framework for optimal trade execution.  In 2014 IEEE Conference on Computational Intelligence for Financial Engineering & Economics (CIFEr), pages 457–464. IEEE, 2014.
  * Huberman and Stanzl [2000] G. Huberman and W. Stanzl.  Optimal liquidity trading.  Yale School of Management Working Papers, 9(2):165–200, 2000.
  * Kakade et al. [2004] Sham M Kakade, Michael Kearns, Yishay Mansour, and Luis E Ortiz.  Competitive algorithms for vwap and limit order trading.  In Proceedings of the 5th ACM conference on Electronic commerce, pages 189–198, 2004.
  * Lin and Beling [2019] Siyu Lin and AP Beling.  Optimal liquidation with deep reinforcement learning.  In 33rd Conference on Neural Information Processing Systems (NeurIPS 2019) Deep Reinforcement Learning Workshop, Vancouver, Canada, 2019\.
  * Lin and Beling [2020] Siyu Lin and Peter A Beling.  An end-to-end optimal trade execution framework based on proximal policy optimization.  In IJCAI, pages 4548–4554, 2020.
  * Mnih et al. [2013] Volodymyr Mnih, Koray Kavukcuoglu, David Silver, Alex Graves, Ioannis Antonoglou, Daan Wierstra, and Martin Riedmiller.  Playing atari with deep reinforcement learning.  2013\.  cite arxiv:1312.5602Comment: NIPS Deep Learning Workshop 2013.
  * Nevmyvaka et al. [2006] Yuriy Nevmyvaka, Yi Feng, and Michael Kearns.  Reinforcement learning for optimized trade execution.  In Proceedings of the 23rd international conference on Machine learning, pages 673–680, 2006.
  * Ning et al. [2018] Brian Ning, Franco Ho Ting Lin, and Sebastian Jaimungal.  Double deep q-learning for optimal execution.  arXiv preprint arXiv:1812.06600, 2018.
  * Pan et al. [2019] Feiyang Pan, Qingpeng Cai, An-Xiang Zeng, Chun-Xiang Pan, Qing Da, Hualin He, Qing He, and Pingzhong Tang.  Policy optimization with model-based explorations.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 33, pages 4675–4682, 2019.
  * Schulman et al. [2017] John Schulman, Filip Wolski, Prafulla Dhariwal, Alec Radford, and Oleg Klimov.  Proximal policy optimization algorithms.  CoRR, abs/1707.06347, 2017.
  * Vaswani et al. [2017] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin.  Attention is all you need.  In Advances in neural information processing systems, pages 5998–6008, 2017.
  * Watkins and Dayan [1992] Cjch Watkins and P. Dayan.  Q-learning.  In Machine Learning, 1992.


$\displaystyle\leq\frac{11}{10}\cdot\frac{1}{10}+\frac{1}{10}=\frac{21}{100}$

(iii) By Lemma C.5, we have $\boldsymbol{\mathit{s}}^{+},s_{gap}^{+}>0$ and

$\tilde{\eta}(\boldsymbol{\mathit{s}}^{+},s_{gap}^{+})\leq\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap}^{\prime})^{2}+\frac{1}{20}\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap}^{\prime})\leq\left(\frac{21}{100}\right)^{2}+\frac{1}{20}\cdot\frac{21}{100}<\frac{1}{10}$

∎

$\boldsymbol{\mathit{x}}=\mathtt{InteriorPoint}(A,\boldsymbol{\mathit{b}},\boldsymbol{\mathit{c}},\boldsymbol{\mathit{y}}^{0},\epsilon)$
• Compute
$(\boldsymbol{\mathit{y}}^{C},z^{C})=\mathtt{FindCentralPath}(A,\boldsymbol{\mathit{b}},\boldsymbol{\mathit{c}},\boldsymbol{\mathit{y}}^{0})$
and $\left[\begin{matrix}\boldsymbol{\mathit{s}}^{C}\\\
s_{gap}^{C}\end{matrix}\right]=\left[\begin{matrix}\boldsymbol{\mathit{c}}-A^{T}\boldsymbol{\mathit{y}}^{C}\\\
\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}^{C}-z^{C}\end{matrix}\right]$
• Set
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap},z):=(\boldsymbol{\mathit{y}}^{C},\boldsymbol{\mathit{s}}^{C},s_{gap}^{C},z^{C})$
• While $s_{gap}>\frac{\epsilon}{3}$: – Set
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap},z):=\mathtt{Shift}(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap},z)$
• Compute
$\boldsymbol{\mathit{v}}=\mathtt{Solve}(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T},\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m},\epsilon_{4})$
$\displaystyle\text{where }\tilde{A}$
$\displaystyle=\left[\begin{matrix}A&-\boldsymbol{\mathit{b}}{\mbox{\boldmath$1$}}_{m}^{T}\end{matrix}\right]$
$\displaystyle\text{and
}\epsilon_{4}=\min\left(1,\frac{s_{min}}{TU}\cdot\frac{m^{1/2}}{n}\right)$
$\displaystyle\text{and }\tilde{\boldsymbol{\mathit{s}}}$
$\displaystyle=\left[\begin{matrix}\boldsymbol{\mathit{s}}\\\
s_{gap}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$ and $s_{min}$ is the
smallest entry of $\tilde{\boldsymbol{\mathit{s}}}$ • Return
$\boldsymbol{\mathit{x}}=\frac{\boldsymbol{\mathit{x}}^{\prime}}{mx_{gap}^{\prime}}$
where $\left[\begin{matrix}\boldsymbol{\mathit{x}}^{\prime}\\\
x^{\prime}_{gap}\end{matrix}\right]=\left[\begin{matrix}S^{-1}{\mbox{\boldmath$1$}}_{m}-S^{-2}A^{T}\boldsymbol{\mathit{v}}\\\
s_{gap}^{-1}+s_{gap}^{-2}\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{v}}\end{matrix}\right]$

Figure 4: Dual path-following interior-point algorithm using an approximate
solver

We now present the complete path-following $\mathtt{InteriorPoint}$ algorithm,
implemented using an approximate solver, in Figure 4. For now we postpone
describing the $\mathtt{FindCentralPath}$ subroutine, which gives an initial
point near the central path. In particular, it produces a $z^{C}<z^{*}$ and
$(\boldsymbol{\mathit{y}}^{C},\boldsymbol{\mathit{s}}^{C},s_{gap}^{C})\in\Omega^{D}_{\boldsymbol{\mathit{b}},z^{C}}$
satisfying
$\tilde{\eta}(\boldsymbol{\mathit{s}}^{C},s_{gap}^{C})\leq\frac{1}{10}$. Once
we have this initial central path point, Lemma C.8 tells us that after each
call to $\mathtt{Shift}$ we have a new value $z<z^{*}$ and new central path
point
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap})\in\mathring{\Omega}^{D}_{\boldsymbol{\mathit{b}},z}$
that satisfies
$\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})\leq\frac{1}{10}$.

Later we will analyze the number of calls to $\mathtt{Shift}$ before the
algortihm terminates. First, let us confirm that the algorithm returns the
correct output:

######  Lemma C.9.

The output of
$\boldsymbol{\mathit{x}}=\mathtt{InteriorPoint}(A,\boldsymbol{\mathit{b}},\boldsymbol{\mathit{c}},\boldsymbol{\mathit{y}}^{0},\epsilon)$
satisfies

  * (i)

$\boldsymbol{\mathit{x}}>0$

  * (ii)

$\left\|A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}\right\|\leq\frac{\epsilon}{12\sqrt{2}\cdot
Tn^{1/2}}$

  * (iii)

$\boldsymbol{\mathit{c}}^{T}\boldsymbol{\mathit{x}}<z^{*}+\epsilon$

###### Proof.

(i) To assist in our proof, let us define
$\tilde{\boldsymbol{\mathit{x}}}^{\prime}=\left[\begin{matrix}\boldsymbol{\mathit{x}}^{\prime}\\\
x^{\prime}_{gap}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$ and note that
$\tilde{\boldsymbol{\mathit{x}}}^{\prime}=\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}-\tilde{S}^{-2}\tilde{A}^{T}\boldsymbol{\mathit{v}}$.
We have

$\displaystyle\left\|\tilde{S}\tilde{\boldsymbol{\mathit{x}}}^{\prime}-{\mbox{\boldmath$1$}}_{2m}\right\|$
$\displaystyle=\left\|\tilde{S}^{-1}\tilde{A}^{T}\boldsymbol{\mathit{v}}\right\|$
$\displaystyle=\left\|\boldsymbol{\mathit{v}}\right\|_{\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}}$
$\displaystyle\leq\left\|(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}}+\left\|\boldsymbol{\mathit{v}}-(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}}$
$\displaystyle\leq(1+\epsilon_{4})\left\|(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}}\qquad\text{(by
guarantee of {Solve})}$
$\displaystyle=(1+\epsilon_{4})\cdot\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})\leq
2\cdot\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})\leq
2\cdot\frac{1}{10}\leq\frac{1}{5}$ (10)

Since $\tilde{\boldsymbol{\mathit{s}}}$ is positive, we conclude that
$\tilde{\boldsymbol{\mathit{x}}}^{\prime}$ must also be positive, and so must
be $\boldsymbol{\mathit{x}}$.

(ii) We have

$\displaystyle\left\|A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}\right\|$
$\displaystyle=\frac{1}{mx_{gap}^{\prime}}\left\|A\boldsymbol{\mathit{x}}^{\prime}-mx_{gap}^{\prime}\boldsymbol{\mathit{b}}\right\|$
$\displaystyle=\frac{1}{mx_{gap}^{\prime}}\left\|\tilde{A}\tilde{\boldsymbol{\mathit{x}}}^{\prime}\right\|$
$\displaystyle=\frac{1}{mx_{gap}^{\prime}}\left\|\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}-\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}\boldsymbol{\mathit{v}}\right\|$
Observe that the largest eigenvalue of the matrix
$\tilde{A}\tilde{A}^{T}=AA^{T}+m\boldsymbol{\mathit{b}}\boldsymbol{\mathit{b}}^{T}$
is less than the trace, which is at most $2nmU^{2}$. Thus, the largest
eigenvalue of $\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}$ is at most
$2nmU^{2}s_{min}^{-2}$. So we proceed
$\displaystyle\leq\frac{(2nmU^{2}s_{min}^{-2})^{1/2}}{mx_{gap}^{\prime}}\left\|\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}-\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}\boldsymbol{\mathit{v}}\right\|_{(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}}$
$\displaystyle=\frac{(2n)^{1/2}U}{m^{1/2}x_{gap}^{\prime}s_{min}}\left\|(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}-\boldsymbol{\mathit{v}}\right\|_{\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}}$
$\displaystyle\leq\frac{(2n)^{1/2}U}{m^{1/2}x_{gap}^{\prime}s_{min}}\cdot\epsilon_{4}\left\|(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}}\qquad\text{(by
guarantee of {Solve})}$
$\displaystyle=\frac{(2n)^{1/2}U}{m^{1/2}x_{gap}^{\prime}s_{min}}\cdot\epsilon_{4}\cdot\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})$
$\displaystyle\leq\frac{(2n)^{1/2}U}{m^{1/2}x_{gap}^{\prime}s_{min}}\cdot\epsilon_{4}\cdot\frac{1}{10}$
$\displaystyle\leq\frac{(2n)^{1/2}U}{m^{1/2}x_{gap}^{\prime}s_{min}}\cdot\frac{m^{1/2}s_{min}}{nTU}\cdot\frac{1}{10}$
$\displaystyle=\frac{1}{5\sqrt{2}\cdot
Tn^{1/2}}\cdot\frac{1}{x_{gap}^{\prime}}$
$\displaystyle\leq\frac{1}{5\sqrt{2}\cdot Tn^{1/2}}\cdot\frac{5}{4}\cdot
s_{gap}\qquad\text{(from equation \ref{eq:sxprime}, we know
$s_{gap}x_{gap}^{\prime}\geq\frac{4}{5}$)}$
$\displaystyle\leq\frac{1}{5\sqrt{2}\cdot
Tn^{1/2}}\cdot\frac{5}{4}\cdot\frac{\epsilon}{3}$
$\displaystyle\leq\frac{\epsilon}{12\sqrt{2}\cdot Tn^{1/2}}$ (11)

(iii) We have

$\displaystyle\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}$
$\displaystyle=\frac{\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}^{\prime}}{mx^{\prime}_{gap}}$
$\displaystyle\geq\frac{5\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}^{\prime}}{6m}\cdot
s_{gap}\qquad\text{(from equation \ref{eq:sxprime}, we know
$s_{gap}x_{gap}^{\prime}\leq\frac{6}{5}$)}$
$\displaystyle=\frac{5\left\|S\boldsymbol{\mathit{x}}^{\prime}\right\|_{1}}{6m}\cdot
s_{gap}$
$\displaystyle\geq\frac{5}{6m}\left(\left\|{\mbox{\boldmath$1$}}_{m}\right\|_{1}-\left\|S\boldsymbol{\mathit{x}}^{\prime}-{\mbox{\boldmath$1$}}_{m}\right\|_{1}\right)\cdot
s_{gap}$
$\displaystyle=\frac{5}{6m}\left(m-\left\|S\boldsymbol{\mathit{x}}^{\prime}-{\mbox{\boldmath$1$}}_{m}\right\|_{1}\right)\cdot
s_{gap}$
$\displaystyle\geq\frac{5}{6m}\left(m-\sqrt{m}\left\|S\boldsymbol{\mathit{x}}^{\prime}-{\mbox{\boldmath$1$}}_{m}\right\|\right)\cdot
s_{gap}$ $\displaystyle\geq\frac{5}{6m}\left(m-\frac{1}{5}\sqrt{m}\right)\cdot
s_{gap}\qquad\text{(by equation \ref{eq:sxprime})}$
$\displaystyle\geq\frac{5}{6m}\cdot\frac{4}{5}\cdot m\cdot s_{gap}$
$\displaystyle=\frac{2}{3}\cdot s_{gap}$ (12)

$\displaystyle\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}$
$\displaystyle=\frac{\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}^{\prime}}{mx^{\prime}_{gap}}$
$\displaystyle\leq\frac{5\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}^{\prime}}{4m}\cdot
s_{gap}\qquad\text{(from equation \ref{eq:sxprime}, we know
$s_{gap}x_{gap}^{\prime}\geq\frac{4}{5}$)}$
$\displaystyle=\frac{5\left\|S\boldsymbol{\mathit{x}}^{\prime}\right\|_{1}}{4m}\cdot
s_{gap}$
$\displaystyle\leq\frac{5}{4m}\left(\left\|S\boldsymbol{\mathit{x}}^{\prime}-{\mbox{\boldmath$1$}}_{m}\right\|_{1}+\left\|{\mbox{\boldmath$1$}}_{m}\right\|_{1}\right)\cdot
s_{gap}$
$\displaystyle=\frac{5}{4m}\left(\left\|S\boldsymbol{\mathit{x}}^{\prime}-{\mbox{\boldmath$1$}}_{m}\right\|_{1}+m\right)\cdot
s_{gap}$
$\displaystyle\leq\frac{5}{4m}\left(\sqrt{m}\left\|S\boldsymbol{\mathit{x}}^{\prime}-{\mbox{\boldmath$1$}}_{m}\right\|+m\right)\cdot
s_{gap}$ $\displaystyle\leq\frac{5}{4m}\left(\frac{1}{5}\sqrt{m}+m\right)\cdot
s_{gap}\qquad\text{(by equation \ref{eq:sxprime})}$
$\displaystyle\leq\frac{5}{4m}\cdot\frac{6}{5}\cdot m\cdot s_{gap}$
$\displaystyle=\frac{3}{2}\cdot s_{gap}$ (13)

We then have

$\displaystyle\boldsymbol{\mathit{c}}^{T}\boldsymbol{\mathit{x}}-z^{*}$
$\displaystyle>\boldsymbol{\mathit{c}}^{T}\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}$
$\displaystyle=(\boldsymbol{\mathit{c}}^{T}-\boldsymbol{\mathit{y}}^{T}A)\boldsymbol{\mathit{x}}+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})$
$\displaystyle=\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})$
$\displaystyle\geq\frac{2}{3}\cdot
s_{gap}+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})\qquad\text{(by
Equation \ref{eq:sxupperbound})}$
$\displaystyle\leq\dots\frac{2}{3}\cdot\epsilon+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})$
$\displaystyle\leq\frac{2}{3}\cdot\epsilon+\left\|\boldsymbol{\mathit{y}}\right\|\left\|A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}\right\|$
$\displaystyle\leq\frac{2}{3}\cdot\epsilon+Tn^{1/2}\left\|A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}\right\|$
$\displaystyle\leq\frac{2}{3}\cdot\epsilon+\frac{1}{12\sqrt{2}}\cdot\epsilon\qquad\text{(by
Lemma \ref{lem:intptoutput}(ii))}$ $\displaystyle<\epsilon$
$\displaystyle\boldsymbol{\mathit{c}}^{T}\boldsymbol{\mathit{x}}-z^{*}$
$\displaystyle<\boldsymbol{\mathit{c}}^{T}\boldsymbol{\mathit{x}}-z$
$\displaystyle=(\boldsymbol{\mathit{c}}^{T}-\boldsymbol{\mathit{y}}^{T}A)\boldsymbol{\mathit{x}}+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})+\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-z$
$\displaystyle=\boldsymbol{\mathit{s}}^{T}\boldsymbol{\mathit{x}}+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})+s_{gap}$
$\displaystyle\leq\frac{5}{2}\cdot
s_{gap}+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})\qquad\text{(by
Equation \ref{eq:sxupperbound})}$
$\displaystyle\leq\frac{5}{6}\cdot\epsilon+\boldsymbol{\mathit{y}}^{T}(A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}})$
$\displaystyle\leq\frac{5}{6}\cdot\epsilon+\left\|\boldsymbol{\mathit{y}}\right\|\left\|A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}\right\|$
$\displaystyle\leq\frac{5}{6}\cdot\epsilon+Tn^{1/2}\left\|A\boldsymbol{\mathit{x}}-\boldsymbol{\mathit{b}}\right\|$
$\displaystyle\leq\frac{5}{6}\cdot\epsilon+\frac{1}{12\sqrt{2}}\cdot\epsilon\qquad\text{(by
Lemma \ref{lem:intptoutput}(ii))}$ $\displaystyle<\epsilon$

∎

Next, we analyze the number of $\mathtt{Shift}$ iterations until the algorithm
terminates. We can measure the progress of the algorithm with the potential
function $B(z)$:

$B(z)=\sum_{j=1}^{m}\log\accentset{*}{\boldsymbol{\mathit{s}}}_{j}+m\log\accentset{*}{s}_{gap}\qquad\text{where
$(\accentset{*}{\boldsymbol{\mathit{y}}},\accentset{*}{\boldsymbol{\mathit{s}}},\accentset{*}{s}_{gap})$
is the analytic center of $\Omega^{D}_{\boldsymbol{\mathit{b}},z}$}$

Soon, we will show how a decrease in $B(z)$ implies that $s_{gap}$ is
decreasing and thus the algorithm is making progress. Let us first show that
the value of $B(z)$ decreases by $\Omega(\sqrt{m})$ after each iteration.

######  Lemma C.10 (compare [Ye97, Lem 4.6]).

Given
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap})\in\mathring{\Omega}^{D}_{\boldsymbol{\mathit{b}},z}$
satisfiying $\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})<\frac{1}{10}$, let
$(\boldsymbol{\mathit{y}}^{+},z^{+})=\mathtt{Shift}(\boldsymbol{\mathit{y}},z)$.
Then $B(z^{+})\leq B(z)-\Theta(\sqrt{m})$.

###### Proof.

Let
$(\accentset{*}{\boldsymbol{\mathit{y}}},\accentset{*}{\boldsymbol{\mathit{s}}},\accentset{*}{s}_{gap})$
and
$(\accentset{*}{\boldsymbol{\mathit{y}}}^{+},\accentset{*}{\boldsymbol{\mathit{s}}}^{+},\accentset{*}{s}_{gap}^{+})$
respectively be the analytic centers of
$\Omega^{D}_{\boldsymbol{\mathit{b}},z}$ and
$\Omega^{D}_{\boldsymbol{\mathit{b}},z^{+}}$.

Following Lemma C.7, we define
$\boldsymbol{\mathit{x}}=\frac{\accentset{*}{s}_{gap}}{m}\accentset{*}{S}^{-1}{\mbox{\boldmath$1$}}_{m}$
that satisfies $A\boldsymbol{\mathit{x}}=\boldsymbol{\mathit{b}}$. We have

$\displaystyle e^{\frac{B(z^{+})-B(z)}{2m}}$
$\displaystyle=\sqrt{\left(\prod_{j=1}^{m}\frac{\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}}{\accentset{*}{\boldsymbol{\mathit{s}}}_{j}}\right)^{\frac{1}{m}}\cdot\frac{\accentset{*}{s}_{gap}^{+}}{\accentset{*}{s}_{gap}}}$
$\displaystyle\leq\frac{1}{2m}\sum_{j=1}^{m}\frac{\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}}{\accentset{*}{\boldsymbol{\mathit{s}}}_{j}}+\frac{1}{2}\frac{\accentset{*}{s}_{gap}^{+}}{\accentset{*}{s}_{gap}}$
$\displaystyle\leq
1+\frac{1}{2m}\sum_{j=1}^{m}\frac{\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}-\accentset{*}{\boldsymbol{\mathit{s}}}_{j}}{\accentset{*}{\boldsymbol{\mathit{s}}}_{j}}+\frac{1}{2}\frac{\accentset{*}{s}_{gap}^{+}-\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}}$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\cdot\frac{\accentset{*}{s}_{gap}}{m}(\accentset{*}{\boldsymbol{\mathit{s}}}^{+}-\accentset{*}{\boldsymbol{\mathit{s}}})\accentset{*}{S}^{-1}{\mbox{\boldmath$1$}}_{m}+\frac{1}{2\accentset{*}{s}_{gap}}(\accentset{*}{s}_{gap}^{+}-\accentset{*}{s}_{gap})$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\left((\accentset{*}{\boldsymbol{\mathit{s}}}^{+}-\accentset{*}{\boldsymbol{\mathit{s}}})^{T}\boldsymbol{\mathit{x}}+(\accentset{*}{s}_{gap}^{+}-\accentset{*}{s}_{gap})\right)$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\left(\left((\boldsymbol{\mathit{c}}-A^{T}\accentset{*}{\boldsymbol{\mathit{y}}}^{+})-(\boldsymbol{\mathit{c}}-A^{T}\accentset{*}{\boldsymbol{\mathit{y}}})\right)^{T}\boldsymbol{\mathit{x}}+(\accentset{*}{s}_{gap}^{+}-\accentset{*}{s}_{gap})\right)$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\left((\accentset{*}{\boldsymbol{\mathit{y}}}-\accentset{*}{\boldsymbol{\mathit{y}}}^{+})^{T}A\boldsymbol{\mathit{x}}+(\accentset{*}{s}_{gap}^{+}-\accentset{*}{s}_{gap})\right)$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\left((\accentset{*}{\boldsymbol{\mathit{y}}}-\accentset{*}{\boldsymbol{\mathit{y}}}^{+})^{T}\boldsymbol{\mathit{b}}+(\accentset{*}{s}_{gap}^{+}-\accentset{*}{s}_{gap})\right)\qquad\text{(by
Lemma \ref{lem:constructx})}$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\left((\boldsymbol{\mathit{b}}^{T}\accentset{*}{\boldsymbol{\mathit{y}}}-\accentset{*}{s}_{gap})-(\boldsymbol{\mathit{b}}^{T}\accentset{*}{\boldsymbol{\mathit{y}}}^{+}-\accentset{*}{s}_{gap}^{+})\right)$
$\displaystyle=1+\frac{1}{2\accentset{*}{s}_{gap}}\left(z-z^{+}\right)$ (14)
$\displaystyle=1-\frac{s_{gap}}{20\sqrt{m}\cdot\accentset{*}{s}_{gap}}$
$\displaystyle\leq
1-\frac{9}{200\sqrt{m}}\qquad\text{($\frac{\accentset{*}{s}_{gap}}{s_{gap}}\leq\frac{10}{9}$
by Lemma \ref{lem:nearcenter})}$ $\displaystyle\leq
e^{-\frac{9}{200\sqrt{m}}}$

We conclude that $B(z^{+})-B(z)\leq-\frac{9}{100}\sqrt{m}$. ∎

Let us now show that a decrease in the potential function $B(z)$ implies a
decrease in the value of $s_{gap}$:

######  Lemma C.11 (compare [Ye97, Prop 4.2]).

Given
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap})\in\mathring{\Omega}^{D}_{\boldsymbol{\mathit{b}},z}$
and
$(\boldsymbol{\mathit{y}}^{+},\boldsymbol{\mathit{s}}^{+},s_{gap}^{+})\in\mathring{\Omega}^{D}_{\boldsymbol{\mathit{b}},z^{+}}$
where $z^{+}>z$ and $\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})\leq\eta$
and $\tilde{\eta}(\boldsymbol{\mathit{s}}^{+},s_{gap}^{+})\leq\eta$ for
$\eta<1$. Then

$\frac{s_{gap}}{s_{gap}^{+}}\leq(1-2\eta)\cdot\left(e^{\frac{B(z_{1})-B(z_{2})}{m}}-1\right)$

###### Proof.

Let
$(\accentset{*}{\boldsymbol{\mathit{y}}},\accentset{*}{\boldsymbol{\mathit{s}}},\accentset{*}{s}_{gap})$
and
$(\accentset{*}{\boldsymbol{\mathit{y}}}^{+},\accentset{*}{\boldsymbol{\mathit{s}}}^{+},\accentset{*}{s}_{gap}^{+})$
respectively be the analytic centers of
$\Omega^{D}_{\boldsymbol{\mathit{b}},z}$ and
$\Omega^{D}_{\boldsymbol{\mathit{b}},z^{+}}$. We define
$\boldsymbol{\mathit{x}}=\frac{\accentset{*}{s}_{gap}}{m}\accentset{*}{S}^{-1}{\mbox{\boldmath$1$}}_{m}$
and
$\boldsymbol{\mathit{x}}^{+}=\frac{\accentset{*}{s}_{gap}^{+}}{m}(\accentset{*}{S}^{+})^{-1}{\mbox{\boldmath$1$}}_{m}$,
which by Lemma C.7 satisfy
$A\boldsymbol{\mathit{x}}=\boldsymbol{\mathit{b}}=A\boldsymbol{\mathit{x}}^{+}$.

We have

$\displaystyle e^{\frac{B(z)-B(z^{+})}{m}}$
$\displaystyle=\left(\prod_{j=1}^{m}\frac{\accentset{*}{\boldsymbol{\mathit{s}}}_{j}}{\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}}\right)^{\frac{1}{m}}\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle\leq\frac{1}{m}\left(\sum_{j=1}^{m}\frac{\accentset{*}{\boldsymbol{\mathit{s}}}_{j}}{\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{m}\sum_{j=1}^{m}\frac{\accentset{*}{\boldsymbol{\mathit{s}}}_{j}-\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}}{\accentset{*}{\boldsymbol{\mathit{s}}}^{+}_{j}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}(\accentset{*}{\boldsymbol{\mathit{s}}}-\accentset{*}{\boldsymbol{\mathit{s}}}^{+})^{T}\boldsymbol{\mathit{x}}^{+}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}\left((\boldsymbol{\mathit{c}}-A^{T}\accentset{*}{\boldsymbol{\mathit{y}}})-(\boldsymbol{\mathit{c}}-A^{T}\accentset{*}{\boldsymbol{\mathit{y}}}^{+})\right)^{T}\boldsymbol{\mathit{x}}^{+}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}(\accentset{*}{\boldsymbol{\mathit{y}}}-\accentset{*}{\boldsymbol{\mathit{y}}}^{+})^{T}A\boldsymbol{\mathit{x}}^{+}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}(\accentset{*}{\boldsymbol{\mathit{y}}}-\accentset{*}{\boldsymbol{\mathit{y}}}^{+})^{T}\boldsymbol{\mathit{b}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}(\accentset{*}{\boldsymbol{\mathit{y}}}-\accentset{*}{\boldsymbol{\mathit{y}}}^{+})^{T}A\boldsymbol{\mathit{x}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}\left((\boldsymbol{\mathit{c}}-A^{T}\accentset{*}{\boldsymbol{\mathit{y}}})-(\boldsymbol{\mathit{c}}-A^{T}\accentset{*}{\boldsymbol{\mathit{y}}}^{+})\right)^{T}\boldsymbol{\mathit{x}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}(\accentset{*}{\boldsymbol{\mathit{s}}}-\accentset{*}{\boldsymbol{\mathit{s}}}^{+})^{T}\boldsymbol{\mathit{x}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle\leq\left(1+\frac{1}{\accentset{*}{s}_{gap}^{+}}\accentset{*}{\boldsymbol{\mathit{s}}}^{T}\boldsymbol{\mathit{x}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle=\left(1+\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}\right)\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle\leq\left(1+\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}\right)^{2}$

So

$\displaystyle\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}$
$\displaystyle\geq e^{\frac{B(z)-B(z^{+})}{2m}}-1$

Using Lemma C.4, we may conclude

$\frac{s_{gap}}{s_{gap}^{+}}=\frac{s_{gap}}{\accentset{*}{s}_{gap}}\cdot\frac{\accentset{*}{s}_{gap}}{\accentset{*}{s}_{gap}^{+}}\cdot\frac{\accentset{*}{s}_{gap}^{+}}{s_{gap}^{+}}\geq\frac{1-\frac{\eta}{1-\eta}}{1+\frac{\eta}{1-\eta}}\cdot\left(e^{\frac{B(z)-B(z^{+})}{2m}}-1\right)$

∎

######  Corollary C.12.

The $\mathtt{InteriorPoint}$ algorithm makes
$\mathcal{O}\left(\sqrt{m}\log\frac{s_{gap}^{C}}{\epsilon}\right)$ calls to
$\mathtt{Shift}$.

###### Proof.

Recall that the algorithm will terminate only when the value of $s_{gap}$ has
decreased from its initial value of $s_{gap}^{C}$ to below
$\frac{\epsilon}{3}$. Thus, Lemma C.11 ensures us that $s_{gap}$ will be
smaller than $\frac{\epsilon}{3}$ once $B(z)$ has decreased by
$\Omega\left(m\log\frac{s_{gap}^{C}}{\epsilon}\right)$. According to Lemma
C.10, this occurs after
$\mathcal{O}\left(\sqrt{m}\log\frac{s_{gap}^{C}}{\epsilon}\right)$
$\mathtt{Shift}$ iterations. ∎

###  C.3 Finding the Central Path

It remains for us to describe how to initialize the path-following algorithm
by finding a point near the central path. Essentially, this is accomplished by
running the path-following algorithm in reverse. Instead of stepping towards
the optimum given by $\boldsymbol{\mathit{b}}$, we step away from the optimum
given by the vector
$\underaccent{\bar}{\bb}=A(S^{0})^{-1}{\mbox{\boldmath$1$}}_{m}$ that depends
on our initial feasible point
$(\boldsymbol{\mathit{y}}^{0},\boldsymbol{\mathit{s}}^{0})\in\mathring{\Omega}^{D}$.

Our analysis parallels that in the previous section. The following function
$\underaccent{\bar}{\eta}$ measures the proximity of a point
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})\in\mathring{\Omega}^{D}_{\underaccent{\bar}{\bb},\underaccent{\bar}{z}}$
to the central path given by $\underaccent{\bar}{\bb}$:

$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})=\eta_{\tilde{\underaccent{\bar}{A}}}(\tilde{\underaccent{\bar}{\ss}})\qquad\text{where
$\tilde{\underaccent{\bar}{A}}=\left[\begin{matrix}A&-\underaccent{\bar}{\bb}{\mbox{\boldmath$1$}}_{m}^{T}\end{matrix}\right]$
and
$\tilde{\underaccent{\bar}{\ss}}=\left[\begin{matrix}\boldsymbol{\mathit{s}}\\\
\underaccent{\bar}{s}_{gap}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$}$

To initialize the algorithm, we observe that
$(\boldsymbol{\mathit{y}}^{0},\boldsymbol{\mathit{s}}^{0},m)\in\Omega^{D}_{\underaccent{\bar}{\bb},\underaccent{\bar}{z}^{0}}$
is on the $\underaccent{\bar}{\bb}$ central path, where we define
$\underaccent{\bar}{z}^{0}=\underaccent{\bar}{\bb}^{T}\boldsymbol{\mathit{y}}^{0}-m$:

######  Lemma C.13.

$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}}^{0},m)=0$

###### Proof.

Defining
$\tilde{\underaccent{\bar}{\ss}}^{0}=\left[\begin{matrix}\boldsymbol{\mathit{s}}^{0}\\\
m{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$, we have

$\tilde{\underaccent{\bar}{A}}(\tilde{\underaccent{\bar}{S}}^{0})^{-1}{\mbox{\boldmath$1$}}_{2m}=\left[\begin{matrix}A&-\underaccent{\bar}{\bb}{\mbox{\boldmath$1$}}_{m}^{T}\end{matrix}\right]\left[\begin{matrix}(S^{0})^{-1}{\mbox{\boldmath$1$}}_{m}\\\
m^{-1}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]=A(S^{0})^{-1}{\mbox{\boldmath$1$}}_{m}-\underaccent{\bar}{\bb}\cdot\frac{{\mbox{\boldmath$1$}}_{m}^{T}{\mbox{\boldmath$1$}}_{m}}{m}=\underaccent{\bar}{\bb}-\underaccent{\bar}{\bb}=0$

Thus,
$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}}^{0},m)=\left\|\tilde{\underaccent{\bar}{A}}(\tilde{\underaccent{\bar}{S}}^{0})^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(\tilde{\underaccent{\bar}{A}}(\tilde{\underaccent{\bar}{S}}^{0})^{-2}\tilde{\underaccent{\bar}{A}}^{T})^{-1}}=0$

∎

$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap},z)=\mathtt{FindCentralPath}(A,\boldsymbol{\mathit{b}},\boldsymbol{\mathit{c}},\boldsymbol{\mathit{y}}^{0})$
• Define $\displaystyle\left[\begin{matrix}\boldsymbol{\mathit{s}}^{0}\\\
\underaccent{\bar}{s}_{gap}^{0}\end{matrix}\right]$
$\displaystyle=\left[\begin{matrix}\boldsymbol{\mathit{c}}-A^{T}\boldsymbol{\mathit{y}}^{0}\\\
m\end{matrix}\right]$ $\displaystyle\underaccent{\bar}{\bb}$
$\displaystyle=A(S^{0})^{-1}{\mbox{\boldmath$1$}}_{m}$
$\displaystyle\underaccent{\bar}{z}^{0}$
$\displaystyle=\underaccent{\bar}{\bb}^{T}\boldsymbol{\mathit{y}}^{0}-\underaccent{\bar}{s}_{gap}^{0}$
• Set
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap},\underaccent{\bar}{z}):=(\boldsymbol{\mathit{y}}^{0},\boldsymbol{\mathit{s}}^{0},\underaccent{\bar}{s}_{gap}^{0},\underaccent{\bar}{z}^{0})$
• While
$\underaccent{\bar}{s}_{gap}<40\lambda_{min}^{-1/2}Tm\left\|\underaccent{\bar}{\bb}\right\|$:
– Set
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap},\underaccent{\bar}{z}):=\mathtt{Unshift}\left(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap},\underaccent{\bar}{z}\right)$
• Return $(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}})$ and
$z=\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-40\lambda_{min}^{-1/2}Tm\left\|\boldsymbol{\mathit{b}}\right\|$
and $s_{gap}=\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-z$
$(\boldsymbol{\mathit{y}}^{+},\boldsymbol{\mathit{s}}^{+},\underaccent{\bar}{s}_{gap}^{+},\underaccent{\bar}{z}^{+})=\mathtt{Unshift}(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap},\underaccent{\bar}{z})$
• Let
$\underaccent{\bar}{z}^{+}=\underaccent{\bar}{z}-\frac{\underaccent{\bar}{s}_{gap}}{10\sqrt{m}}$
• Let
$\boldsymbol{\mathit{y}}^{+}=\mathtt{NewtonStep}(\tilde{\underaccent{\bar}{A}},\tilde{\underaccent{\bar}{\cc}},\boldsymbol{\mathit{y}})$
where
$\tilde{\underaccent{\bar}{A}}=\left[\begin{matrix}A&-\underaccent{\bar}{\bb}{\mbox{\boldmath$1$}}_{m}^{T}\end{matrix}\right]$
and
$\tilde{\underaccent{\bar}{\cc}}=\left[\begin{matrix}\boldsymbol{\mathit{c}}\\\
-\underaccent{\bar}{z}^{+}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$ • Let
$\left[\begin{matrix}\boldsymbol{\mathit{s}}^{+}\\\
\underaccent{\bar}{s}_{gap}^{+}\end{matrix}\right]=\left[\begin{matrix}\boldsymbol{\mathit{c}}-A^{T}\boldsymbol{\mathit{y}}^{+}\\\
\underaccent{\bar}{\bb}^{T}\boldsymbol{\mathit{y}}^{+}-\underaccent{\bar}{z}^{+}\end{matrix}\right]$

Figure 5: Algorithm for finding point near central path given feasible
interior point

We present the $\mathtt{FindCentralPath}$ algorithm in Figure 5. Starting with
$\underaccent{\bar}{z}=\underaccent{\bar}{z}^{0}$, we take steps along the
$\underaccent{\bar}{\bb}$ central path, decreasing $\underaccent{\bar}{z}$
until it is sufficiently small that the analytic center of
$\Omega^{D}_{\underaccent{\bar}{\bb},\underaccent{\bar}{z}}$ is close to the
analytic center of $\Omega^{D}$, and therfore also close to the analytic
center of $\Omega^{D}_{\boldsymbol{\mathit{b}},z}$ for some sufficiently small
$z$.

Let us show that the $\mathtt{Unshift}$ procedure indeed takes steps near the
$\underaccent{\bar}{\bb}$ central path:

######  Lemma C.14 (compare Lemmas C.8).

Given
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})\in\mathring{\Omega}^{D}_{\underaccent{\bar}{\bb},z}$
satisfying
$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})\leq\frac{1}{40}$.
Let
$\underaccent{\bar}{s}_{gap}^{\prime}=\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-\underaccent{\bar}{z}^{+}$
and
$(\boldsymbol{\mathit{y}}^{+},\boldsymbol{\mathit{s}}^{+},\underaccent{\bar}{s}_{gap}^{+},\underaccent{\bar}{z}^{+})=\mathtt{Unshift}(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap},\underaccent{\bar}{z})$
. Then

  1. (i)

$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}^{\prime}_{gap})\leq\frac{51}{400}$

  2. (ii)

$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}}^{+},\underaccent{\bar}{s}^{+}_{gap})\leq\frac{1}{40}$.

###### Proof of C.14(i).

Following the proof of Lemma C.8(i) through equation 9, we have

$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap}^{\prime})\leq\frac{11}{10}\cdot\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})+\frac{1}{10}\leq\frac{11}{10}\cdot\frac{1}{40}+\frac{1}{10}=\frac{51}{400}$

∎

###### Proof of C.14(ii).

By Lemma C.5, we have

$\tilde{\eta}(\boldsymbol{\mathit{s}}^{+},\underaccent{\bar}{s}_{gap}^{+})\leq\tilde{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap}^{\prime})^{2}+\frac{1}{20}\tilde{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap}^{\prime})\leq\left(\frac{51}{400}\right)^{2}+\frac{1}{20}\cdot\frac{51}{400}<\frac{1}{40}$

∎

Next, let us prove that the point returned by $\mathtt{FindCentralPath}$ is
indeed near the original central path (i.e. the path given by
$\boldsymbol{\mathit{b}}$):

######  Lemma C.15.

For $\boldsymbol{\mathit{y}}^{0}$ satisfying
$A^{T}\boldsymbol{\mathit{y}}^{0}<\boldsymbol{\mathit{c}}$, let
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap},z)=\mathtt{FindCentralPath}(A,\boldsymbol{\mathit{b}},\boldsymbol{\mathit{c}},\boldsymbol{\mathit{y}}^{0})$.

Then
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap})\in\mathring{\Omega}^{D}_{\boldsymbol{\mathit{b}},z}$
and $\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})\leq\frac{1}{10}$.

###### Proof.

Using the values at the end of the algorithm, we write
$\tilde{\boldsymbol{\mathit{s}}}=\left[\begin{matrix}\boldsymbol{\mathit{s}}\\\
s_{gap}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$ and
$\tilde{\underaccent{\bar}{\ss}}=\left[\begin{matrix}\boldsymbol{\mathit{s}}\\\
\underaccent{\bar}{s}_{gap}{\mbox{\boldmath$1$}}_{m}\end{matrix}\right]$.

To begin, we note

$\displaystyle\underaccent{\bar}{s}_{gap}$ $\displaystyle\geq
40\lambda_{min}^{-1/2}Tm\left\|\underaccent{\bar}{\bb}\right\|$
$\displaystyle=40m(T^{-2}\lambda_{min})^{-1/2}\left\|\underaccent{\bar}{\bb}\right\|$
$\displaystyle\geq
40m\left\|\underaccent{\bar}{\bb}\right\|_{(AS^{-2}A^{T})^{-1}}$ (15)

where the last inequality follows because the smallest eigenvalue of
$AS^{-2}A^{T}$ is at least $T^{-2}\lambda_{min}$.

Similarly,

$s_{gap}=40\lambda_{min}^{-1/2}Tm\left\|\boldsymbol{\mathit{b}}\right\|\geq
40m\left\|\boldsymbol{\mathit{b}}\right\|_{(AS^{-2}A^{T})^{-1}}$ (16)

We have

$\displaystyle\tilde{\eta}(\boldsymbol{\mathit{s}},s_{gap})$
$\displaystyle=\left\|\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(\tilde{A}\tilde{S}^{-2}\tilde{A}^{T})^{-1}}$
$\displaystyle\leq\left\|\tilde{A}\tilde{S}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(AS^{-2}A^{T})^{-1}}$
(because
$\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}-AS^{-2}A^{T}=ms_{gap}^{-2}\boldsymbol{\mathit{b}}\boldsymbol{\mathit{b}}^{T}$
is positive semidefinite)
$\displaystyle=\left\|\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-1}{\mbox{\boldmath$1$}}_{2m}-ms_{gap}^{-1}\boldsymbol{\mathit{b}}+m\underaccent{\bar}{s}_{gap}^{-1}\underaccent{\bar}{\bb}\right\|_{(AS^{-2}A^{T})^{-1}}$
$\displaystyle\leq\left\|\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(AS^{-2}A^{T})^{-1}}+m\underaccent{\bar}{s}_{gap}^{-1}\left\|\underaccent{\bar}{\bb}\right\|_{(AS^{-2}A^{T})^{-1}}+ms_{gap}^{-1}\left\|\boldsymbol{\mathit{b}}\right\|_{(AS^{-2}A^{T})^{-1}}$
$\displaystyle\leq\left\|\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(AS^{-2}A^{T})^{-1}}+\frac{1}{40}+\frac{1}{40}\qquad\text{(by
equations \ref{eq:unshift-hatsgap} and \ref{eq:unshift-sgap})}$
$\displaystyle\leq\left(1+m\underaccent{\bar}{s}_{gap}^{-1}\left\|\underaccent{\bar}{\bb}\right\|_{(AS^{-2}A^{T})^{-1}}\right)^{1/2}\left\|\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-2}\tilde{\underaccent{\bar}{A}}^{T})^{-1}}+\frac{1}{40}+\frac{1}{40}$
(by Lemma B.3, using the fact that
$\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-2}\tilde{\underaccent{\bar}{A}}^{T}-AS^{-2}A^{T}=m\underaccent{\bar}{s}_{gap}^{-2}\underaccent{\bar}{\bb}\underaccent{\bar}{\bb}^{T}$)
$\displaystyle\leq\left(1+\frac{1}{40}\right)^{1/2}\left\|\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-1}{\mbox{\boldmath$1$}}_{2m}\right\|_{(\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-2}\tilde{\underaccent{\bar}{A}}^{T})^{-1}}+\frac{1}{40}+\frac{1}{40}\qquad\text{(by
equation \ref{eq:unshift-hatsgap})}$
$\displaystyle=\left(1+\frac{1}{40}\right)^{1/2}\cdot\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})+\frac{1}{40}+\frac{1}{40}$
$\displaystyle\leq
2\cdot\underaccent{\bar}{\eta}(\tilde{\underaccent{\bar}{\ss}},\underaccent{\bar}{s}_{gap})+\frac{1}{20}$
$\displaystyle\leq 2\cdot\frac{1}{40}+\frac{1}{20}\qquad\text{(by Lemma
\ref{lem:unshift}(ii))}$ $\displaystyle=\frac{1}{10}$

∎

To measure the progress of the $\mathtt{FindCentralPath}$ algorithm, we define
$\underaccent{\bar}{B}(\underaccent{\bar}{z})$:

$\underaccent{\bar}{B}(\underaccent{\bar}{z})=\sum_{j=1}^{m}\log\accentset{*}{\boldsymbol{\mathit{s}}}_{j}+m\log\accentset{*}{\underaccent{\bar}{s}}_{gap}\qquad\text{where
$(\accentset{*}{\boldsymbol{\mathit{y}}},\accentset{*}{\boldsymbol{\mathit{s}}},\accentset{*}{\underaccent{\bar}{s}}_{gap})$
is the analytic center of
$\Omega^{D}_{\underaccent{\bar}{\bb},\underaccent{\bar}{z}}$}$

######  Lemma C.16 (compare Lemma C.10).

Given
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})\in\mathring{\Omega}^{D}_{\underaccent{\bar}{\bb},z}$
satisfying
$\underaccent{\bar}{\eta}(\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap})\leq\frac{1}{40}$.
Let
$(\boldsymbol{\mathit{y}}^{+},\boldsymbol{\mathit{s}}^{+},\underaccent{\bar}{s}_{gap}^{+},\underaccent{\bar}{z}^{+})=\mathtt{Unshift}(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},\underaccent{\bar}{s}_{gap},\underaccent{\bar}{z})$.
Then
$\underaccent{\bar}{B}(\underaccent{\bar}{z}^{+})\geq\underaccent{\bar}{B}(\underaccent{\bar}{z})+\Theta(\sqrt{m})$.

###### Proof.

We will follow the proof of Lemma C.10, with some minor changes.

Before we proceed, let us recall the definition
$\underaccent{\bar}{z}^{+}=\underaccent{\bar}{z}-\frac{\underaccent{\bar}{s}_{gap}}{10\sqrt{m}}$
to note that

$\underaccent{\bar}{s}_{gap}^{\prime}=\underaccent{\bar}{s}_{gap}+\underaccent{\bar}{z}-\underaccent{\bar}{z}^{+}=10\sqrt{m}(\underaccent{\bar}{z}-\underaccent{\bar}{z}^{+})+\underaccent{\bar}{z}-\underaccent{\bar}{z}^{+}\leq
11\sqrt{m}(\underaccent{\bar}{z}-\underaccent{\bar}{z}^{+})$ (17)

Now, we switch the places of $z$ and $z^{+}$, and follow the proof of Lemma
C.10 up to Equation 14:

$\displaystyle
e^{\frac{\underaccent{\bar}{B}(\underaccent{\bar}{z})-\underaccent{\bar}{B}(\underaccent{\bar}{z}^{+})}{2m}}$
$\displaystyle\leq
1+\frac{1}{2\accentset{*}{\underaccent{\bar}{s}}_{gap}^{+}}\cdot(\underaccent{\bar}{z}^{+}-\underaccent{\bar}{z})$
We continue: $\displaystyle\leq
1-\frac{\underaccent{\bar}{s}_{gap}^{\prime}}{22\sqrt{m}\cdot\accentset{*}{\underaccent{\bar}{s}}_{gap}^{+}}\qquad\text{(by
Equation \ref{eq:stosprime})}$ $\displaystyle\leq
1-\frac{1}{22\sqrt{m}}\cdot\frac{349}{400}\qquad\text{(by Lemmas
\ref{lem:unshift} and \ref{lem:nearcenter})}$ $\displaystyle\leq
e^{-\frac{169}{4400\sqrt{m}}}$

∎

######  Corollary C.17.

The $\mathtt{FindCentralPath}$ algorithm makes
$O\left(\sqrt{m}\log\frac{TUm}{\lambda_{min}s^{0}_{min}}\right)$ calls to
$\mathtt{Unshift}$, where $s^{0}_{min}$ is the smallest entry of
$\boldsymbol{\mathit{s}}^{0}=\boldsymbol{\mathit{c}}-A^{T}\boldsymbol{\mathit{y}}^{0}$.

###### Proof.

Recall that the algorithm will terminate only when the value of
$\underaccent{\bar}{s}_{gap}$ has increased from its initial value of $m$ to
at least $40\lambda_{min}^{-1/2}Tm\left\|\underaccent{\bar}{\bb}\right\|$. So,
by Lemma C.11, this will have happened once $\underaccent{\bar}{B}(z)$ has
increased by
$\Omega\left(m\log\left(\lambda_{min}^{-1/2}T\left\|\underaccent{\bar}{\bb}\right\|\right)\right)$
.

According to Lemma C.16, this occurs after
$O\left(\sqrt{m}\log\left(\lambda_{min}^{-1/2}T\left\|\underaccent{\bar}{\bb}\right\|\right)\right)$
iterations.

To complete the proof, we note that

$\left\|\underaccent{\bar}{\bb}\right\|=\left\|A(S^{0})^{-1}{\mbox{\boldmath$1$}}_{m}\right\|\leq\frac{n^{1/2}mU}{s^{0}_{min}}$

∎

###  C.4 Calls to the Solver

In each call to $\mathtt{Unshift}$, we solve one system in a matrix of the
form

$\tilde{\underaccent{\bar}{A}}\tilde{\underaccent{\bar}{S}}^{-2}\tilde{\underaccent{\bar}{A}}^{T}=AS^{-2}A^{T}+m\underaccent{\bar}{s}_{gap}^{-2}\underaccent{\bar}{\bb}\underaccent{\bar}{\bb}^{T}$

and in each call to $\mathtt{Shift}$, we solve one system in a matrix of the
form

$\tilde{A}\tilde{S}^{-2}\tilde{A}^{T}=AS^{-2}A^{T}+ms_{gap}^{-2}\boldsymbol{\mathit{b}}\boldsymbol{\mathit{b}}^{T}$

At the end of the interior-point algorithm we have one final call of the
latter form.

In order to say something about the condition number of the above matrices, we
must bound the slack vector $\boldsymbol{\mathit{s}}$. We are given an upper
bound of $T$ on the elements of $\boldsymbol{\mathit{s}}$, so it remains to
prove a lower bound:

######  Lemma C.18.

Throughout the $\mathtt{InteriorPoint}$ algorithm,
$\boldsymbol{\mathit{s}}\geq\frac{\epsilon}{48nmTU}\boldsymbol{\mathit{s}}^{0}$

###### Proof.

At all times duing the algorithm, we know from Lemma C.4 that the elements of
$\boldsymbol{\mathit{s}}$ are bounded by a constant factor from the slacks at
the current central path point $\accentset{*}{\boldsymbol{\mathit{s}}}$. In
particular, taking into account Lemmas C.8 and C.14, we surely have
$\boldsymbol{\mathit{s}}\geq\frac{1}{2}\accentset{*}{\boldsymbol{\mathit{s}}}$.
So let us bound from below the elements of
$\accentset{*}{\boldsymbol{\mathit{s}}}$.

During the $\mathtt{FindCentralPath}$ subroutine, as we decrease
$\underaccent{\bar}{z}$ and expand the polytope
$\Omega^{D}_{\underaccent{\bar}{\bb},\underaccent{\bar}{z}}$, clearly the
initial point $\boldsymbol{\mathit{s}}^{0}$ remains in the interior of
$\Omega^{D}_{\underaccent{\bar}{\bb},\underaccent{\bar}{z}}$ throughout. Thus,
by Lemma C.3, we have
$\accentset{*}{\boldsymbol{\mathit{s}}}\geq\frac{1}{2m}\boldsymbol{\mathit{s}}^{0}$,
and so
$\boldsymbol{\mathit{s}}\geq\frac{1}{2}\accentset{*}{\boldsymbol{\mathit{s}}}\geq\frac{1}{4m}\boldsymbol{\mathit{s}}^{0}$.

Unfortunately, during the main part of the algorithm, as we increase $z$ and
shrink the polytope $\Omega^{D}_{\boldsymbol{\mathit{b}},z}$, the initial
point may not remain inside the polytope. In particular, once we have
$z\geq\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}^{0}$, the initial
point is no longer in $\Omega^{D}_{\boldsymbol{\mathit{b}},z}$, but we may
define a related point
$(\boldsymbol{\mathit{y}}^{z},\boldsymbol{\mathit{s}}^{z},s_{gap}^{z})$ that
is in $\Omega^{D}_{\boldsymbol{\mathit{b}},z}$.

Given our current point
$(\boldsymbol{\mathit{y}},\boldsymbol{\mathit{s}},s_{gap})\in\Omega^{D}_{\boldsymbol{\mathit{b}},z}$
for $z\geq\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}^{0}$, let us
define
$r=\frac{\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-z}{2(\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}^{0})}$
and note that $0<r<\frac{1}{2}$. We then define

$\displaystyle\boldsymbol{\mathit{y}}^{z}$
$\displaystyle=r\boldsymbol{\mathit{y}}^{0}+(1-r)\boldsymbol{\mathit{y}}$
$\displaystyle\left[\begin{matrix}\boldsymbol{\mathit{s}}^{z}\\\
s_{gap}^{z}\end{matrix}\right]$
$\displaystyle=\left[\begin{matrix}\boldsymbol{\mathit{c}}-A^{T}\boldsymbol{\mathit{y}}^{z}\\\
\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}^{z}-z\end{matrix}\right]=\left[\begin{matrix}r\boldsymbol{\mathit{s}}^{0}+(1-r)\boldsymbol{\mathit{s}}\\\
\frac{1}{2}(\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-z)\end{matrix}\right]>0$

Therefore Lemma C.3 gives

$\accentset{*}{\boldsymbol{\mathit{s}}}\geq\frac{1}{2m}\boldsymbol{\mathit{s}}^{z}=\frac{r\boldsymbol{\mathit{s}}^{0}+(1-r)\boldsymbol{\mathit{s}}}{2m}\geq\frac{r}{2m}\boldsymbol{\mathit{s}}^{0}$

We then find

$r=\frac{s_{gap}}{2(\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}-\boldsymbol{\mathit{b}}^{T}\boldsymbol{\mathit{y}}^{0})}\geq\frac{s_{gap}}{4nTU}\geq\frac{\epsilon}{24nTU}$

The last inequality follows because, when $s_{gap}$ decreased below
$\frac{\epsilon}{3}$ on the final step, using Lemma C.4 we find that it
certainly could not have decreased by more than a factor of $\frac{1}{2}$.

We conclude
$\boldsymbol{\mathit{s}}\geq\frac{1}{2}\accentset{*}{\boldsymbol{\mathit{s}}}\geq\frac{r}{4m}\boldsymbol{\mathit{s}}^{0}\geq\frac{\epsilon}{48nmTU}\boldsymbol{\mathit{s}}^{0}$
∎

We may now summarize the calls to the solver as follows:

######  Theorem C.19.

The
$\mathtt{InteriorPoint}(A,\boldsymbol{\mathit{b}},\boldsymbol{\mathit{c}},\boldsymbol{\mathit{y}}^{0},\epsilon)$
algorithm makes
$\mathcal{O}\left(\sqrt{m}\log\frac{TUm}{\lambda_{min}s^{0}_{min}\epsilon}\right)$
calls to the approximate solver, of the form

$\mathtt{Solve}\left(AS^{-2}A^{T}+\boldsymbol{\mathit{v}}\boldsymbol{\mathit{v}}^{T},\cdot,\Theta(m^{-1/2})\right)$

and one call of the form

$\mathtt{Solve}\left(AS^{-2}A^{T}+\boldsymbol{\mathit{v}}\boldsymbol{\mathit{v}}^{T},\cdot,\Omega\left(\frac{s^{0}_{min}\epsilon}{m^{1/2}n^{2}T^{2}U^{2}}\right)\right)$

where $S$ is a positive diagonal matrix with condition number
$\mathcal{O}\left(\frac{T^{2}Um^{2}}{\epsilon}\right)$, and
$\boldsymbol{\mathit{v}}$ satisfies

$\left\|\boldsymbol{\mathit{v}}\right\|=\mathcal{O}\left(\frac{U(mn)^{1/2}}{s^{0}_{min}\epsilon}\right)$

###### Proof.

From Lemmas C.17 and C.12, the total number of solves is
$\mathcal{O}\left(\sqrt{m}\left(\log\frac{TUm}{\lambda_{min}s^{0}_{min}}+\log\frac{s_{gap}^{C}}{\epsilon}\right)\right)$,
where we know from the $\mathtt{FindCentralPath}$ algorithm that
$s_{gap}^{C}=40\frac{Tm\left\|\boldsymbol{\mathit{b}}\right\|}{\lambda_{min}^{1/2}}=\mathcal{O}\left(\frac{TUmn^{1/2}}{\lambda_{min}^{1/2}}\right)$

As we noted above, all solves are in matrices that take the form
$AS^{-2}A^{T}+\boldsymbol{\mathit{v}}\boldsymbol{\mathit{v}}^{T}$, where

$\displaystyle\boldsymbol{\mathit{v}}$
$\displaystyle=m^{1/2}s_{gap}^{-1}\boldsymbol{\mathit{b}}$ or
$\displaystyle\boldsymbol{\mathit{v}}$
$\displaystyle=m^{1/2}\underaccent{\bar}{s}_{gap}^{-1}A(S^{0})^{-1}{\mbox{\boldmath$1$}}_{m}$
We know that $s_{gap}=\Omega(\epsilon)$ and
$\underaccent{\bar}{s}_{gap}=\Omega(m)$, so we obtain the respective bounds
$\displaystyle\left\|\boldsymbol{\mathit{v}}\right\|$
$\displaystyle=\mathcal{O}\left(\frac{U(mn)^{1/2}}{\epsilon}\right)$
$\displaystyle\left\|\boldsymbol{\mathit{v}}\right\|$
$\displaystyle=\mathcal{O}\left(\frac{U(mn)^{1/2}}{s^{0}_{min}}\right)$

The condition number of $S$ comes from Lemma C.18 and the upper bound of $T$
on the slacks.

The error parameter for the solver is $\Theta\left(m^{-1/2}\right)$ from the
all $\mathtt{NewtonStep}$ calls. In the final solve, the error parameter is
$\frac{s_{min}m^{1/2}}{TUn}\geq\frac{m^{1/2}}{TUn}\cdot\frac{s^{0}_{min}\epsilon}{48nmTU}$,
again invoking Lemma C.18. ∎