New challenges in fleet deployment considering EU oil sanctions

Yiwei Wu; Yao Lu; Shuaian Wang; Lu Zhen; Yiwei Wu; Yao Lu; Shuaian Wang; Lu Zhen

doi:10.3934/era.2023230

Electronic Research Archive

2023, Volume 31, Issue 8: 4507-4529. doi: 10.3934/era.2023230

Previous Article Next Article

Research article Special Issues

New challenges in fleet deployment considering EU oil sanctions

1.
Department of Logistics and Maritime Studies, Faculty of Business, The Hong Kong Polytechnic University, Hung Hom, Kowloon, Hong Kong 999077, China
2.
School of Maritime Economics and Management, Dalian Maritime University, Dalian 116026, China
3.
School of Management, Shanghai University, Shanghai 200436, China

Received: 29 April 2023 Revised: 28 May 2023 Accepted: 30 May 2023 Published: 16 June 2023

Due to European Union (EU) oil sanctions, tanker shipping companies need to redeploy their tankers by moving tankers between ship routes with the consideration of flag states of tankers, but the literature lacks quantitative methods for this problem. To fill this research gap, this paper studies an integrated problem of fleet deployment, fleet repositioning, round trip completion, and speed optimization with the consideration of flag states of tankers. The problem is formulated as a nonlinear integer programming model to minimize the total cost, including the fleet repositioning cost, the mismatch cost, and the fuel cost, during the planning period while satisfying the total crude oil transportation demand of each voyage and the minimum shipping frequency. Some linearization methods are used to transform the nonlinear model to a linear one which can be directly solved by Gurobi. The average solving time required for 17 computational instances is 4.5 minutes, which validates the effectiveness of the proposed model. Sensitivity analyses, including the influences of the unit fuel price, the total crude oil transportation demand, the mismatch cost of completing a round trip by a deployed tanker, and the repositioning cost for each deployed tanker, on operations decisions, are conducted to obtain managerial insights.

Keywords:

Citation: Yiwei Wu, Yao Lu, Shuaian Wang, Lu Zhen. New challenges in fleet deployment considering EU oil sanctions[J]. Electronic Research Archive, 2023, 31(8): 4507-4529. doi: 10.3934/era.2023230

Related Papers:

[1]	Liping Yang, Hu Li . A hybrid algorithm based on parareal and Schwarz waveform relaxation. Electronic Research Archive, 2022, 30(11): 4086-4107. doi: 10.3934/era.2022207
[2]	Chunhua Chu, Yijun Chen, Qun Zhang, Ying Luo . Joint linear array structure and waveform design for MIMO radar under practical constraints. Electronic Research Archive, 2022, 30(9): 3249-3265. doi: 10.3934/era.2022165
[3]	Xuerong Hu, Yuxiang Han, Junyan Lu, Linxiang Wang . Modeling and experimental investigation of quasi-zero stiffness vibration isolator using shape memory alloy springs. Electronic Research Archive, 2025, 33(2): 768-790. doi: 10.3934/era.2025035
[4]	Qian He, Wenxin Du, Feng Shi, Jiaping Yu . A fast method for solving time-dependent nonlinear convection diffusion problems. Electronic Research Archive, 2022, 30(6): 2165-2182. doi: 10.3934/era.2022109
[5]	Li Tian, Ziqiang Wang, Junying Cao . A high-order numerical scheme for right Caputo fractional differential equations with uniform accuracy. Electronic Research Archive, 2022, 30(10): 3825-3854. doi: 10.3934/era.2022195
[6]	Yi Wei . The Riccati-Bernoulli subsidiary ordinary differential equation method to the coupled Higgs field equation. Electronic Research Archive, 2023, 31(11): 6790-6802. doi: 10.3934/era.2023342
[7]	Adrian Nicolae Branga . Some conditions for the existence and uniqueness of monotonic and positive solutions for nonlinear systems of ordinary differential equations. Electronic Research Archive, 2022, 30(6): 1999-2017. doi: 10.3934/era.2022101
[8]	Junseok Kim . Maximum principle preserving the unconditionally stable method for the Allen–Cahn equation with a high-order potential. Electronic Research Archive, 2025, 33(1): 433-446. doi: 10.3934/era.2025021
[9]	Martin Bohner, Alexander Domoshnitsky, Oleg Kupervasser, Alex Sitkin . Floquet theory for first-order delay equations and an application to height stabilization of a drone's flight. Electronic Research Archive, 2025, 33(5): 2840-2861. doi: 10.3934/era.2025125
[10]	David Cheban, Zhenxin Liu . Averaging principle on infinite intervals for stochastic ordinary differential equations. Electronic Research Archive, 2021, 29(4): 2791-2817. doi: 10.3934/era.2021014

Abstract

1. Introduction

Waveform relaxation (WR) is an iterative method for systems of ordinary differential equations (ODEs). It is introduced for the first time by Lelarasmee et al. ^[1] for the time domain analysis of large-scale nonlinear dynamical systems. A lot of studies have been done on WR methods and most of them focus on convergence (see ^{[2,3,4,5,6,7,8,9]}). Zero-stability is the basic property of numerical methods of ODEs^[10]. However, to our best knowledge, so far there exists little work on zero-stability of WR methods.

We first of all propose two examples to show the fact that tiny perturbation may lead to huge change of WR approximation solutions, at the same time, to confirm that it is necessary to study zero-stability of WR methods. For this, it is enough to consider a simple method, and suppose only the initial value has tiny perturbations.

Example 1. Consider the WR method

$\begin{equation} x^{(k+1)}_{n+2}-(1+\alpha)x^{(k+1)}_{n+1}+\alpha x^{(k+1)}_{n} = -(1-\alpha)h(x^{(k+1)}_{n}+x^{(k)}_{n}) \end{equation}$

(1.1)

of solving numerically the initial value problem

$\dot{x}(t) = -2x(t), t\in[0, T]; x(0) = x_0,$

which is gotten by applying the two step method

$x_{n+2}-(1+\alpha)x_{n+1}+\alpha x_n = -2(1-\alpha)hx_n$

that is consistent for any $\alpha$ and zero-stable for $-1 < \alpha < 1$ to the iterative scheme

$\dot{x}^{(k+1)}(t) = -x^{(k+1)}(t)-x^{(k)}(t).$

Suppose that

$\begin{equation} N \;{\rm{ is\;a\;positive\;integer\;and }}\; h = \dfrac{T}{N}, \end{equation}$

(1.2)

$\begin{equation} x^{(0)}_n = x_0, n = 0, 1, \cdots, N, x^{(k)}_0 = x_0, k = 0, 1, \cdots \end{equation}$

(1.3)

and

$\begin{equation} x^{(k+1)}_1 = x^{(k+1)}_0-h (x^{(k+1)}_0+ x^{(k)}_0), k = 0, 1, \cdots. \end{equation}$

(1.4)

We can get the approximate solution $\{\{x^{(k)}_n, n = 0, 1, \cdots, N\}, k = 0, 1, \cdots\}$ by using Eqs (1.1)–(1.4). Suppose that tiny perturbations denoted by $\delta$ lead to the initial values become

$\begin{equation} \tilde{x}^{(0)}_n = x_0+\delta, n = 0, 1, \cdots, N, \tilde{x}^{(k)}_0 = x_0+\delta, k = 0, 1, \cdots \end{equation}$

(1.5)

and the resulting perturbed solution of Eq (1.1) is denoted by $\{\tilde{x}^{(k)}_n\}$ . Let $\delta = 1$ for simplicity and denote $\tilde{x}^{(k)}_n-x^{(k)}_n$ with $e^{(k)}_n$ . Then we have

$\begin{equation} \begin{aligned} &e^{(k+1)}_{n+2}-(1+\alpha)e^{(k+1)}_{n+1}+\alpha e^{(k+1)}_{n} = -(1-\alpha)h(e^{(k+1)}_{n}+e^{(k)}_{n})\\ &\qquad \qquad \qquad (n = 0, 1, \cdots, N-2, k = 0, 1, \cdots), \\ &e^{(k+1)}_0 = 1, e^{(k+1)}_1 = (1-2h)e^{(k+1)}_0(k = 0, 1, \cdots), e^{(0)}_n = 1(n = 0, 1, \cdots, N).\\ \end{aligned} \end{equation}$

(1.6)

By virtue of Eq (1.6) the following numerical results are obtained.

The data in and show that errors of the approximate solution of Eq (1.1) brought by perturbations of initial values are controllable for $\alpha = 0.5$ , but they are unbounded for $\alpha = 2$ when $h$ tends to zero.

Table 1. Results obtained by (1.6) with

$\alpha = 0.5$ .

$N$	10	100	250	500	750	1000
$e^{(5)}_N$	0.10621	0.1319	0.13343	0.13388	0.13402	0.13409
$e^{(10)}_N$	0.10621	0.13238	0.13422	0.13478	0.13497	0.13506
$e^{(100)}_N$	0.10621	0.13238	0.13422	0.13478	0.13497	0.13506

| Show Table

DownLoad: CSV

Table 2. Results obtained by (1.6) with

$\alpha = 2$ .

$N$	10	100	250	500	1000
$e^{(5)}_N$	268.74	6.2863e+028	3.7892e+073	3.4922e+148	5.7695e+298
$e^{(10)}_N$	268.74	6.2867e+028	3.7897e+073	3.4927e+148	5.7704e+298
$e^{(100)}_N$	268.74	6.2867e+028	3.7897e+073	3.4927e+148	5.7704e+298

| Show Table

DownLoad: CSV

In this numerical examples, taking $\delta = 1$ instead of tiny perturbations seem to be unreasonable. The main reason for this is that we are usually interested in the ratio $\frac{e^{(k)}_N}{\delta}$ instead of $\delta$ itself. In fact we have also performed some other experiments using smaller $\delta$ , in which the similar ratios are obtained for different $\delta$ .

Example 2. Consider the one-step WR method

$\begin{equation} \begin{aligned} &x^{(k+1)}_{n+1} = x^{(k+1)}_n+h\left(500(x^{(k+1)}_n)^{1.001}+500(x^{(k)}_n)^{1.001}\right), n = 0, 1, \cdots, \dfrac{2}{h}, \\ &x^{(k+1)}_0 = 0;\ x^{(0)}_n = 0, \;{ \rm{ for\;all } }\;n, \end{aligned} \end{equation}$

(1.7)

where $2/h$ is an integer. Suppose that the initial value $0$ becomes $\delta$ by the influence of perturbations and the solution $x_n^{(k)}$ becomes $\tilde{x}_n^{(k)}$ by the influence of $\delta$ , which satisfies

$\begin{equation} \begin{aligned} &\tilde{x}^{(k+1)}_{n+1} = \tilde{x}^{(k+1)}_n+h\left(500(\tilde{x}^{(k+1)}_n)^{1.001}+500(\tilde{x}^{(k)}_n)^{1.001}\right), n = 0, 1, \cdots, \dfrac{2}{h}, \\ &\tilde{x}^{(k+1)}_0 = \delta;\ \tilde{x}^{(0)}_n = \delta, \;{ \rm{ for\;all } }\;n. \end{aligned} \end{equation}$

(1.8)

Clearly the numerical solution generated by Eq (1.7) is the constant zero and Eq (1.8) can be regarded as the WR method of solving the following initial value problem

$\begin{equation} \dot{x}(t) = 1000(x(t))^{1.001}, t\in[0, 2]; x(0) = \delta. \end{equation}$

(1.9)

It is easy to derive that $t_0 = \dfrac{1}{\delta^{0.001}}$ is the blow-up time of the solution of Eq (1.9) and $t_0\in [0, 2]$ whenever $\delta > 10^{-301}$ . Consequently, $\sup\limits_{0 < t < 2}x(t) = +\infty$ for any $\delta > 10^{-301}$ . If the "2" in Eqs (1.8) and (1.9) is replaced by $t_0-\varepsilon$ , where $\varepsilon > 0$ is small enough, then Eq (1.8) converges to Eq (1.9), that is,

$\lim\limits_{h\rightarrow 0, k\rightarrow \infty}\sup\limits_{0\leq n\leq \frac{t_0-\varepsilon}{h}}|\tilde{x}^{(k)}_n-x(nh)| = 0,$

see ^[4]. Thus

$\lim\limits_{h\rightarrow 0, k\rightarrow \infty}\sup\limits_{0\leq n\leq \frac{2}{h}}|\tilde{x}^{(k)}_n| = +\infty, \;{ \rm{ for\;any } }\;\delta > 10^{-301},$

that is, the change of the solution of Eq (1.8) brought by the tiny change of initial values is unbounded as $h$ tends to zero.

The above two examples show that tiny perturbations can lead to huge change of WR approximation solutions, although the methods used in the examples are too simple to be practically useful.

A numerical method is said to be zero-stable if the change of its solution brought by the tiny perturbations is controllable when $h$ tends to zero (see ^[10]). In this paper we will explore what conditions guarantee the zero-stability of WR methods.

2. Preliminary

Consider the initial problem

$\dot{x}(t) = f(t, x(t)), t\in[0, T]; x(0) = x_0 (x(t)\in\mathbb{R}^d).$

Taking the splitting function $F(t, x, x) = f(t, x)$ we can construct the iterative scheme

$\begin{equation} \begin{aligned} &\dot{x}^{(k+1)}(t) = F(t, x^{(k+1)}(t), x^{(k)}(t)), t\in[0, T], k = 0, 1, \cdots;\\ & x^{(0)}(t)\equiv x_0, x^{(k+1)}(0) = x_0. \end{aligned} \end{equation}$

(2.1)

Using one-step methods with variable step-size to discrete Eq (2.1) one arrives at

$\begin{equation} \begin{aligned} &x^{(k+1)}_{n+1} = x^{(k+1)}_n+h_n \Phi(t_{n}, h_n, x^{(k+1)}_{n+1}, x^{(k+1)}_n, x^{(k)}_{n+1}, x^{(k)}_n);x^{(k+1)}_0 = x_0, \end{aligned} \end{equation}$

(2.2)

where $n = 0, 1, \cdots, N-1, k = 0, 1, \cdots$ , $\sum\limits_{n = 0}^{N-1}h_n = T$ and $x^{(0)}_n = x_0$ for all $n$ . Applying a class of multi-step methods with fixed step-size to Eq (2.1) we get

$\begin{equation} x^{(k+1)}_{n+1} = x^{(k+1)}_n+h\Psi(t_n, h, x^{(k+1)}_{n+1}, x^{(k+1)}_{n}, \cdots, x^{(k+1)}_{n-p}, x^{(k)}_{n+1}, x^{(k)}_{n}, \cdots, x^{(k)}_{n-p}), \end{equation}$

(2.3)

where $n = p, p+1, \cdots, N-1, k = 0, 1, \cdots$ , $h = T/N$ , the initial values except $x^{(0)}_n(n = 0, 1, \cdots, N)$ and $x^{(k)}_0(k = 0, 1, \cdots)$ are generated by a suitable one-step method. Here we take the initial values as

$\begin{equation} x^{(k)}_n = c^{(k)}_n, n = 0, 1, \cdots, p, k = 0, 1, \cdots(c^{(k)}_n = x_0, \;{\rm{ if }}\;n\cdot k = 0). \end{equation}$

(2.4)

A classical example of Eq (2.3) is Admas-type linear multi-step methods:

$x^{(k+1)}_{n+1} = x^{(k+1)}_n+h\sum\limits_{j = 0}^{p+1} \beta_j F(t_{n+1-j}, x^{(k+1)}_{n+1-j}, x^{(k)}_{n+1-j}), n = p, p+1, \cdots, N-1,$

where $\beta_j, j = 0, 1, \cdots, p+1$ are constants. The linear multi-step methods use generally the fixed step-size and therefore we only consider the Eq (2.3) with the fixed step-size.

Let $_ix^{(k)}$ denote the component of $x^{(k)}$ satisfying

$x^{(k)} = \left(\left(_1x^{(k)}\right)^{\rm{T}}, \left(_2x^{(k)}\right)^{\rm{T}}, \cdots, \left(_mx^{(k)}\right)^{\rm{T}}\right)^{\rm{T}},$

where $_ix^{(k)}\in \mathbb{R}^{d_i}$ , $x^{(k)}\in \mathbb{R}^{d}$ and $d_1+d_2+\cdots+d_m = d$ . The Eq (2.1) can be rewritten as

$_i\dot{x}^{(k+1)} = F_i(t, \sideset{_1}{}{\mathop{x}^{(k+1)}}, \sideset{_2}{}{\mathop{x}^{(k+1)}}, \cdots, \sideset{_m}{}{\mathop{x}^{(k+1)}}, \sideset{_1}{}{\mathop{x}^{(k)}}, \sideset{_2}{}{\mathop{x}^{(k)}}, \cdots, \sideset{_m}{}{\mathop{x}^{(k)}}),$

with $i = 1, 2, \cdots, m$ , that is, the large Eq (2.1) is divided into $m$ subsystems. When taking

$F_i(t, \cdots) = f_i(t, \sideset{_1}{}{\mathop{x}^{(k)}}, \cdots, \sideset{_{i-1}}{}{\mathop{x}^{(k)}}, \sideset{_i}{}{\mathop{x}^{(k+1)}}, \sideset{_{i+1}}{}{\mathop{x}^{(k)}}, \cdots, \sideset{_m}{}{\mathop{x}^{(k)}}), i = 1, 2, \cdots, m$

one arrives at Guass-Jacobi WR method used frequently in actual computation, which consists of $m$ independent subsystems and is hence parallel in nature(see ^[2,6,9]). In Eq (2.2) the unique mesh $0 = t_0 < t_1 < \cdots < t_{N} = T$ is applied to all subsystems. However the subsystems of Gauss-Jocobi WR methods are independent and may have distinct behaviors, and consequently it is better to use different meshes for the subsystems with different behaviors. Applying the mesh

$0 = \sideset{_i}{_0}{\mathop{t}} < \sideset{_i}{_1}{\mathop{t}} < \cdots < \sideset{_i}{_{N_i}}{\mathop{t}} = T, \sideset{_i}{_j}{\mathop{h}} = \sideset{_i}{_{j+1}}{\mathop{t}}-\sideset{_i}{_{j}}{\mathop{t}}, j = 0, 1, \cdots, N_i-1$

to $i$ th subsystem yields the following multi-rate WR method

$\begin{equation} \begin{aligned} \sideset{_i}{^{(k+1)}_{n+1}}{\mathop{x}} = &\sideset{_i}{^{(k+1)}_{n}}{\mathop{x}}+\sideset{_i}{_n}{\mathop{h}}\Phi_i\left(\sideset{_i}{_n}{\mathop{t}}, \sideset{_i}{_n}{\mathop{h}}, \sideset{_1}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \cdots, \sideset{_{i-1}}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \sideset{_i}{^{(k+1)}_{n+1}}{\mathop{x}}, \right.\\ &\quad \sideset{_{i+1}}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \cdots, \sideset{_m}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \sideset{_1}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n}}{\mathop{t}}), \cdots, \sideset{_{i-1}}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n}}{\mathop{t}}), \\ &\quad\left. \sideset{_i}{^{(k+1)}_{n}}{\mathop{x}}, \sideset{_{i+1}}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n}}{\mathop{t}}), \cdots, \sideset{_m}{^{(k)}}{\mathop{y}}(\sideset{_i}{_{n}}{\mathop{t}})\right), n = 0, 1, \cdots, N_i-1, \end{aligned} \end{equation}$

(2.5)

where $i = 1, 2, \cdots, m, k = 0, 1, \cdots$ . Here the initial values $\sideset{_i}{^{(k+1)}_0}{\mathop{x}} = \sideset{_i}{^{(0)}_n}{\mathop{x}} = \sideset{_i}{_0}{\mathop{x}}$ for all $k$ and $n$ , $\left(\sideset{_1}{_0}{\mathop{x}}^{\rm{T}}, \sideset{_2}{_0}{\mathop{x}}^{\rm{T}}, \cdots, \sideset{_m}{_0}{\mathop{x}}^{\rm{T}} \right)^{\rm{T}} = x_0$ and $\sideset{_j}{^{(k)}}{\mathop{y}}(t)$ is the interpolation function satisfying $\sideset{_j}{^{(k)}}{\mathop{y}}(\sideset{_j}{_n}{\mathop{t}}) = \sideset{_j}{^{(k)}_n}{\mathop{x}}$ for $n = 0, 1, \cdots, N_j$ .

Because of the presence of errors in actual computing, it is necessary to consider the following perturbed systems of Eqs (2.2), (2.3) and (2.5) generated by the perturbations $\{\delta^{(k)}_n, n = 0, 1, \cdots, N, k = 0, 1, \cdots\}$

$\begin{equation} \begin{aligned} &&\tilde{x}^{(k+1)}_{n+1} = \tilde{x}^{(k+1)}_n+h_n \Phi(t_{n}, h_n, \tilde{x}^{(k+1)}_{n+1}, \tilde{x}^{(k+1)}_n, \tilde{x}^{(k)}_{n+1}, \tilde{x}^{(k)}_n)+h_n\delta^{(k+1)}_{n+1}, \\ &&n = 0, 1, \cdots, N-1;\\ &&\tilde{x}^{(k+1)}_0 = x_0+\delta^{(k+1)}_{0}, \tilde{x}^{(0)}_n = x_0+\delta^{(0)}_{n}(n = 0, 1, \cdots, N), \end{aligned} \end{equation}$

(2.6)

$\begin{equation} \begin{aligned} &\tilde{x}^{(k+1)}_{n+1} = \tilde{x}^{(k+1)}_n+h \Psi\left(t_n, h, \tilde{x}^{(k+1)}_{n+1}, \tilde{x}^{(k+1)}_{n}, \cdots, \tilde{x}^{(k+1)}_{n-p}, \right.\\ &\qquad \qquad \qquad \qquad \left. \tilde{x}^{(k)}_{n+1}, \tilde{x}^{(k)}_{n}, \cdots, \tilde{x}^{(k)}_{n-p}\right)+h\delta^{(k+1)}_{n+1}, n = p, p+1, \cdots, N-1;\\ &\tilde{x}^{(k+1)}_n = c^{(k+1)}_n+\delta^{(k+1)}_n, n = 0, 1, \cdots, p , \tilde{x}^{(0)}_n = x_0+\delta^{(0)}_{n}(n = 0, 1, \cdots, N), \end{aligned} \end{equation}$

(2.7)

and

$\begin{equation} \begin{aligned} &\sideset{_i}{^{(k+1)}_{n+1}}{\mathop{\tilde{x}}} = \sideset{_i}{^{(k+1)}_{n}}{\mathop{\tilde{x}}}+\sideset{_i}{_n}{\mathop{h}}\Phi_i\left(\sideset{_i}{_n}{\mathop{t}}, \sideset{_i}{_n}{\mathop{h}}, \sideset{_1}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \cdots, \sideset{_{i-1}}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \sideset{_i}{^{(k+1)}_{n+1}}{\mathop{\tilde{x}}}, \right.\\ &\qquad \qquad \quad \sideset{_{i+1}}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \cdots, \sideset{_m}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n+1}}{\mathop{t}}), \sideset{_1}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n}}{\mathop{t}}), \cdots, \sideset{_{i-1}}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n}}{\mathop{t}}), \\ &\qquad \qquad \quad \left.\sideset{_i}{^{(k+1)}_{n}}{\mathop{\tilde{x}}}, \sideset{_{i+1}}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n}}{\mathop{t}}), \cdots, \sideset{_m}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_i}{_{n}}{\mathop{t}})\right)+\sideset{_i}{_n}{\mathop{h}}\cdot\sideset{_i}{^{(k+1)}_{n+1}}{\mathop{\delta}}, \\ &\qquad \qquad\qquad \qquad \qquad \qquad \qquad\qquad \qquad \qquad n = 0, 1, \cdots, N_i-1;\\ &\sideset{_i}{^{(k+1)}_0}{\mathop{\tilde{x}}} = \sideset{_i}{_0}{\mathop{x}}+\sideset{_i}{^{(k+1)}_0}{\mathop{\delta}}, \sideset{_i}{^{(0)}_n}{\mathop{\tilde{x}}} = \sideset{_i}{_0}{\mathop{x}}+\sideset{_i}{^{(0)}_n}{\mathop{\delta}}(n = 0, 1, \cdots, N_i), \end{aligned} \end{equation}$

(2.8)

where $i = 1, 2, \cdots, m$ and $\sideset{_j}{^{(k)}}{\mathop{\tilde{y}}}(t)$ is the interpolation function satisfying $\sideset{_j}{^{(k)}}{\mathop{\tilde{y}}}(\sideset{_j}{_n}{\mathop{t}}) = \sideset{_j}{^{(k)}_n}{\mathop{\tilde{x}}}$ for $n = 0, 1, \cdots, N_j, j = 0, 1, \cdots, m$ .

Definition 2.1. (Page 32 of ^[10]) Let $\{\delta^{(k)}_n, n = 0, 1, \cdots, N, k = 0, 1, \cdots\}$ be any perturbation of Eq (2.2) or Eq (2.3) and Eq (2.4), and let $\{\tilde{x}^{(k)}_n, n = 0, 1, \cdots, N, k = 0, 1, \cdots\}$ be the solution of the resulting perturbed system Eq (2.6) or Eq (2.7). Then if there exist constants $C$ and $h_0$ such that, for all $h\in(0, h_0]$ ( $h = \max\limits_n h_n$ for Eq (2.2))

$\max\limits_{\begin{subarray}{c}0\leq k, 0\leq n\leq N\end{subarray}}\|\tilde{x}^{(k)}_n-x^{(k)}_n\|\leq C \varepsilon,$

whenever

$\max\limits_{\begin{subarray}{c}0\leq k, 0\leq n\leq N\end{subarray}}\|\delta^{(k)}_n\|\leq \varepsilon,$

we say that the Eq (2.2) or Eqs (2.3) and Eq (2.4) is zero-stable.

Here and hereafter $\|\cdot\|$ denotes any norm defined in $\mathbb{R}^d$ .

Definition 2.2. Let $\{_i\delta^{(k)}_n, i = 1, 2, \cdots, m, n = 0, 1, \cdots, N, k = 0, 1, \cdots\}$ be any perturbation of Eq (2.5), and let $\{_i\tilde{x}^{(k)}_n, i = 1, 2, \cdots, m, n = 0, 1, \cdots, N, k = 0, 1, \cdots\}$ be the solution of the resulting perturbed Eq (2.8). Then if there exist constants $h_0$ and $C_k$ depended on $k$ such that

$\max\limits_{\begin{subarray}{c}1\leq i\leq m\\0\leq l\leq k\end{subarray}} \max\limits_{0\leq n\leq N_i}\|_i\tilde{x}^{(l)}_n-\sideset{_i}{^{(l)}_n}{\mathop{x}}\|\leq C_k \varepsilon$

for all $h = \max\limits_{\begin{subarray}{c}1\leq i\leq m\end{subarray}} \max\limits_{0\leq n\leq N_i}\sideset{_i}{_n}{\mathop{h}}\in(0, h_0]$ , whenever

$\max\limits_{\begin{subarray}{c}1\leq i\leq m\\0\leq l\leq k\end{subarray}}\max\limits_{0\leq n\leq N_i}\|_i\delta^{(l)}_n\| \leq \varepsilon,$

we say that the Eq (2.5) is weakly zero-stable.

The following several lemmas are useful for studying zero-stability of the WR methods mentioned above.

Lemma 2.3. Let $a, b$ and $T$ be nonnegative constants, and $N$ be a positive integer. Suppose that the sets $\{h_n > 0, n = 0, 1, \cdots, N-1\}$ and $\{u_n > 0, n = 0, 1, \cdots, N\}$ satisfy

$\begin{equation} u_{n+1}\leq (1+ah_n)u_n+bh_n, n = 0, 1, \cdots, N-1. \end{equation}$

(2.9)

Then

$\max\limits_{0\leq n\leq N} u_n\leq \mathrm{e} ^{aT}u_0+\dfrac{b}{a}(\mathrm{e}^{aT}-1)$

provided that $\sum\limits_{n = 0}^{N-1}h_n = T$ .

Proof. By using Eq (2.9) repeatedly we have

$\begin{aligned} u_{n+1}\leq &(1+ah_n)(1+ah_{n-1})\cdots (1+ah_0)u_0\\ &+(1+ah_n)(1+ah_{n-1})\cdots (1+ah_1)bh_0\\ &+(1+ah_n)(1+ah_{n-1})\cdots (1+ah_2)bh_1\\ &+\cdots+(1+ah_n)bh_{n-1}+bh_n, \end{aligned}$

where $n = 0, 1, \cdots, N-1$ . Let $t_0 = 0, t_{n+1} = t_n+h_n, n = 0, 1, \cdots, N-1$ , then $t_N = T$ . Noting that $1+x < {\mathrm{e}}^x$ for $x > 0$ we can derive easily from above inequality

$\begin{aligned} u_{n+1}\leq &{\mathrm{e}}^{a(t_{n+1}-t_0)}u_0+{\mathrm{e}}^{a(t_{n+1}-t_1)}b(t_1-t_0)+{\mathrm{e}}^{a(t_{n+1}-t_2)}b(t_2-t_1)\\ &+\cdots+{\mathrm{e}}^{a(t_{n+1}-t_n)}b(t_n-t_{n-1})+{\mathrm{e}}^{a(t_{n+1}-t_{n+1})}b(t_{n+1}-t_{n})\\ \leq &{\mathrm{e}}^{a(T-t_0)}u_0+{\mathrm{e}}^{aT}b\left({\mathrm{e}}^{-at_1}(t_1-t_0)+\cdots+{\mathrm{e}}^{-at_n}(t_n-t_{n-1}) \right.\\ &\qquad \qquad \qquad \qquad \left.+\cdots+{\mathrm{e}}^{-at_N}(t_N-t_{N-1})\right) \end{aligned}$

for $n = 0, 1, \cdots, N-1$ . This together with monotonicity of the function ${\mathrm{e}}^{-at}$ imply

$\begin{aligned} u_{n+1}\leq {\mathrm{e}}^{a(T-t_0)}u_0+{\mathrm{e}}^{aT}b\int_{t_0}^T{\mathrm{e}}^{-at}\mathrm{d}t = {\mathrm{e}}^{aT}u_0+\dfrac{b}{a}\left({\mathrm{e}}^{aT}-1\right) \end{aligned}$

for all $n = 0, 1, \cdots, N-1$ .

Lemma 2.4. Let $a, b, c, d$ and $T$ be nonnegative constants, $N$ be a positive integer and $h_n, n = 0, 1, \cdots, N-1,$ be positive real numbers.

Suppose that the sequence of positive numbers $\{e_{n}^{(k)} > 0, n = 0, 1, \cdots N, k = 0, 1, \cdots\}$ satisfies

$\begin{equation} e^{(k+1)}_{n+1}\leq (1+ah_n)e^{(k+1)}_{n}+bh_ne^{(k)}_{n}+ch_ne^{(k)}_{n+1}+dh_n, n = 0, 1, \cdots N-1. \end{equation}$

(2.10)

Then for any $0 < h = \max\limits_{0\leq n\leq N-1} h_n < \dfrac{1}{2c}$

$\begin{equation} \max\limits_{\begin{subarray}{c} 0\leq k, 0\leq n\leq N\end{subarray}}e^{(k)}_n\leq \mathrm{e}^{2(a+b+c)T}\max\left\{\max\limits_{0\leq k}e^{(k)}_0, \max\limits_{0\leq n\leq N}e^{(0)}_n\right\}+ \dfrac{d}{a+b+c}(\mathrm{e}^{2(a+b+c)T}-1) \end{equation}$

(2.11)

provided that $\sum\limits_{n = 0}^{N-1}h_n = T$ .

Proof. By Eq (2.10) we have

$\begin{equation} \max\limits_{0\leq k}e^{(k+1)}_{n+1}\leq (1+ah_n)\max\limits_{0\leq k}e^{(k+1)}_{n}+bh_n\max\limits_{0\leq k}e^{(k)}_{n} +ch_n\max\limits_{0\leq k}e^{(k)}_{n+1}+dh_n \end{equation}$

(2.12)

for $n = 0, 1, \cdots, N-1$ . Let $\varepsilon = \max\limits_{0\leq n\leq N}e^{(0)}_n$ and $u_n = \max\left\{\max\limits_{0\leq k}e^{(k)}_n, \varepsilon\right\}, n = 0, 1, \cdots, N$ . Clearly $\max\left\{\max\limits_{0\leq k}e^{(k+1)}_n, \varepsilon\right\} = u_n$ . Hence we can get by Eq (2.12)

$\begin{aligned} u_{n+1}&\leq (1+ah_n)u_n+bh_nu_n+ch_nu_{n+1}+dh_n\\ &\leq \dfrac{1+ah_n+bh_n}{1-ch_n}u_n+\dfrac{dh_n}{1-ch_n}, n = 0, 1, \cdots, N-1. \end{aligned}$

Suppose that $0 < h < \dfrac{1}{2c}$ . Consequently, $\dfrac{1}{1-ch_n} < 2$ and the above inequality thus yield

$u_{n+1}\leq (1+2(a+b+c)h_n)u_n+2dh_n, n = 0, 1, \cdots, N-1.$

This together with Lemma 2.3 yield

$\max\limits_{0\leq n\leq N}u_n\leq \mathrm{e}^{2(a+b+c)T}u_0+\dfrac{d}{a+b+c}(\mathrm{e}^{2(a+b+c)T}-1).$

Thus Eq (2.11) holds true. The proof is complete.

Lemma 2.5. Let $a_j, b_j(j = 0, 1, \cdots, p), c, d$ and $T$ be nonnegative constants. Let $N$ be a positive integer and $h = T/N$ .Suppose that the sequence of positive numbers $\{e_{n}^{(k)} > 0, n = 0, 1, \cdots N, k = 0, 1, \cdots\}$ satisfies

$\begin{equation} e^{(k+1)}_{n+1}\leq e^{(k+1)}_{n}+h\sum\limits_{j = 0}^p \left(a_je^{(k+1)}_{n-j}+b_je^{(k)}_{n-j}\right)+che^{(k)}_{n+1}+dh, n = p, p+1, \cdots, N-1. \end{equation}$

(2.13)

Then for any $0 < h < \dfrac{1}{2c}$

$\max\limits_{\begin{subarray}{c} 0\leq k, p+1\leq n\leq N\end{subarray}}e^{(k)}_n\leq \mathrm{e}^{2(A+c)T}\varepsilon+\dfrac{d}{A+c}\left(\mathrm{e}^{2(A+c)T}-1\right),$

where $\varepsilon = \max\left\{ \max\limits_{0\leq n\leq N}e^{(0)}_n, \max\limits_{\begin{subarray}{c} 0\leq k, 0\leq n\leq p \end{subarray}}e^{(k)}_n\right\}$ , $A = \sum\limits_{j = 0}^p (a_j+b_j)$ .

Proof. Note that $e^{(k)}_{n}\leq \varepsilon$ for all $n\leq p$ and all $k\geq 0$ . By using Eq (2.13) repeatedly we have

$\begin{equation} \begin{aligned} e^{(k+1)}_{p+1}\leq &(1+ch+\cdots+(ch)^k)(1+Ah)\varepsilon +(ch)^{k+1}\varepsilon\\ &+(1+ch+\cdots+(ch)^k)dh, k = 0, 1, \cdots. \end{aligned} \end{equation}$

(2.14)

Let $0 < h < \dfrac{1}{2c}$ . With the notation $\alpha_h = \dfrac{1+Ah}{1-ch}$ and $\beta_h = \dfrac{dh}{1-ch}$ we can derive from Eq (2.14)

$e^{(k)}_{p+1}\leq \alpha_h\varepsilon+\beta_h, k = 0, 1, \cdots.$

Noting that $\varepsilon\leq \alpha_h\varepsilon+\beta_h$ and therefore $e^{(k)}_{n}\leq \alpha_h\varepsilon+\beta_h$ for all $n\leq p+1$ and all $k\geq 0$ we thus have by repeating the above process

$e^{(k)}_{p+2}\leq \alpha_h(\alpha_h\varepsilon+\beta_h)+\beta_h = (\alpha_h)^2\varepsilon+\beta_h(1+\alpha_h), k = 0, 1, \cdots.$

Consequently, we get

$e^{(k)}_{p+m}\leq (\alpha_h)^m\varepsilon+\beta_h(1+\alpha_h+(\alpha_h)^2+\cdots+(\alpha_h)^{m-1}), k = 0, 1, \cdots$

hold true for $m = 1, 2, \cdots, N-p$ . Note that $\alpha_h\geq 1$ . Thus

$e^{(k)}_n \leq (\alpha_h)^N\varepsilon+\beta_h(1+\alpha_h+(\alpha_h)^2+\cdots+(\alpha_h)^{N-1})$

for $k = 0, 1, \cdots$ , $n = 0, 1, \cdots, N$ . This together with $0 < ch < 1/2$ yield

$\begin{aligned} e^{(k)}_n \leq & (\alpha_h)^N\varepsilon+\dfrac{(\alpha_h)^N-1}{\alpha_h-1}\beta_h\\ \leq &(1+2(A+c)h)^{T/h}\varepsilon+\dfrac{d}{A+c}\left((1+2(A+c)h)^{T/h}-1\right)\\ \leq &\mathrm{e}^{2(A+c)T}\varepsilon+\dfrac{d}{A+c}\left(\mathrm{e}^{2(A+c)T}-1\right) \end{aligned}$

for $k = 0, 1, \cdots, n = 0, 1, \cdots, N$ . The proof is complete.

3. Main results

Theorem 3.1. Suppose that there exist constants $L_1$ , $L_2$ , $L_3$ and $L_4$ such that

$\begin{equation} \begin{aligned} &\|\Phi(t, h, x_1, x_2, x_3, x_4)-\Phi(t, h, y_1, y_2, y_3, y_4)\|\leq \sum\limits_{i = 1}^{4} L_i\|x_i-y_i\| \end{aligned} \end{equation}$

(3.1)

for any real numbers $t, h\in \mathbb{R}, x_i, y_i\in \mathbb{R}^{d}, i = 1, 2, 3, 4$ . Then the WR Eq (2.2) is zero-stable.

Proof. Let $\varepsilon^{(k)}_n = \|\tilde{x}^{(k)}_n-x^{(k)}_n\|$ and for any $\varepsilon > 0$ the perturbations in Eq (2.6) $\{\delta^{(k)}_n, k = 0, 1, \cdots, n = 0, 1, \cdots, N\}$ satisfies $\max\limits_{\begin{subarray}{c}0\leq k, 0\leq n\leq N\end{subarray}}\|\delta^{(k)}_n\| < \varepsilon$ . Then by virtue of Eqs (2.2) and (2.6) and the condition (3.1) we have

$\begin{equation} \varepsilon^{(k+1)}_{n+1}\leq \varepsilon^{(k+1)}_n +h_n\left( L_1\varepsilon^{(k+1)}_{n+1}+L_2 \varepsilon^{(k+1)}_n+L_3 \varepsilon^{(k)}_{n+1}+L_4\varepsilon^{(k)}_n\right)+h_n\varepsilon \end{equation}$

(3.2)

for $n = 0, 1, \cdots, N-1, k = 0, 1, \cdots$ and

$\begin{equation} \varepsilon^{(k+1)}_0\leq \varepsilon, k = 0, 1, \cdots, \varepsilon^{(0)}_n\leq \varepsilon, n = 0, 1, \cdots, N. \end{equation}$

(3.3)

Let $h_n < \dfrac{1}{2L_1}$ . Then $\dfrac{1}{1-h_nL_1} < 2$ . We therefore derive from Eq (3.2)

$\begin{equation} \varepsilon^{(k+1)}_{n+1}\leq (1+2h_n(L_1+L_2))\varepsilon^{(k+1)}_{n}+2h_nL_3\varepsilon^{(k)}_{n+1}+2h_nL_4\varepsilon^{(k)}_{n}+2h_n\varepsilon \end{equation}$

(3.4)

for $n = 0, 1, \cdots, N-1, k = 0, 1, \cdots$ . Consequently, by Eqa (3.3), (3.4) and Lemma 2.4 we have

$\varepsilon^{(k)}_{n}\leq \mathrm{e}^{4(L_1+L_2+L_3+L_4)T}\varepsilon+\dfrac{\varepsilon}{L_1+L_2+L_3+L_4}\left(\mathrm{e}^{4(L_1+L_2+L_3+L_4)T}-1\right)$

for any $k = 0, 1, \cdots, n = 0, 1, \cdots, N$ if $0 < h < \min\left\{\dfrac{1}{2L_1}, \dfrac{1}{4L_3}\right\}$ . The proof is therefore complete.

Theorem 3.2. Suppose that there exist constants $L_i, i = 1, 2, \cdots, 2p+4$ such that

$\begin{equation} \begin{aligned} &\|\Psi(t, h, x_1, x_2, \cdots, x_{2p+4})-\Psi(t, h, y_1, y_2, \cdots, y_{2p+4})\|\leq \sum\limits_{i = 1}^{2p+4}L_i\|x_i-y_i\| \end{aligned} \end{equation}$

(3.5)

for any real numbers $t, h\in \mathbb{R}, x_i, y_i\in \mathbb{R}^d, i = 1, 2, \cdots, 2p+4$ . Then the WR Eq (2.3) and Eq (2.4) is zero-stable.

Proof. Let $\varepsilon^{(k)}_n$ denote $\|\tilde{x}^{(k)}_n-x^{(k)}_n\|$ and the perturbations $\{\delta^{(k)}_n, k = 0, 1, \cdots, n = 0, 1, \cdots, N\}$ in Eq (2.7) satisfy $\max\limits_{\begin{subarray}{c}0\leq k, 0\leq n\leq N\end{subarray}}\|\delta^{(k)}_n\| < \varepsilon$ , where $\varepsilon$ is any positive number. Then by virtue of Eqs (2.3), (2.4), (2.7) and the condition (3.5) we have

$\begin{equation} \begin{aligned} \varepsilon^{(k+1)}_{n+1}\leq & \varepsilon^{(k+1)}_n+h(L_1\varepsilon^{(k+1)}_{n+1}+L_2\varepsilon^{(k+1)}_{n}+\cdots+L_{p+2}\varepsilon^{(k+1)}_{n-p}\\ &+L_{p+3}\varepsilon^{(k)}_{n+1}+L_{p+4}\varepsilon^{(k)}_{n}+\cdots+L_{2p+4}\varepsilon^{(k)}_{n-p})+h\varepsilon, \\ &\qquad\qquad\qquad\qquad\qquad\quad n = p, p+1, \cdots, N-1, \\ \varepsilon^{(k+1)}_{n}\leq & \varepsilon, n = 0, 1, \cdots, p, \varepsilon^{(0)}_n\leq \varepsilon, n = 0, 1, \cdots, N, \end{aligned} \end{equation}$

(3.6)

where $k = 0, 1, \cdots$ . Let $h < \dfrac{1}{2L_1}$ . Consequently, $\dfrac{1}{1-hL_1} < 2$ . We therefore derive from Eq (3.6)

$\begin{equation} \begin{aligned} \varepsilon^{(k+1)}_{n+1}\leq & (1+2h(L_1+L_2))\varepsilon^{(k+1)}_n+2h(L_3\varepsilon^{(k+1)}_{n-1}+\cdots+L_{p+2}\varepsilon^{(k+1)}_{n-p}+L_{p+3}\varepsilon^{(k)}_{n+1}\\ &+L_{p+4}\varepsilon^{(k)}_{n}+\cdots+L_{2p+4}\varepsilon^{(k)}_{n-p})+2h\varepsilon, n = p, p+1, \cdots, N-1, \\ \varepsilon^{(k+1)}_{n}\leq & \varepsilon, n = 0, 1, \cdots, p, \varepsilon^{(0)}_n\leq \varepsilon, n = 0, 1, \cdots, N \end{aligned} \end{equation}$

(3.7)

for $n = 0, 1, \cdots, N-1, k = 0, 1, \cdots$ . By Eq (3.7) and Lemma 2.5 we have

$\varepsilon^{(k)}_{n}\leq \mathrm{e}^{4(L_1+L_2+\cdots+L_{2p+4})T}\varepsilon+\dfrac{\varepsilon}{L_1+L_2+\cdots+L_{2p+4}}\left(\mathrm{e}^{4(L_1+L_2+\cdots+L_{2p+4})T}-1\right)$

for any $k = 0, 1, \cdots, n = 0, 1, \cdots, N$ if $0 < h < \min\left\{\dfrac{1}{2L_1}, \dfrac{1}{4L_{p+3}}\right\}$ . The proof is complete.

Theorem 3.3. For Eqs (2.5) and (2.8) let $\Phi = (\Phi_1^{\mathrm{T}}, \Phi_2^{\mathrm{T}}, \cdots, \Phi_m^{\mathrm{T}})^{\mathrm{T}}$ satisfy the Lipschitz condition:there exist constants $L_1$ and $L_2$ such that

$\begin{equation} \|\Phi(t, h, x_1, x_2)-\Phi(t, h, y_1, y_2)\|\leq \sum\limits_{i = 1}^{2} L_i\|x_i-y_i\|, \;for\;any\; x_1, x_2, y_1, y_2\in \mathbb{R}^d, \end{equation}$

(3.8)

and the interpolation functions $\sideset{_j}{^{(k)}}{\mathop{\tilde{y}}}(t)$ and $\sideset{_j}{^{(k)}}{\mathop{y}}(t)$ satisfy that

$\begin{equation} \sup\limits_{0\leq t\leq T}\|\sideset{_j}{^{(k)}}{\mathop{\tilde{y}}}(t)- \sideset{_j}{^{(k)}}{\mathop{y}}(t)\|\leq C\max\limits_{0\leq n\leq N_j}\|\sideset{_j}{^{(k)}_n}{\mathop{\tilde{x}}}-\sideset{_j}{^{(k)}_n}{\mathop{x}}\|, j = 1, 2, \cdots, m, \end{equation}$

(3.9)

where $C$ is the positive constant. Then the multi-rate Eq (2.5) is weakly zero-stable.

Proof. Define $\|x\|_a = \max\limits_{1\leq i\leq m}\|_ix\|$ , $\|x\|_b = \sum\limits_{i = 1}^m\|_ix\|$ for $x = (\sideset{_1}{^{\rm{T}}}{\mathop{x}}, \sideset{_2}{^{\rm{T}}}{\mathop{x}}, \cdots, \sideset{_m}{^{\rm{T}}}{\mathop{x}})^{\rm{T}}\in \mathbb{R}^d$ , where $_ix\in \mathbb{R}^{d_i}, i = 1, 2, \cdots, m, d_1+d_2+\cdots+d_m = d$ . It is easy to prove that $\|\cdot\|_a$ and $\|\cdot\|_b$ are the norm. Consequently, there exist constants $C_1$ and $C_2$ such that

$\begin{equation} C_1\|x\|_a\leq \|x\|\leq C_2\|x\|_b \;{ \rm{ for\;any } }\;x\in \mathbb{R}^d. \end{equation}$

(3.10)

Eqs (3.8) and (3.10) thus yield

$\begin{aligned} &C_1\|\Phi(t, h, x_1, x_2)- \Phi(t, h, y_1, y_2)\|_a \leq \|\Phi(t, h, x_1, x_2)- \Phi(t, h, y_1, y_2)\| \\ \leq &\sum\limits_{j = 1}^2 L_j\|x_j-y_j\|\leq \sum\limits_{j = 1}^2 L_jC_2\|x_j-y_j\|_b = \sum\limits_{j = 1}^2 L_jC_2\sum\limits_{i = 1}^m\|_ix_j-\sideset{_i}{_j}{\mathop{y}}\| \end{aligned}$

for any $x_1, x_2, y_1, y_2\in \mathbb{R}^d$ and therefore

$\begin{equation} \|\Phi_i(t, h, x_1, x_2)-\Phi_i(t, h, y_1, y_2)\|\leq \sum\limits_{j = 1}^2L_j\dfrac{C_2}{C_1}\sum\limits_{i = 1}^m\|_ix_j-\sideset{_i}{_j}{\mathop{y}}\|. \end{equation}$

(3.11)

With the notation $_ie^{(k)}_n = \|_i\tilde{x}^{(k)}_n-\sideset{_i}{^{(k)}_n}{\mathop{x}}\|$ we can derive from Eqs (2.5), (2.8), (3.9) and (3.11)

$\begin{equation} \begin{aligned} _ie^{(k+1)}_{n+1}\leq &\sideset{_i}{^{(k+1)}_n}{\mathop{e}}+_ih_n \dfrac{L_1C_2C}{C_1}\left(\max\limits_{0\leq j\leq N} \sideset{_1}{^{(k)}_j}{\mathop{e}}+\cdots+ \max\limits_{0\leq j\leq N_i}\sideset{_{i-1}}{^{(k)}_j}{\mathop{e}}+\sideset{_i}{^{(k+1)}_{n+1}}{\mathop{e}}\right.\\ &+\left.\max\limits_{0\leq j\leq N_i} \sideset{_{i+1}}{^{(k)}_j}{\mathop{e}}+\cdots+\max\limits_{0\leq j\leq N_i}\sideset{_m}{^{(k)}_j}{\mathop{e}}\right)+ \sideset{_i}{_n}{\mathop{h}} \dfrac{L_2C_2C}{C_1}\left(\max\limits_{0\leq j\leq N_i} \sideset{_1}{^{(k)}_j}{\mathop{e}}\right.\\ &\left.+\cdots+\max\limits_{0\leq j\leq N_i}\sideset{_{i-1}}{^{(k)}_j}{\mathop{e}}+\sideset{_i}{^{(k+1)}_{n}} {\mathop{e}}+\max\limits_{0\leq j\leq N_i} \sideset{_{i+1}}{^{(k)}_j}{\mathop{e}}+\cdots+\max\limits_{0\leq j\leq N_i}\sideset{_m}{^{(k)}_j}{\mathop{e}}\right)\\ &+\sideset{_i}{_n}{\mathop{h}} \|\sideset{_i}{^{(k+1)}_{n+1}}{\mathop{\delta}}\|. \end{aligned} \end{equation}$

(3.12)

Let $h = \max\limits_{\begin{subarray}{c}1\leq i\leq m \end{subarray}}\max\limits_{0\leq n\leq N_i} \sideset{_i}{_n}{\mathop{h}} < \dfrac{C_1}{2L_1C_2C}$ . Consequently, $_ih_n\dfrac{L_1C_2C}{C_1} < \dfrac12 < 1-\sideset{_i}{_n}{\mathop{h}}\dfrac{L_1C_2C}{C_1}$ which together with Eq (3.12) imply that

$\begin{equation} \begin{aligned} _ie^{(k+1)}_{n+1}\leq &\left(1+a \sideset{_i}{_n}{\mathop{h}}\right) \sideset{_i}{^{(k+1)}_n}{\mathop{e}}+\left(a\sum\limits_{\begin{subarray}{c}l = 1\\ l\neq i\end{subarray}}^m\max\limits_{0\leq j\leq N_i} \sideset{_l}{^{(k)}_j}{\mathop{e}}+2\|\sideset{_i}{^{(k+1)}_{n+1}}{\mathop{\delta}}\|\right)\sideset{_i}{_n}{\mathop{h}}\\ \leq &(1+a \sideset{_i}{_n}{\mathop{h}}) \sideset{_i}{^{(k+1)}_n}{\mathop{e}}+ \left(a\sum\limits_{i = 1}^m\max\limits_{0\leq j\leq N_i} \sideset{_i}{_j^{(k)}}{\mathop{e}}+2\max\limits_{\begin{subarray}{c}1\leq i\leq m \end{subarray}} \max\limits_{0\leq j\leq N_i}\|\sideset{_i}{_j^{(k+1)}}{\mathop{\delta}}\| \right)\sideset{_i}{_n}{\mathop{h}} \end{aligned} \end{equation}$

(3.13)

with $a = \dfrac{2(L_1+L_2)C_2C}{C_1}$ . By using Eq (3.13) with $k = 0$ and Lemma 2.3 we get

$\begin{aligned} _ie^{(1)}_n\leq &\mathrm{e}^{aT}\sideset{_i}{^{(1)}_0}{\mathop{e}}+\dfrac{1}{a}({\mathrm{e}}^{aT}-1) \left(a\sum\limits_{i = 1}^m\max\limits_{0\leq j\leq N_i} \sideset{_i}{_j^{(0)}}{\mathop{e}}+2\max\limits_{\begin{subarray}{c}1\leq i\leq m\end{subarray}} \max\limits_{0\leq j\leq N_i}\|\sideset{_i}{_j^{(1)}}{\mathop{\delta}}\| \right). \end{aligned}$

Thus

$\begin{aligned} \max\limits_{\begin{subarray}{c}1\leq i\leq m\end{subarray}}\max\limits_{0\leq n\leq N_i}\sideset{_i}{^{(1)}_n}{\mathop{e}} \leq a_1\max\left(\max\limits_{1\leq i\leq m} \sideset{_i}{^{(1)}_0}{\mathop{e}}, \sum\limits_{i = 1}^m\max\limits_{0\leq j\leq N_i}\sideset{_i}{_j^{(0)}}{\mathop{e}}, \max\limits_{\begin{subarray}{c}1\leq i\leq m\end{subarray}}\max\limits_{0\leq j\leq N_i}\|\sideset{_i}{_j^{(1)}}{\mathop{\delta}}\| \right) \end{aligned}$

with $a_1 = {\mathrm{e}}^{aT}+\dfrac{1}{a}({\mathrm{e}}^{aT}-1)(a+2)$ . Similarly, using the above inequality, Eq (3.13) with $k = 1$ and Lemma 2.3 we have

$\begin{aligned} \max\limits_{\begin{subarray}{c}1\leq i\leq m\end{subarray}}\max\limits_{0\leq n\leq N_i}\sideset{_i}{^{(2)}_n}{\mathop{e}} \leq &a_2\max\left(\max\limits_{ \begin{subarray}{c}1\leq i\leq m\\ 1\leq l\leq 2\end{subarray}} \sideset{_i}{^{(l)}_0}{\mathop{e}}, \sum\limits_{i = 1}^m\max\limits_{0\leq j\leq N_i}\sideset{_i}{_j^{(0)}}{\mathop{e}}, \max\limits_{\begin{subarray}{c}1\leq i\leq m\\1\leq l\leq 2\end{subarray}} \max\limits_{0\leq j\leq N_i}\|\sideset{_i}{_j^{(l)}}{\mathop{\delta}}\| \right) \end{aligned}$

with $a_2 = {\mathrm{e}}^{aT}+({\mathrm{e}}^{aT}-1)(ma_1+\dfrac{2}{a})$ . Repeating the above process yields that for any $k$ there exists the constant $a_k$ such that

$\begin{aligned} &\max\limits_{\begin{subarray}{c}1\leq i\leq m\\ 0\leq l\leq k \end{subarray}} \max\limits_{0\leq n\leq N_i}\|\sideset{_i}{^{(l)}_n}{\mathop{\tilde{x}}}-\sideset{_i}{^{(l)}_n}{\mathop{x}}\| \\ \leq & a_k\max\left(\max\limits_{\begin{subarray}{c}1\leq i\leq m\\ 0\leq l\leq k\\ \end{subarray}}\|\sideset{_i}{^{(l)}_0}{\mathop{\delta}}\|, \sum\limits_{i = 1}^m \max\limits_{0\leq n\leq N_i}\|\sideset{_i}{^{(0)}_n}{\mathop{\delta}}\|, \max\limits_{\begin{subarray}{c}1\leq i\leq m\\0\leq l\leq k \end{subarray}}\max\limits_{0\leq n\leq N_i}\|\sideset{_i}{^{(l)}_n}{\mathop{\delta}}\| \right) \end{aligned}$

for any mesh $\{_it_n, n = 0, 1, \cdots, N, i = 1, 2, \cdots, m\}$ . The proof is therefore complete.

Remark 1. In the above discussions we assume that the perturbation at one step is $\delta_n^{(k)} h$ , which depends on step-size $h$ . This assumption is necessary for the proof of the theorems mentioned above. However it seems to be more reasonable to regard $\delta_n^{(k)}$ , a small number, as the perturbation, if the step-size is very small. In this case the tiny perturbation at each step may bring out the huge change of the numerical solution. We will explore these in next section by some computer experiments.

Remark 2. Note that the convergence property of the WR methods is not be used in the proof of Theorem 3.1-3.3. Hence, a divergent WR method may be zero-stable. In other words, zero-stability does not imply convergence.

4. Numerical experiments

In this section some numerical experiments are performed to verify the theorems obtained in Section 3.

Consider the equation

$\begin{equation} x' = Ax = \left( \begin{array}{ccc} T_1 & I & 0 \\ I & T_2 & I \\ 0 & I & T_3 \\ \end{array} \right)x, t\in [0, 1], x(0) = x_0, \end{equation}$

(4.1)

where $x\in \mathbb{R}^{15}$ , $T_i = \left(\begin{array}{cccc} -2^i & 1 & & \\ 1 & -2^i & \ddots & \\ & \ddots & \ddots & 1\\ & & 1 & -2^i\\ \end{array} \right)\in \mathbb{R}^{5\times 5}$ and $I$ is the identity matrix. The special case of such system has been examined by Burrage ^[11], which is obtained by discretizing the heat equation in two variables.

Let $A_1 = \left(\begin{array}{ccc} T_1 & & \\ & T_2 & \\ & & T_3 \\ \end{array} \right)$ and $A_2 = A-A_1$ . First, choose the following WR methods of solving the Eq (4.1):

(i) Euler method with variable step-size

$\begin{aligned} &x^{(k+1)}_{n+1} = x^{(k+1)}_n+h_n(A_1x^{(k+1)}_n+A_2x^{(k)}_n), n = 0, 1, \cdots, N-1, k = 0, 1, \cdots, \\ &x^{(k+1)}_0 = x_0\;{ \rm{ for\;all } }\;k, x^{(0)}_n = x_0\;{ \rm{ for\;all } }\;n. \end{aligned}$

(ii) Two step methods with fixed step-size

$\begin{aligned} &x^{(k+1)}_{n+1} = x^{(k+1)}_n+\dfrac{h}{2}\left[3(A_1x^{(k+1)}_n+A_2x^{(k)}_n)-(A_1x^{(k+1)}_{n-1}+A_2x^{(k)}_{n-1})\right], \\ &\qquad \qquad \qquad \qquad \qquad \qquad n = 1, 2, \cdots, N-1, k = 0, 1, \cdots, \\ &x^{(k+1)}_0 = x_0\;{ \rm{ for\;all } }\;k, x^{(0)}_n = x_0\;{ \rm{ for\;all } }\;n, \end{aligned}$

where $x^{(k+1)}_1$ is obtained by a suitable one-step method. Here let

$x^{(k+1)}_1 = x_0+hAx_0+\dfrac{h^2}{2}A^2x_0.$

(iii) Multi-rate Jacobi methods based on Euler method

$\left. \begin{aligned} &&_ix^{(k+1)}_{n+1} = \sideset{_i}{_{n}^{(k+1)}}{\mathop{x}}+ \sideset{_i}{_{n}}{\mathop{h}}\left(T_i\cdot \sideset{_i}{_{n+1}^{(k+1)}}{\mathop{x}} +\sideset{_{i-1}}{_{n+1}^{(k)}}{\mathop{x}}+\sideset{_{i+1}}{_{n+1}^{(k)}}{\mathop{x}}\right), \\ && n = 0, 1, \cdots, N_i-1, \\ &&\sideset{_i}{_{0}^{(k+1)}}{\mathop{x}} = \sideset{_i}{_{n}^{(0)}}{\mathop{x}} = \sideset{_i}{_{0}}{\mathop{x}}, \sideset{_0}{_{n}^{(k)}}{\mathop{x}} = \sideset{_4}{_{n}^{(k)}}{\mathop{x}} = (0, 0, 0, 0, 0)^{\rm{T}}\;{ \rm{ for\;all} }\;k\;{ \rm{ and } }\;n \end{aligned}\right\} i = 1, 2, 3.$

Second, assume the initial values and the approximation solution at each calculation step have tiny perturbations as in Eqs (2.6)–(2.8) and let $\delta^{(k)}_n\equiv \Delta = \delta\cdot (1, 1, \cdots, 1)^{\rm{T}}\in \mathbb{R}^{15}$ , $_i\delta^{(k)}_n\equiv \Delta_p = \delta\cdot (1, 1, 1, 1, 1)^{\rm{T}}\in \mathbb{R}^5$ , $\tilde{x}^{(k)}_n$ denote the resulting perturbed solution and $u^{(k)}_n$ denote $\tilde{x}^{(k)}_n-x^{(k)}_n$ , which yields

$\begin{equation} \begin{aligned} &u^{(k+1)}_{n+1} = u^{(k+1)}_n+h_n(A_1u^{(k+1)}_n+A_2u^{(k)}_n)+h_n\Delta, n = 0, 1, \cdots, N-1, k = 0, 1, \cdots, \\ &u^{(k+1)}_0 = \Delta\;{ \rm{ for\;all } }\;k, u^{(0)}_n = \Delta\;{ \rm{ for\;all } }\;n, \end{aligned} \end{equation}$

(4.2)

$\begin{equation} \begin{aligned} &u^{(k+1)}_{n+1} = u^{(k+1)}_n+\dfrac{h}{2}\left[3(A_1u^{(k+1)}_n+A_2u^{(k)}_n)-(A_1u^{(k+1)}_{n-1}+A_2u^{(k)}_{n-1})\right]+\dfrac{h}{2}\Delta, \\ &\qquad \qquad \qquad \qquad \qquad \qquad n = 1, 2, \cdots, N-1, k = 0, 1, \cdots, \\ &u^{(k+1)}_0 = \Delta, u^{(k+1)}_1 = (I+hA+\dfrac{h^2}{2}A^2)\Delta, \;{ \rm{ for\;all } }\;k, u^{(0)}_n = \Delta\;{ \rm{ for\;all } }\;n, \end{aligned} \end{equation}$

(4.3)

and

$\begin{array}{l} \left. \begin{array}{r} _iu^{(k+1)}_{n+1} = \sideset{_i}{_{n}^{(k+1)}}{\mathop{u}}+ \sideset{_i}{_{n}}{\mathop{h}}\left(T_i\cdot \sideset{_i}{_{n+1}^{(k+1)}}{\mathop{u}} +\sideset{_{i-1}}{_{n+1}^{(k)}}{\mathop{u}}+\sideset{_{i+1}}{_{n+1}^{(k)}}{\mathop{u}}+\Delta_p\right), \\ n = 0, 1, \cdots, N_i-1, \\ \sideset{_i}{_{0}^{(k+1)}}{\mathop{u}} = \sideset{_i}{_{n}^{(0)}}{\mathop{u}} = \Delta_p, \sideset{_0}{_{n}^{(k)}}{\mathop{u}} = \sideset{_4}{_{n}^{(k)}}{\mathop{u}} = (0, 0, 0, 0, 0)^{\rm{T}}, \forall k, \forall n \end{array}\right\} i = 1, 2, 3. \end{array}$

(4.4)

Lastly, reach a reasonable conclusion by analysing data generated by Eqs (4.2)–(4.4).

We choose $A$ , $A_1$ , $A_2$ and $\delta^{(k)}_n$ as above to guarantee that $u^{(k)}_n$ increase with $k$ and $n$ that is the worst case of error propagation.

In our experiments the step-sizes $h_n$ in Eq (4.2) and $_ih_n$ in Eq (4.4) are taken as

$h_n = \min\left\{10h, \max\left\{0.1h, |\xi_n|h\right\}\right\}$

and

$_ih_n = \min\{10^ih, \max\{10^{i-2}h, 10^{i-1}|\xi_n|h\}\},$

where $h > 0$ , $\{\xi_0, \xi_1, \cdots\}$ is a sequence of independent and identically distributed random variables satisfying $\xi_i\sim N(0, 1)$ .

The data in , and show clearly the errors resulted by the perturbations $h_n\Delta(h\Delta/2)$ in Eq (4.2) (Eq (4.3) are controllable as $k$ and $n$ tend to infinity which support the theorems developed in Section 3). However the data in and show that the errors increase with $k$ and $1/h$ when the perturbations $h_n\Delta (h\Delta/2)$ are replaced by the constant perturbation $\Delta(\Delta/2)$ .

Table 3. The computing results of Eq (4.2) with

$\delta = 0.0001$ (

$\varepsilon = \max\limits_{0\leq l\leq k}\max\limits_{0\leq n\leq N}u^{(l)}_n$ ).

$h$	0.1		0.001
$k$	5	10	10	15
$\varepsilon$	3.1483e-04	3.1483e-04	3.1605e-04	3.1605e-04
$h$	0.00001
$k$	10	15
$\varepsilon$	3.1606e-04	3.1606e-04

| Show Table

DownLoad: CSV

Table 4. The computing results of Eq (4.2) with

$\delta = 0.0001$ (

$\varepsilon = \max\limits_{0\leq l\leq k}\max\limits_{0\leq n\leq N}u^{(l)}_n$ ) for the case of replacing perturbation

$h_n \Delta$ with the constant perturbation

$\Delta$ .

$h$	0.1		0.001		0.00001
$k$	5	10	10	15	10	15
$\varepsilon$	0.0018	0.0018	0.1824	0.1824	17.6723	17.6723
$h$	0.000001
$k$	10	15
$\varepsilon$	176.8079	176.8079

| Show Table

DownLoad: CSV

Table 5. The computing results of Eq (4.3) with

$\delta = 0.0001$ (

$\varepsilon = \max\limits_{0\leq l\leq k}\max\limits_{0\leq n\leq N}u^{(l)}_n$ ).

$h$	0.1		0.001
$k$	10	15	10	15
$\varepsilon$	3.1568e-04	3.1568e-04	3.1606e-04	3.1606e-04
$h$	0.00001
$k$	10	15
$\varepsilon$	3.1606e-04	3.1606e-04

| Show Table

DownLoad: CSV

Table 6. The computing results of Eq (4.3) with

$\delta = 0.0001$ (

$\varepsilon = \max\limits_{0\leq l\leq k}\max\limits_{0\leq n\leq N}u^{(l)}_n$ ) for the case of replacing perturbation

$h \Delta/2$ with

$\Delta/2$ .

$h$	0.1		0.001		0.00001
$k$	10	15	10	15	10	15
$\varepsilon$	0.1782	0.1782	0.1435	0.1435	14.1724	14.1724
$h$	0.000001
$k$	10	15
$\varepsilon$	141.7252	141.7252

| Show Table

DownLoad: CSV

Table 7. The computing results of Eq (4.4) with

$\delta = 0.0001$ (

$\varepsilon = \max\limits_{0\leq l\leq k}\max\limits_{0\leq n\leq N}u^{(l)}_n$ ).

$h$	0.1		0.001
$k$	5	10	10	15
$\varepsilon$	3.1085e-04	3.1085e-04	3.1567e-04	3.1567e-04
$h$	0.00001
$k$	10	15
$\varepsilon$	3.1606e-04	3.1606e-04

| Show Table

DownLoad: CSV

5. Further work

It is well known that a linear multi-step method with fixed step-size is convergent if, and only if, it is both consistent and zero-stable. In this paper we only present the sufficient conditions of zero-stability for some special methods. We will explore the relationship between convergence and zero-stability for WR methods based on general linear multi-step methods in the future. Moreover, we shall apply our methods to the numerical study of some inverse problems ^[12,13,14].

Acknowledgments

Research is supported by the Natural Science Foundation of Fujian Province, China (grant number 2021J011031).

Conflict of interest

The author declares that he has no conflicts of interest regarding the publication of this paper.

References

[1]	Statista, Transport volume of crude oil in seaborne trade worldwide from 2010 to 2021, 2023. Available from: https://www.statista.com/statistics/264013/transport-volume-of-crude-oil-in-seaborne-trade/.
[2]	Statista, Oil production worldwide from 1998 to 2021, 2023. Available from: https://www.statista.com/statistics/265229/global-oil-production-in-million-metric-tons/.
[3]	J. Qi, S. Wang, H. Psaraftis, Bi-level optimization model applications in managing air emissions from ships: a review, Commun. Transp. Res., 1 (2021), 100020. https://doi.org/10.1016/j.commtr.2021.100020 doi: 10.1016/j.commtr.2021.100020
[4]	R. Yan, S. Wang, L. Zhen, G. Laporte, Emerging approaches applied to maritime transport research: Past and future, Commun. Transp. Res., 1 (2021), 100011. https://doi.org/10.1016/j.commtr.2021.100011 doi: 10.1016/j.commtr.2021.100011
[5]	U.S. Energy Information Administration (USEIA), Oil and petroleum products explained, 2022. Available from: https://www.eia.gov/energyexplained/oil-and-petroleum-products/where-our-oil-comes-from.php.
[6]	Statista, Largest importers of Russian crude oil and petroleum motor fuels in the European Union and the United Kingdom in 2021, by country, 2023. Available from: https://www.statista.com/statistics/1295165/oil-import-value-from-russia-to-eu-countries/.
[7]	European Council (EC), EU sanctions against Russia explained, 2023. Available from: https://www.consilium.europa.eu/en/policies/sanctions/restrictive-measures-against-russia-over-ukraine/sanctions-against-russia-explained.
[8]	A. Almendral, Russia is amassing a shadow fleet of tankers to avoid EU oil sanctions, 2022. Available from: https://qz.com/russia-is-amassing-a-shadow-fleet-of-tankers-to-avoid-e-1849853809.
[9]	X. Bai, X. Zhang, K. X. Li, Y. Zhou, K. F. Yuen, Research topics and trends in the maritime transport: A structural topic model, Transp. Policy, 102 (2021), 11–24. https://doi.org/10.1016/j.tranpol.2020.12.013 doi: 10.1016/j.tranpol.2020.12.013
[10]	H. Zhao, Q. Meng, Y. Wang, Robust container slot allocation with uncertain demand for liner shipping services, Flexible Serv. Manuf. J., 34 (2022), 551–579. https://doi.org/10.1007/s10696-021-09420-z doi: 10.1007/s10696-021-09420-z
[11]	Q. Meng, S. Wang, H. Andersson, K. Thun, Containership routing and scheduling in liner shipping: overview and future research directions, Transp. Sci., 48 (2014), 265–280. https://doi.org/10.1287/trsc.2013.0461 doi: 10.1287/trsc.2013.0461
[12]	S. Wang, Q. Meng, Container liner fleet deployment: a systematic overview, Transp. Res. Part C Emerg. Technol., 77 (2017), 389–404. https://doi.org/10.1016/j.trc.2017.02.010 doi: 10.1016/j.trc.2017.02.010
[13]	A. Baykasoğlu, K. Subulan, A. S. Taşan, N. Dudaklı, A review of fleet planning problems in single and multimodal transportation systems, Transportmetrica A : Transport Sci., 15 (2019), 631–697. https://doi.org/10.1080/23249935.2018.1523249 doi: 10.1080/23249935.2018.1523249
[14]	M. Christiansen, E. Hellsten, D. Pisinger, D. Sacramento, C. Vilhelmsen, Liner shipping network design, Eur. J. Oper. Res., 286 (2020), 1–20. https://doi.org/10.1016/j.ejor.2019.09.057 doi: 10.1016/j.ejor.2019.09.057
[15]	M. A. Dulebenets, J. Pasha, O. F. Abioye, M. Kavoosi, Vessel scheduling in liner shipping: a critical literature review and future research needs, Flexible Serv. Manuf. J., 33 (2021), 43–106. https://doi.org/10.1007/s10696-019-09367-2 doi: 10.1007/s10696-019-09367-2
[16]	D. Kizilay, D. T. Eliiyi, A comprehensive review of quay crane scheduling, yard operations and integrations thereof in container terminals, Flexible Serv. Manuf. J., 33 (2021), 1–42. https://doi.org/10.1007/s10696-020-09385-5 doi: 10.1007/s10696-020-09385-5
[17]	Q. Meng, T. Wang, A chance constrained programming model for short-term liner ship fleet planning problems, Marit. Policy Manage., 37 (2010), 329–346. https://doi.org/10.1080/03088839.2010.486635 doi: 10.1080/03088839.2010.486635
[18]	M. Ng, D. Y. Lin, Fleet deployment in liner shipping with incomplete demand information, Transp. Res. Part E Logist. Transp. Rev., 116 (2018) 184–189. https://doi.org/10.1016/j.tre.2018.06.004 doi: 10.1016/j.tre.2018.06.004
[19]	S. Wang, Q. Meng, Liner ship fleet deployment with container transshipment operations, Transp. Res. Part E Logist. Transp. Rev., 48 (2012), 470–484. https://doi.org/10.1016/j.tre.2011.10.011 doi: 10.1016/j.tre.2011.10.011
[20]	Q. Meng, T. Wang, S. Wang, Short-term liner ship fleet planning with container transshipment and uncertain container shipment demand, Eur. J. Oper. Res., 223 (2012), 96–105. https://doi.org/10.1016/j.ejor.2012.06.025 doi: 10.1016/j.ejor.2012.06.025
[21]	T. Wang, Q. Meng, S. Wang, Z. Tan, Risk management in liner ship fleet deployment: a joint chance constrained programming model, Transp. Res. Part E Logist. Transp. Rev., 60 (2013), 1–12. https://doi.org/10.1016/j.tre.2013.09.001 doi: 10.1016/j.tre.2013.09.001
[22]	M. Ng, Distribution-free vessel deployment for liner shipping, Eur. J. Oper. Res., 238 (2014), 858–862. https://doi.org/10.1016/j.ejor.2014.04.019 doi: 10.1016/j.ejor.2014.04.019
[23]	M. Ng, Container vessel fleet deployment for liner shipping with stochastic dependencies in shipping demand, Transp. Res. Part B Methodol., 74 (2015), 79–87. https://doi.org/10.1016/j.trb.2015.01.004 doi: 10.1016/j.trb.2015.01.004
[24]	S. Wang, Z. Liu, X. Qu, Minimax regret model for liner shipping fleet deployment with uncertain demand, Transp. Res. Rec., 2549 (2016), 45–53. https://doi.org/10.3141/2549-06 doi: 10.3141/2549-06
[25]	Q. Meng, S. Wang, Optimal operating strategy for a long-haul liner service route, Eur. J. Oper. Res., 215 (2011), 105–114. https://doi.org/10.1016/j.ejor.2011.05.057 doi: 10.1016/j.ejor.2011.05.057
[26]	L. Zhen, Y. Hu, S. Wang, G. Laporte, Y. Wu, Fleet deployment and demand fulfillment for container shipping liners, Transp. Res. Part B Methodol., 120 (2019), 15–32. https://doi.org/10.1016/j.trb.2018.11.011 doi: 10.1016/j.trb.2018.11.011
[27]	L. Zhen, Y. Wu, S. Wang, G. Laporte, Green technology adoption for fleet deployment in a shipping network, Transp. Res. Part B Methodol., 139 (2020), 388–410. https://doi.org/10.1016/j.trb.2020.06.004 doi: 10.1016/j.trb.2020.06.004
[28]	C. F. Gao, Z. H. Hu, Speed optimization for container ship fleet deployment considering fuel consumption, Sustainability, 13 (2021), 5242. https://doi.org/10.3390/su13095242 doi: 10.3390/su13095242
[29]	J. Pasha, M. A. Dulebenets, A. M. Fathollahi-Fard, G. Tian, Y. Y. Lau, P. Singh, et al., An integrated optimization method for tactical-level planning in liner shipping with heterogeneous ship fleet and environmental considerations, Adv. Eng. Inf., 48 (2021), 101299. https://doi.org/10.1016/j.aei.2021.101299 doi: 10.1016/j.aei.2021.101299
[30]	X. Lai, L. Wu, K. Wang, F. Wang, Robust ship fleet deployment with shipping revenue management, Transp. Res. Part B Methodol., 161 (2022), 169–196. https://doi.org/10.1016/j.trb.2022.05.005 doi: 10.1016/j.trb.2022.05.005
[31]	K. Tierney, B. Áskelsdóttir, R. M. Jensen, D. Pisinger, Solving the liner shipping fleet repositioning problem with cargo flows, Transp. Sci., 49 (2015), 652–674. https://doi.org/10.1287/trsc.2013.0515 doi: 10.1287/trsc.2013.0515
[32]	K. Tierney, A. Coles, A. Coles, C. Kroer, A. Britt, R. Jensen, Automated planning for liner shipping fleet repositioning, in Proceedings of the International Conference on Automated Planning and Scheduling, 22 (2012), 279–287. https://doi.org/10.1609/icaps.v22i1.13500 doi: 10.1609/icaps.v22i1.13500
[33]	Y. F. Huang, J. K. Hu, B. Yang, Liner services network design and fleet deployment with empty container repositioning, Comput. Ind. Eng., 89 (2015), 116–124. https://doi.org/10.1016/j.cie.2015.01.021 doi: 10.1016/j.cie.2015.01.021
[34]	D. Müller, K. Tierney, Decision support and data visualization for liner shipping fleet repositioning, Inf. Technol. Manage., 18 (2017), 203–221. https://doi.org/10.1007/s10799-016-0259-3 doi: 10.1007/s10799-016-0259-3
[35]	D. Wetzel, K. Tierney, Integrating fleet deployment into liner shipping vessel repositioning, Transp. Res. Part E Logist. Transp. Rev., 143 (2020), 102101. https://doi.org/10.1016/j.tre.2020.102101 doi: 10.1016/j.tre.2020.102101
[36]	S. Kuhlemann, J. Ksciuk, K. Tierney, A. Koberstein, The stochastic liner shipping fleet repositioning problem with uncertain container demands and travel times, EURO J. Transp. Logist., 10 (2021), 100052. https://doi.org/10.1016/j.ejtl.2021.100052 doi: 10.1016/j.ejtl.2021.100052
[37]	X. Bai, L. Cheng, Ç. Iris, Data-driven financial and operational risk management: empirical evidence from the global tramp shipping industry, Transp. Res. Part E Logist. Transp. Rev., 158 (2022), 102617. https://doi.org/10.1016/j.tre.2022.102617 doi: 10.1016/j.tre.2022.102617
[38]	P. Balcombe, J. Brierley, C. Lewis, L. Skatvedt, J. Speirs, A. Hawkes, et al., How to decarbonise international shipping: options for fuels, technologies and policies, Energy Convers. Manag., 182 (2019), 72–88. https://doi.org/10.1016/j.enconman.2018.12.080 doi: 10.1016/j.enconman.2018.12.080
[39]	L. Zhang, L. Guan, D. Z. Long, H. Shen, H. Tang, Who is better off by selling extended warranties in the supply chain: the manufacturer, the retailer, or both?, Ann. Oper. Res., (2020), 1–27. https://doi.org/10.1007/s10479-020-03728-z doi: 10.1007/s10479-020-03728-z
[40]	L. Zhen, S. Wang, G. Laporte, Y. Hu, Integrated planning of ship deployment, service schedule and container routing, Comput. Oper. Res., 104 (2019), 304–318. https://doi.org/10.1016/j.cor.2018.12.022 doi: 10.1016/j.cor.2018.12.022
[41]	X. Xin, X. Wang, X. Tian, Z. Chen, K. Chen, Green scheduling model of shuttle tanker fleet considering carbon tax and variable speed factor, J. Clean. Prod., 234 (2019), 1134–1143. https://doi.org/10.1016/j.jclepro.2019.06.275 doi: 10.1016/j.jclepro.2019.06.275
[42]	B. D. Brouer, J. F. Alvarez, C. E. M. Plum, D. Pisinger, M. M. Sigurd, A base integer programming model and benchmark suite for liner-shipping network design, Transp. Sci., 48 (2013), 281–312. https://doi.org/10.1287/trsc.2013.0471 doi: 10.1287/trsc.2013.0471
[43]	Ship & Bunker (S & B), World bunker prices, 2022. Available from: https://shipandbunker.com/prices/av/global/av-g20-global-20-ports-average.
[44]	Y. Wu, Y. Huang, H. Wang, L. Zhen, Joint planning of fleet deployment, ship refueling, and speed optimization for dual-fuel ships considering methane slip, J. Mar. Sci. Eng., 10 (2022), 1690. https://doi.org/10.3390/jmse10111690 doi: 10.3390/jmse10111690
[45]	R. Greiner, Ship operating costs: current and future trends, 2017. Available from: http://greece.moorestephens.com/MediaLibsAndFiles/media/greeceweb.moorestephens.com/Documents/1-Richard-Greiner.pdf.
[46]	Y. Wu, Y. Huang, H. Wang, L. Zhen, Nonlinear programming for fleet deployment, voyage planning and speed optimization in sustainable liner shipping, Electron. Res. Arch., 31 (2023), 147–168. https://doi.org/10.3934/era.2023008 doi: 10.3934/era.2023008
[47]	S. Wang, Q. Meng, Robust schedule design for liner shipping services, Transp. Res. Part E Logist. Transp. Rev., 48 (2012), 1093–1106. https://doi.org/10.1016/j.tre.2012.04.007 doi: 10.1016/j.tre.2012.04.007
[48]	N. Rasmussen, Recorded webinar: tanker shipping market overview & outlook Q1 2023: the stars align to create the strongest market in 15 years, 2023. Available from: https://www.bimco.org/news/market_analysis/2023/20230228-smoo-tanker.
[49]	Moore Greece, Moore maritime index 2022, 2022. Available from: https://www.moore-greece.gr/el-gr/insights.
[50]	D. Sheppard, C. Cook, Russia assembles 'shadow fleet' of tankers to help blunt oil sanctions, 2022. Available from: https://www.ft.com/content/cdef936b-852e-43d8-ae55-33bcbbb82eb6.
[51]	Ç. Iris, J. S. L. Lam, Optimal energy management and operations planning in seaports with smart grid while harnessing renewable energy under uncertainty, Omega, 103 (2021), 102445. https://doi.org/10.1016/j.omega.2021.102445 doi: 10.1016/j.omega.2021.102445
[52]	W. Wang, Y. Wu, Is uncertainty always bad for the performance of transportation systems?, Commun. Transp. Res, 1 (2021), 100021. https://doi.org/10.1016/j.commtr.2021.100021 doi: 10.1016/j.commtr.2021.100021
[53]	B. Liu, Z. C. Li, Y. Wang, A two-stage stochastic programming model for seaport berth and channel planning with uncertainties in ship arrival and handling times, Transp. Res. Part E Logist. Transp. Rev., 167 (2022), 102919. https://doi.org/10.1016/j.tre.2022.102919 doi: 10.1016/j.tre.2022.102919
[54]	M. Lashgari, A. A. Akbari, S. Nasersarraf, A new model for simultaneously optimizing ship route, sailing speed, and fuel consumption in a shipping problem under different price scenarios, Appl. Ocean Res., 113 (2021), 102725. https://doi.org/10.1016/j.apor.2021.102725 doi: 10.1016/j.apor.2021.102725
[55]	J. Zhang, D. Z. Long, R. Wang, C. Xie, Impact of penalty cost on customers' booking decisions, Prod. Oper. Manage., 30 (2021), 1603–1614. https://doi.org/10.1111/poms.13297 doi: 10.1111/poms.13297
[56]	Y. Wang, J. Sarkis, Emerging digitalisation technologies in freight transport and logistics: current trends and future directions, Transp. Res. Part E Logist. Transp. Rev., 148 (2021), 102291. https://doi.org/10.1016/j.tre.2021.102291 doi: 10.1016/j.tre.2021.102291
[57]	H. Yalcin, T. U. Daim, Logistics, supply chain management and technology research: an analysis on the axis of technology mining, Transp. Res. Part E Logist. Transp. Rev., 168 (2022), 102943. https://doi.org/10.1016/j.tre.2022.102943 doi: 10.1016/j.tre.2022.102943
[58]	Y. Li, S. E. Li, X. Jia, S. Zeng, Y. Wang, FPGA accelerated model predictive control for autonomous driving, J. Intell. Connected Veh., 5 (2022), 63–71. https://doi.org/10.1108/JICV-03-2021-0002 doi: 10.1108/JICV-03-2021-0002
[59]	N. Lyu, Y. Wang, C. Wu, L. Peng, A. F. Thomas, Using naturalistic driving data to identify driving style based on longitudinal driving operation conditions, J. Intell. Connected Veh., 5 (2022), 17–35. https://doi.org/10.1108/JICV-07-2021-0008 doi: 10.1108/JICV-07-2021-0008
[60]	L. A. H. Hassan, M. Hewitt, H. S. Mahmassani, Daily load planning under different autonomous truck deployment scenarios, Transp. Res. Part E Logist. Transp. Rev., 166 (2022), 102885. https://doi.org/10.1016/j.tre.2022.102885 doi: 10.1016/j.tre.2022.102885
[61]	Y. Yang, B. Jia, X. Y. Yan, J. Li, Z. Yang, Z. Gao, Identifying intercity freight trip ends of heavy trucks from GPS data, Transp. Res. Part E Logist. Transp. Rev., 157 (2022), 102590. https://doi.org/10.1016/j.tre.2021.102590 doi: 10.1016/j.tre.2021.102590
[62]	M. G. Demissie, L. Kattan, Estimation of truck origin-destination flows using GPS data, Transp. Res. Part E Logist. Transp. Rev., 159 (2022), 102621. https://doi.org/10.1016/j.tre.2022.102621 doi: 10.1016/j.tre.2022.102621
[63]	R. Yan, S. Wang, Integrating prediction with optimization: models and applications in transportation management, Multimodal Transp., 1 (2022), 100018. https://doi.org/10.1016/j.multra.2022.100018 doi: 10.1016/j.multra.2022.100018
[64]	S. Wang, X. Tian, R. Yan, Y. Liu, A deficiency of prescriptive analytics—no perfect predicted value or predicted distribution exists, Electron. Res. Arch., 30 (2022), 3586–3594. https://doi.org/10.3934/era.2022183 doi: 10.3934/era.2022183
[65]	B. Niu, Z. Mu, B. Cao, J. Gao, Should multinational firms implement blockchain to provide quality verification?, Transp. Res. Part E Logist. Transp. Rev., 145 (2021), 102121. https://doi.org/10.1016/j.tre.2020.102121 doi: 10.1016/j.tre.2020.102121
[66]	S. Liu, G. Hua, Y. Kang, T. E. Cheng, Y. Xu, What value does blockchain bring to the imported fresh food supply chain?, Transp. Res. Part E Logist. Transp. Rev., 165 (2022), 102859. https://doi.org/10.1016/j.tre.2022.102859 doi: 10.1016/j.tre.2022.102859
[67]	T. M. Choi, T. Siqin, Blockchain in logistics and production from Blockchain 1.0 to Blockchain 5.0: an intra-inter-organizational framework, Transp. Res. Part E Logist. Transp. Rev., 160 (2022), 102653. https://doi.org/10.1016/j.tre.2022.102653 doi: 10.1016/j.tre.2022.102653
[68]	W. Yi, S. Wu, L. Zhen, G. Chawynski, Bi-level programming subsidy design for promoting sustainable prefabricated product logistics, Cleaner Logist. Supply Chain, 2021 (2021), 100005. https://doi.org/10.1016/j.clscn.2021.100005 doi: 10.1016/j.clscn.2021.100005
[69]	W. Yi, L. Zhen, Y. Jin, Stackelberg game analysis of government subsidy on sustainable off-site construction and low-carbon logistics, Cleaner Logist. Supply Chain, 2021 (2021), 100013. https://doi.org/10.1016/j.clscn.2021.100013 doi: 10.1016/j.clscn.2021.100013
[70]	S. Wang, R. Yan, A global method from predictive to prescriptive analytics considering prediction error for "predict, then optimize" with an example of low-carbon logistics, Cleaner Logist. Supply Chain, 4 (2022), 100062. https://doi.org/10.1016/j.clscn.2022.100062 doi: 10.1016/j.clscn.2022.100062
[71]	S. Wang, R. Yan, "Predict, then optimize" with quantile regression: a global method from predictive to prescriptive analytics and applications to multimodal transportation, Multimodal Transp., 1 (2022), 100035. https://doi.org/10.1016/j.multra.2022.100035 doi: 10.1016/j.multra.2022.100035
[72]	D. Huang, S. Wang, A two-stage stochastic programming model of coordinated electric bus charging scheduling for a hybrid charging scheme, Multimodal Transp., 1 (2022), 100006. https://doi.org/10.1016/j.multra.2022.100006 doi: 10.1016/j.multra.2022.100006

Reader Comments

Your name:*

Email:*
© 2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)