1 Introduction

We are living in a historic, very challenging era with the COVID-19 pandemic and the escalating negative impacts of climate change, with the latter fueling an increasing number of various disasters from hurricanes to heat waves and wildfires [1]. In addition to such “sudden-onset” disasters, climate change is also exacerbating droughts and rainfall in different areas and affecting agriculture and the food supply [2]. Furthermore, the increasing emissions due to the use of fossil fuels are reducing the quality of life and the health of people [3]. In our highly connected planet, viruses can spread easily through travel [4] and smoke, due to wildfires, can result in unhealthy air thousands of miles from the origin locations (cf. [5]).

Labor has been significantly and adversely impacted in this new era of COVID-19, coupled with climate change, with many firms and organizations realizing that they have to invest in labor for health and safety reasons [6]. Not everyone can work from home in the pandemic and that includes essential workers from various manufacturers, farmers and food processors, freight service providers, technicians, first responders, as well as healthcare workers, among others. Furthermore, working conditions, plus higher wages, have been critical in attracting and keeping workers throughout the pandemic; even now, as more and more people are vaccinated and as economies begin to gradually open up [7]. It is important to realize that climate change and COVID-19 are actually linked, in that it has been found that smoke from wildfires in the Pacific Northwest may have led to thousand more cases of COVID-19 and more deaths among those who tested positive for the coronavirus [8].

Early on in the pandemic, the large number of cases of COVID-19, as well as deaths of workers, including those in the food sector, from workers at meat processing plants in the USA to migrant farm labor in many countries captured the attention of many, including consumers who dealt with empty store shelves, a reduced portfolio of products, as well as higher prices of grocery items (see [9] and the references therein). Companies realized that they needed to invest in enhanced coronavirus mitigation procedures in their workplaces, from social distancing and the use of PPEs, to the implementation of various barriers, enhanced ventilation, etc. [10]. Such investments were seen in various grocery stores, manufacturing plants, restaurants, various retail outlets, and other work venues (cf. [11]).

The importance of labor to the production, transportation, storage, and distribution of products in supply chain networks has resonated in the COVID-19 pandemic. Firms are noting that keeping workers healthy and safe affects profits positively as well as a company’s reputation, since the negative fallout from widely publicized illnesses at food processing plants, factories, distribution centers, and warehouses can have lasting effects. Similarly, with global warming and climate change, workers, as well as firms, are needing to adapt since higher temperatures are affecting the productivity of labor (see [12]). De Lima et al. [13] note that labor supply and labor productivity are sensitive to increasing heat stress expected under climate change, especially in economic sectors that depend heavily on outdoor work; notably, agriculture [14, 15]. Being subjected to heat reduces the capability for physical activity in a wide range of work environments, and heat kills farm workers every year [16].

2 Literature Review and Contributions

In this paper, we construct a multiperiod supply chain network optimization model in which a firm seeks to determine its optimal enhancements in labor productivity associated with its supply chain network activities, at a cost, along with its optimal product path flows to the demand markets, under profit maximization, and applying a net present value (NPV) approach. We allow for the operational costs to depend not only on the product flow but also on the productivity enhancements, since there may be maintenance costs encumbered over the time horizon. We also allow the demand price functions to depend not only on the product demand at the demand markets but also on the productivity enhancements. Hence, our optimization model captures that consumers may be more (or less) sensitive to investments by a company in labor productivity. This feature has some analogues in the sustainability literature where consumers may be willing to pay higher prices for products that encumber lower emissions, that is, are more “green” (see, e.g., [17] and the references therein). Examples of enhancements associated with labor productivity can include better temperature-controlled work environments, regular access to water and restroom breaks, provision of PPEs and other type of sanitary supplies, better lighting, ventilation, and work area design, among others. Here we are specifically interested in labor productivity enhancements of the “hardware” variety, rather than education and teaching of new skills, since the focus is on the pandemic healthcare disaster as well as the climate change one (see [18]) and impacts on labor.

This work is the first to construct a supply chain network optimization model with labor and multiple time periods, with a focus on enhancing labor productivity, which is a topic of great current relevance in the pandemic and climate change, as well. There is only limited published work on the integration of labor as a primary resource into supply chain networks from a modeling perspective, and with an eye towards garnering insights that have emerged in the pandemic and are being further amplified because of climate change. Nagurney [19] constructed a supply chain network optimization model for perishable products, with a focus on food, in which there were bounds on the availability of labor on each link. This allowed for the quantification of the impacts of disruptions on labor availability in different supply chain network activities on product flows, consumer prices, and the profit of the firm. Nagurney [20] proposed another supply chain network optimization model with labor and three distinct sets of constraints, with the first consisting of upper bounds on labor availability on links; the second set consisting of an upper bound on labor availability on each tier of the supply chain network (production, transportation, storage, and distribution); and the third set of constraints consisting of a single labor bound in the supply chain network. Nagurney [21] then extended the latter paper to competing multiple firms using the Nash equilibrium (cf. [22, 23]) for the model under the first set of constraints and generalized Nash equilibrium as a concept for the second two sets of constraints, where firms competed for a limited amount of labor resources. All the above-noted publications, however, were single-period models, whereas the model in this paper is a multiperiod one. Nagurney [24] recently proposed a single-period supply chain game theory model with labor but with fixed amounts of labor on the links. In that model, production output increased as a function of wages and there were no investments in enhanced labor productivity. Nagurney [25], also in a single-period model, proposed investments in labor productivity, but the investments did not affect the demand price functions and the financial outlay had a specific form, unlike the general construct in this paper with respect to total cost associated with productivity enhancements, focusing on health and safety of workers. Addressing the health and safety of laborers is now of paramount importance with the double major stressors of the COVID-19 pandemic and climate change. A recent article [26] reveals how outdoor laborers involved in picking produce, who already have hard jobs, are suffering and even dying because of the extreme heat conditions, including in the Pacific Northwest of the USA.

Furthermore, the new model in this paper, in addition to being a multiperiod one, differs from the above noted published works in the following ways:

  1. 1.

    A supernetwork (cf. [27]) is used to abstract the multiperiod supply chain network optimization problem in order to precisely define the paths for the product path flows and to also capture the managerial activity associated with the management of the supply chains.

  2. 2.

    We consider labor availability to be wage-dependent in that more workers will be willing to work at higher wages, which is justifiable from economics and is being used as a strategy in the pandemic now with labor shortages in many sectors. The wages, hence, are not fixed, as was the case in the above-noted published papers.

  3. 3.

    We allow for investments on link labor productivity enhancements, which are associated with the health and safety of the laborers in their activities of production, transportation, storage, and distribution.

  4. 4.

    The link operational cost functions are dependent not only on the product flows but also on the link labor productivity enhancements made in the first time period, as are the demand price functions. The former capture the possible maintenance needed in subsequent time periods, whereas the latter reflects that consumers (to greater or lesser degrees) are willing to pay more for a product if workers have enhanced protection, which enables greater productivity.

It is important to emphasize that the economics literature is rich with respect to labor issues (cf. [28,29,30,31]) but not in terms of supply chains and the integration of labor in supply chains, which we focus on here.

This paper also adds to the disaster literature. For collections of interesting papers on many disaster themes as well as phases of disaster management, see the edited volumes by [32,33,34].

The paper is organized as follows: In Sect. 3, the multiperiod supply chain network model with link labor productivity enhancements is constructed and the supernetwork topology identified. We provide the variational inequality formulation of the problem, which is used for algorithmic purposes. Section 4 presents the algorithm and a series of numerical examples for which the complete input and output results are reported. Section 4 includes additional results in the form of sensitivity analysis exercises in order to ascertain quantitatively the impact on the firm’s NPV of changes in the discount rate, modifications to the managerial link productivity factors (which are not subject to productivity enhancements), and changes in the demand price functions to reflect greater sensitivity in the demand price functions to the labor productivity enhancements. Section 5 presents a summary of the results in this paper and the conclusions.

3 The Multiperiod Supply Chain Network Optimization Model with Investments in Labor Productivity Enhancements

We consider a single firm, which seeks to maximize the NPV of its investments in labor productivity enhancements over a finite period planning horizon. We consider a discretized time horizon T where \(t=1,\ldots ,T\). In each time period t, the firm has \(n_M\) production sites available, \(n_{DC}\) distribution centers at which storage can take place and serves \(n_R\) demand markets. The multiperiod supply chain network topology is depicted in Fig. 1, as a supernetwork, with the firm denoted as the super node 0. We make use of a supernetwork (cf. [27]) since it makes the definitions of paths precise and allows for the inclusion of managerial links. The firm in period t is then denoted by node (0, t) with lower tiered nodes corresponding, in time period t, respectively, to: production (\(M_{1,t},\ldots , M_{n_M,t}\)), distribution (\(DC_{1,t},\ldots , DC_{n_{DC},t}\)), and the demand markets (\({1,t},\ldots , {n_R,t}\)).

Fig. 1
figure 1

The Supernetwork Structure of the Multiperiod Supply Chain Network Model with Labor Productivity Enhancements

Let \(L_t\) denote the set of directed links of the firm in its supply chain network, at time period t, where \(t=1,\ldots ,T\), with the corresponding link from node 0 to node (0, t) also included. Let \(\hat{L}_t\) then denote the set of directed links of the firm in its supply chain network, at time period t, in the supernetwork in Fig. 1 not including the top-most links joining firm node 0 to node (0, t). The links emanating from the node 0 are grouped into the set \(L^0\). We denote a typical link by a. The full supply chain supernetwork in Fig. 1 is represented by the graph \(G=\left[ N,L\right]\), where N corresponds to the set of nodes and L to the set of links. \(\hat{L}\), in turn, consists of all links \(\hat{L}_t\); \(t=1,\ldots ,T\). The links in \(\hat{L}\) represent the supply chain network operational links, whereas the links in \(L^0\) are the managerial links. The managerial links represent the management of the supply chain networks that is needed in the different time periods. Investments can take place in any / all the operational links but do not take place on the managerial links. The managerial links, however, encumber wages to be paid, as do the operational links for labor provision. We assume that investment in the supply chain links takes place in the first time period with those investments then affecting operational link costs and demand price functions in the first and subsequent time periods.

Let \(P_{k,t}\) denote the set of paths joining the firm node 0 to demand market node (kt) representing demand market k in period t, where \(t=1,\ldots ,T\); \(k=1,\ldots ,n_R\). The set of paths P, with the number of elements equal to \(n_P\), consists of all the paths \(P_{k,t}\) for all k and t. Observe that the paths are acyclic and each path consists of a sequence of links representing the supply chain network activities of production, transportation, and distribution to a demand market, within a time period, as well as storage, if the link is over two time periods. Each path also includes a managerial link.

We assume, in the model, that the production sites, storage sites, transportation and distribution links, as well as the demand markets, are fixed from time period to time period. Hence, the number of production sites, the number distribution centers, and the number of demand markets as in time period 1 remain throughout the time horizon as in time period 1. Furthermore, we assume that investments are made only in the first time period. Similar assumptions in terms of supply chain investments, but not concerning labor, in multiperiod supply chain network models have also been made by [17, 35,36,37].

The additional notation for the model is given in Table 1. All vectors are column vectors.

Table 1 Notation for the Multiperiod Supply Chain Network Model with Labor Productivity Enhancements

The firm seeks to maximize the NPV of its investments in enhancing labor productivity over the planning horizon with the objective function denoted by U; that is:

$$\begin{aligned} \text{ Maximize }\quad U=& \ \sum _{t=1}^T {{1}\over {(1+r)^t}}\{\sum _{k=1}^{n_R}\rho _{k,t}(d_t,v^1)d_{(k,t)}-\sum _{a\in \hat{L}_t}\hat{c}_a(f_a,v_a^1)-\sum _{a\in \hat{L}_t}w_al_a\}\\ & -\sum _{a\in L^0}w_al_a-TC(v^1) \end{aligned}$$
(1)

subject to:

$$\begin{aligned} \sum _{p\in P_{k,t}}x_p=d_{k,t}, \quad k=1,\ldots ,n_R; t=1,\ldots ,T, \end{aligned}$$
(2)
$$\begin{aligned} f_a=\sum _{p\in P}x_p\delta _{ap}, \quad \forall a\in L, \end{aligned}$$
(3)
$$\begin{aligned} x_p\ge 0, \quad \forall p\in P, \end{aligned}$$
(4)
$$\begin{aligned} 0\le v_a^1\le v_a^{max}, \quad \forall a\in \hat{L}_1, \end{aligned}$$
(5)
$$\begin{aligned} f_a=(\alpha _a+v_a^1)l_a, \quad \forall a\in \hat{L}, \end{aligned}$$
(6)
$$\begin{aligned} f_a=\alpha _al_a, \quad \forall a\in L^0, \end{aligned}$$
(7)
$$\begin{aligned} l_a=\gamma _aw_a, \quad \forall a\in \hat{L}. \end{aligned}$$
(8)

Constraints (2), (3), and (4) are the conservation of product flow equations. Equation (2) guarantees that the demand for the product is met at each demand market in each time period. Equation (3) makes sure that the product flow on each link is equal to the sum of the product flows on paths that use that link. Equation (4) represents the nonnegativity assumption on the product path flows. The constraint (5) is the bounds on the labor productivity enhancements on the links, which takes place in time period 1 and then is sustained on the identical link in each subsequent time period. There are upper bounds since there may be physical limitations as to what productivity enhancements can be achieved. The expressions (6) and (7), in turn, reflect that the product output is equal to the link productivity function times the labor in hours available for that link activity (see also, [9, 19, 20]. Note that investments in link productivity enhancements are not allowed on managerial links. Note also that by (6) is also meant that the link productivity function in each time period on a given link is the same over all time periods after time period 1 for the same link (cf. Fig. 1) since the labor productivity link enhancement is determined in time period 1 for each link with the investment and then sustained throughout the planning horizon. In addition, here, according to (9), we assume that labor availability is linear with respect to the wage, which is consistent with the work of the Nobel Laureate Angus Deaton; cf. [38], where many references can be found as to such a frequently used structure by economists. Note that total financial outlay, in turn, that is, the total investments in link labor productivity enhancements, is captured in the function \(TC(v^1)\) in (1).

We now show that objective function (1) can be expressed solely in terms of product path flows and link labor productivity enhancements. Indeed, in view of (3), we can define link total operational cost \(\tilde{c}_a(x,v_a^1)\equiv \hat{c}_a(f_a,v_a^1)\), for all links \(a\in \hat{L}\). Also, in view of (2), we can define demand price function \(\tilde{\rho }_{k,t}(x, v^1)\equiv \rho _{k,t}(d, v^1)\), for all k and all t. Applying (6), (7), and (8) and, subsequently, (4), we conclude that

$$\begin{aligned}w_al_a=\left( {{\sum _{p\in P}x_p\delta _{ap}}\over {(\alpha _a+v_a^1)}}\right) ^2{{1}\over {\gamma _a}}, \quad \forall a \in \hat{L}, \end{aligned}$$
(9)
$$\begin{aligned}w_al_a=\left( {{\sum _{p\in P}x_p\delta _{ap}}\over {\alpha _a}}\right) ^2{{1}\over {\gamma _a}}, \quad \forall a \in L^0.\end{aligned}$$
(10)

Letting \(\tilde{U}(x,v^1)\equiv U\), we can now rewrite (1) as:

$$\begin{aligned} \tilde{U}(x,v^1)=\sum _{t=1}^T{{1}\over {(1+r)^t}}\{\sum _{k=1}^{n_R}\tilde{\rho }_{k,t}(x, v^1)\sum _{p\in P_{k,t}} x_p-\sum _{a\in \hat{L}_t} \tilde{c}_a(x,v_a^1)-\sum _{a\in \hat{L}_t}\left( {{\sum _{p\in P}x_p\delta _{ap}}\over {(\alpha _a+v_a^1)}}\right) ^2{{1}\over {\gamma _a}}\} \\ \sum _{a\in L^0}\left( {{\sum _{p\in P}x_p\delta _{ap}}\over {\alpha _a}}\right) ^2{{1}\over {\gamma _a}}-TC(v^1),\end{aligned}$$
(11)

which is subject to constraints (4) and (5).

The firm’s goal, hence, is to maximize (11) subject to constraints (4) and (5). We define the feasible set \(K\equiv \{(x,v^1)\vert \, (4) \,\text{ and }\, (5)\, \text{ hold }\}\). The feasible set K is convex. Also, we assume that the NPV function \(U(x,v^1)\) is continuously differentiable and concave. Then, it follows that an optimal solution to the above supply chain network optimization problem coincides with the solution of the variational inequality problem (cf. [39]): determine \((x^*,v^{1*})\in K\) such that

$$\begin{aligned} -\sum _{p\in P}{{\partial \tilde{U}(x^*,v^{1*})}\over {\partial x_p}}\times (x_p-x_p^*)-\sum _{a\in \hat{L}}{{\partial \tilde{U}(x^*,v^{1*})}\over {\partial v_a^1}}\times (v_a^1-v_a^{1*})\ge 0, \quad \forall (x,v^1)\in K. \end{aligned}$$
(12)

Variational inequality (12) can be put into standard form (cf. [39]), VI\((F,\mathcal{K})\), where we seek to determine a vector \(X^* \in \mathcal{K} \subset R^{\mathcal{N}}\), such that

$$\begin{aligned} \langle F(X^*),X-X^*\rangle \ge 0,\quad \forall X \in \mathcal{K}, \end{aligned}$$
(13)

where F is a given continuous function from \(\mathcal{K}\) to \(R^{\mathcal{N}}\), \(\mathcal{K}\) is a given closed, convex set, and \(\langle \cdot ,\cdot \rangle\) denotes the inner product in \(\mathcal{N}\)-dimensional Euclidean space.

In our supernetwork model \(X\equiv (x,v^1)\) and \(F(X)\equiv (F^1(X), F^2(X))\), where the p-th element of \(F^1(X)\) is: \(-{{\partial \tilde{U}(x,v^{1})}\over {\partial x_p}}\) and the \(a-\)th element of \(F^2(X)\) is: \(-{{\partial \tilde{U}(x,v^{1})}\over {\partial v_a^1}}\). In our model \(\mathcal{N}=n_P+n_{\hat{L}_1}\) and \(\mathcal{K}=K\).

We now proceed to describe the algorithm and to present the computed solutions to several numerical examples.

4 The Algorithm and Numerical Examples

The algorithm that we implemented for the solution of numerical examples is the modified projection method of [40]. It is guaranteed to converge to the solution of a variational inequality as in (13), provided that the function F(X) that enters the variational inequality is monotone and Lipschitz continuous and that a solution exists.

For completeness and convenience, we now recall the algorithmic steps.

4.1 The Modified Projection Method

Step 0: Initialization

Initialize with \(X^0 \in \mathcal{K}\). Set the iteration counter \(\tau :=1\) and let \(\zeta\) be a scalar such that \(0 < \zeta \le \frac{1}{\bar{L}}\), where \(\bar{L}\) is the Lipschitz constant.

Step 1: Computation

Compute \(\bar{X}^{\tau }\) by solving the variational inequality subproblem:

$$\begin{aligned} \langle \bar{X}^{\tau } + \zeta F({X}^{\tau -1}) - {X}^{\tau -1} , X - \bar{X}^{\tau } \rangle \ge 0, \quad \forall X \in \mathcal{K}. \end{aligned}$$
(14)

Step 2: Adaptation

Compute \(X^\tau\) by solving the variational inequality subproblem:

$$\begin{aligned} \langle X^\tau + \zeta F(\bar{X}^{\tau }) - {X}^{\tau -1}, X - X^\tau \rangle \ge 0, \quad \forall X \in \mathcal{K}. \end{aligned}$$
(15)

Step 3: Convergence Verification

If \(\vert X^\tau - {X}^{\tau -1} \vert \le \epsilon\), with \(\epsilon >0\), a pre-specified tolerance, then stop; otherwise, set \(\tau := \tau +1\) and go to Step 1.

4.2 Numerical Examples

The modified projection method was implemented in FORTRAN and a Linux system at the University of Massachusetts Amherst used for the computation of solutions to the numerical examples. The modified projection method was initialized as follows. The link productivity enhancements were set to 0.00. The initial demand at each market was set to 40 with the demand equally distributed among the paths terminating in each demand market. The convergence tolerance was \({10}^{-7}\). This means that the modified projection method was considered to have converged when the absolute value of the difference between each of the variables (the path flows and the link labor productivity enhancements) at two successive iterations differed by no more than this value. The parameter \(\zeta\) was set to .01 for each of the numerical examples, unless noted otherwise.

4.2.1 Examples 1, 2, and 3

In this set of examples, there are two time periods and the firm has two production facilities, a single distribution center for storage, and serves two demand markets. The supernetwork topology of this multiperiod supply chain network problem is depicted in Fig. 2. We consider a generic product for which we can draw insights based on the numerical computations and results. The purpose is to demonstrate that the model is computable and provides a breadth of information. We report the complete example input data as well as solutions for Examples 1, 2, and 3, followed by results of sensitivity analysis.

In Example 1, the link operational costs and the demand price functions are the same for each period and the demand price functions do not depend on the link productivity enhancements. The latter reflects a scenario in which consumers do not care about the work environment of the workers associated with the supply chain network for the product. Example 2 then has the identical data to the data in Example 1 but with the consumers now being responsive to the link labor productivity enhancements, which we consider to being such that make workers more comfortable, in terms of enhanced safety and health, so they can be more productive. Example 3, in turn, has the identical data to the data in Example 2 but now we consider the scenario that the consumers at the two demand markets in the second period are willing to pay a much higher price for the product than they did in the first time period. Hence, the demand price functions reflect this. The scenario in Example 3 could reflect a holiday purchasing scenario for a product or a product that is in greater demand in the pandemic.

Fig. 2
figure 2

Supernetwork Topology of the Multiperiod Supply Chain Network Numerical Examples

Example 1 – Same Data in Time Period 2 as in Time Period 1; Demand Price Functions Independent of Link Labor Productivity Enhancements

The total operational link cost functions are:

$$\begin{aligned} \hat{c}_a(\, f_a)&=2f_a^2+.05{v_a^1}^2,\quad \hat{c}_b(\,f_b)=2f_b^2+.1{v_b^1}^2,\quad \hat{c}_c(\,f_c)=.5f_c^2+.05{v_c^1}^2,\quad \hat{c}_d(\,f_d)=.5f_d^2+.05{v_d^1}^2,\\ \hat{c}_e(\,f_e)&=f_e^2+2f_e+.1{v_e^1}^2,\quad \hat{c}_f(\,f_f)=.5f_f^2+.1{v_f^1}^2, \quad \hat{c}_g(\,f_g)=.5f_g^2+.1{v_g^1}^2, \quad \hat{c}_h(\,f_h)=2f_h^2+.05{v_a^1}^2,\\ \hat{c}_i(\,f_i)&=2f_i^2+.1{v_b^1}^2,\quad \hat{c}_j(\,f_j)=.5f_j^2+.05{v_c^1}^2, \quad \hat{c}_k(\,f_k)=.5f_k^2+.05{v_d^1}^2, \quad \hat{c}_l(\,f_l)=2f_l^2+.1{v_e^1}^2,\\ \hat{c}_m(\,f_m)&=.5f_m^2+.1{v_f^1}^2,\quad \hat{c}_n(\,f_n)=0.00, \quad \hat{c}_o(\,f_o)=0.00. \end{aligned}$$

The demand price functions are:

$$\begin{aligned} \rho _{1,1}(d)=-5d_{1,1}+800,\quad \rho _{1,2}(d)=-5d_{1,2}+850. \end{aligned}$$
$$\begin{aligned} \rho _{2,1}(d)=-5d_{2,1}+800,\quad \rho _{2,2}(d)=-5d_{2,2}+850. \end{aligned}$$

The alpha link parameters are:

$$\begin{aligned} \alpha _a&=55,\quad \alpha _b=50,\quad \alpha _c=35,\quad \alpha _d=35,\quad \alpha _e=60,\quad \alpha _f=38,\quad \alpha _g=36,\\ \alpha _h&=55,\quad \alpha _i=50,\quad \alpha _j=35,\quad \alpha _k=35,\quad \alpha _l=60,\quad \alpha _m=38,\quad \alpha _n=10,\ \alpha _o=10, \end{aligned}$$

and the gamma link parameters are:

$$\begin{aligned} \gamma _a&=.1,\quad \gamma _b=.1,\quad \gamma _c=.09,\quad \gamma _d=.07,\quad \gamma _e=.08\quad \gamma _f=.06,\quad \gamma _g=.08,\\ \gamma _h&=.1,\quad \gamma _i=.1,\quad \gamma _j=.09,\quad \gamma _k=.07,\quad \gamma _l=.08\quad \gamma _m=.06,\quad \gamma _n=.1,\\ \gamma _o&=.1. \end{aligned}$$

The total investment cost function in productivity enhancements on the link is: \(TC=\sum _{a\in \hat{L}_1}{(v_a^1)}^2\).

The bounds on the link productivity enhancements are: \(v_a^{max}=\ldots =v_g^{max}=5\).

The discount rate \(r=.03\) (see [17] who also has used this value and provides a justification).

The paths are defined as: path \(p_1=(n,a,c,e)\), path \(p_2=(n,b,d,e)\), path \(p_3=(n,a,c,f)\), path \(p_4=(n,b,d,f)\), path \(p_5=(n,a,c,g,l)\), path \(p_6=(n,b,d,g,l)\), path \(p_7=(o,h,j,l)\), path \(p_8=(o,i,k,l)\), path \(p_9=(n,a,c,g,m)\), path \(p_{10}=(n,b,d,g,m)\), path \(p_{11}=(o,i,k,m)\), and path \(p_{12}=(o,i,k,m)\).

The modified projection method yields the following equilibrium product path flow pattern:

$$\begin{aligned} x_{p_1}^*&=22.21,\quad x_{p_2}^*=22.17,\quad x_{p_3}^*=26.54, \quad x_{p_4}^*=26.50,\\ x_{p_5}^*&=0.00,\quad x_{p_6}^*=0.00,\quad x_{p_7}^*=22.19, \quad x_{p_8}^*=22.15,\\ x_{p_9}^*&=0.00,\quad x_{p_{10}}^*=0.00,\quad x_{p_{11}}^*=26.52, \quad x_{p_{12}}^*=26.48. \end{aligned}$$

The equilibrium link flows and labor values are reported in Table 2, whereas the equilibrium link labor productivity enhancements and hourly wages are reported in Table 3.

Table 2 Equilibrium Link Flows and Labor Values for Examples 1, 2, and 3
Table 3 Equilibrium Link Productivity Enhancements and Hourly Wages for Examples 1, 2, and 3

The demand price at the first demand market in time period 1 is: 578.09 and at the second demand market the price is: 584.79; The demand price at the first demand market in time period 2 is: 578.31 and at the second demand market the price is: 585.02 with the corresponding equilibrium demands of: 44.38, 53.04, 44.34, and 53.00, respectively.

The value of the firm’s NPV (the objective function) is: 77,010.94.

There is no inventorying from time period 1 to time period 2 since \(f_g^*=0\). The equilibrium link flows in time period 2 are lower than their respective values in time period 1. The wages also decrease on each of the operational supply chain network links in time period 2 as compared to the respective values in time period 1. All the equilibrium link labor productivity enhancements are positive except for the inventorying link one which has \(v_g^{1*}=0\), which makes sense since the inventorying link is not used / needed in this example.

Example 2 – Same Data as in Example 1 but the Demand Price Functions Now Depend on the Link Productivity Enhancements

The data for Example 2 are identical to those in Example 1 but with the demand price functions at the demand markets in the two time periods expanded to include dependency on the link productivity enhancements. This is to model the possible scenario of consumers caring about the working conditions of the laborers in terms of health and safety, for example.

The demand price functions are now:

$$\begin{aligned} \rho _{1,1}(d, v^1)=-5d_{1,1}+10\sum _{a\in \hat{L}_1}v_a^1 +800,\quad \rho _{1,2}(d,v^1)=-5d_{1,2}+10\sum _{a\in \hat{L}_1}v_a^1+850. \end{aligned}$$
$$\begin{aligned} \rho _{2,1}(d,v^1)=-5d_{2,1}+10\sum _{a\in \hat{L}_1}v_a^1+800,\quad \rho _{2,2}(d,v^1)=-5d_{2,2}+10\sum _{a\in \hat{L}_1}v_a^1+850. \end{aligned}$$

The modified projection method yields the following equilibrium product path flow pattern:

$$\begin{aligned} x_{p_1}^*&=22.46,\quad x_{p_2}^*=22.42,\quad x_{p_3}^*=26.81, \quad x_{p_4}^*=26.78,\\ x_{p_5}^*&=0.00,\quad x_{p_6}^*=0.00,\quad x_{p_7}^*=22.44, \quad x_{p_8}^*=22.40,\\ x_{p_9}^*&=0.00,\quad x_{p_{10}}^*=0.00,\quad x_{p_{11}}^*=26.79, \quad x_{p_{12}}^*=26.75. \end{aligned}$$

The equilibrium link flows and labor values are reported in Table 2, whereas the equilibrium link labor productivity enhancements and hourly wages are reported in Table 3.

The demand price at the first demand market in time period 1 is now: 584.39 and at the second demand market the price is: 617.26. The demand price at the first demand market in time period 2 is: 584.61 and at the second demand market the price is: 591.10 with the corresponding equilibrium demands of: 44.81, 53.51, 44.77, and 53.46, respectively.

The equilibrium link product flows, the link labor values, as well as the link wages are all higher in Example 2 than in Example 1. The equilibrium demand market prices that consumers now pay for the product are greater in time period 2 than in time period 1, and the demands are higher. The link productivity enhancements are positive on links e and f in this example.

The firm enjoys an NPV of 80,030.05, which exceeds the NPV in Example 1. The results for this example strongly suggest that letting consumers know that a firm cares about laborers and is willing to invest in their health and safety in the firm’s supply chain network, and having consumers be responsive to such enhancements, can have a positive impact on the firm’s value of its objective function.

Example 3 – Same Data as in Example 2 but with the Demand Price Function Intercepts Increased in the Second Time Period

The data in Example 3 are identical to those in Example 2 but with the demand price function intercepts increased in time period 2; reflecting an increase in the prices the consumers are willing to pay at the demand markets in the second time period.

The demand price functions in the second time period are now:

$$\begin{aligned} \rho _{2,1}(d,v^1)=-5d_{2,1}+10\sum _{a\in \hat{L}_1}v_a^1+1600,\quad \rho _{2,2}(d,v^1)=-5d_{2,2}+10\sum _{a\in \hat{L}_1}v_a^1+1650. \end{aligned}$$

The modified projection method yields the following equilibrium product path flow pattern:

$$\begin{aligned} x_{p_1}^*&=18.57,\quad x_{p_2}^*=18.55,\quad x_{p_3}^*=22.58, \quad x_{p_4}^*=22.56,\\ x_{p_5}^*&=11.33,\quad x_{p_6}^*=11.30,\quad x_{p_7}^*=37.84, \quad x_{p_8}^*=37.79,\\ x_{p_9}^*&=14.69,\quad x_{p_{10}}^*=14.67,\quad x_{p_{11}}^*=41.21, \quad x_{p_{12}}^*=41.15. \end{aligned}$$

The equilibrium link flows and labor values are reported in Table 2, and the equilibrium link labor productivity enhancements and hourly wages are reported in Table 3. In this example, we set \(\zeta =.001\).

The demand price at the first demand market in time period 1 is: 627.34 and at the second demand market the price is: 676.08. The demand price at the first demand market in time period 2 is: 1121.63 and at the second demand market the price is: 1104.30 with the corresponding equilibrium demands of: 37.12, 45.13, 98.26, and 111.73, respectively.

The firm’s value of its NPV at the optimal solution in Example 3 is equal to: 197,428.52. In Example 3, the NPV is more than double the NPV in Example 2. Inventorying now takes place with \(f_g^*=51.99\). Wages now increase on all links. As in the previous two examples, those involved in the managerial links n and o earn the highest wages. Three links now have positive equilibrium productivity enhancement values at the equilibrium with the inventorying link g having the highest value with \(v^{1*}_g=0.66\).

4.2.2 Sensitivity Analysis

We now proceed to conduct several sensitivity analysis exercises. In the first exercise, we evaluate the effects of varying the discount rate on the NPV of the firm. In the second exercise, we explore the impact on the NPV of the firm of changing the \(\alpha _n\) and \(\alpha _o\) factors on the managerial links. In the final sensitivity analysis exercise, we modify the factor preceding the \(\sum _{a\in \hat{L}_1}v_a^1\) term in the demand price functions and examine the impact on the NPV of the firm. Example 3 is used as the baseline for these sensitivity analysis exercises with the data changes for each sensitivity analysis data point as noted on the x-axis of Figs. 3, 4, and 5.

Fig. 3
figure 3

Sensitivity Analysis for Different Discount Rates r and Effects on the Firm’s NPV

Fig. 4
figure 4

Sensitivity Analysis for Different Values of \(\alpha\) for Managerial Links o and n and Effects on the Firm’s NPV

Fig. 5
figure 5

Sensitivity Analysis for Different Values of Factor Preceding \(\sum _{a\in \hat{L}_1} v_a^1\) Expression in All the Demand Price Functions and Effects on the Firm’s NPV

The results of the first sensitivity analysis are reported in Fig. 3.

From Fig. 3, it can be seen that the NPV of the firm decreases linearly as the discount rate is raised. The fact that the firm’s NPV decreases as the discount rate increases is behavior that one would expect from corporate finance.

The results of the second sensitivity analysis are reported in Fig. 4. As has been emphasized earlier in this paper, one of the novel contributions herein is the use of a supernetwork for representing the multiperiod supply chain network optimization problem with link labor productivity enhancements through investments. The use of a supernetwork, with the addition of links to represent managerial activities, provides a richer framework for supply chain management since the wages of management with their labor hourly values and also productivity are now included. Also, recall that for the managerial links there are no link labor productivity enhancements and, therefore, the productivity factors on the managerial links n and o in our numerical examples consist solely of the factors \(\alpha _o\) and \(\alpha _n\). In Fig. 4, we see that as the values of these managerial link factors, which are the same on the managerial link in time period 1 and in time period 2, increase, then the NPV of the firm also increases, but at a decreasing rate over the range that this sensitivity analysis exercise was conducted.

The final sensitivity analysis results are reported in Fig. 5. In this exercise, all the demand price function coefficients preceding the \(\sum _{a\in \hat{L}_1}v_a^1\) term were the same and as delineated in the x axis in Fig. 5. Note that the result for the NPV for 10 is precisely that for Example 3. Here we see that the greater the value of this coefficient, the greater the NPV of the firm. This means that the firm should make sure that consumers respond to the investments of the firm in labor productivity enhancement and are made aware of the investments.

5 Summary and Conclusions

Supply chain networks serve as the critical infrastructure for the production, transportation, storage, and distribution of products globally. The COVID-19 pandemic, and the increase in climate-related disasters and their severity have demonstrated the need to maintain the functionality of supply chains. Essential to the functionality of supply chains is labor. Indeed, without labor, no product can be produced, transported, stored, and distributed. The productivity of labor, in turn, depends on workers being healthy and safe.

In this paper, we add to the literature on the inclusion of labor into supply chain networks, which was initiated in the COVID-19 pandemic and was inspired by news of workers falling ill and many perishing in the pandemic. There have been serious disruptions to supply chains for this and related reasons and, therefore, having rigorous mathematical frameworks that can capture intricacies of supply chain networks, along with labor and various features, including that of labor productivity, is valuable for firms as well as decision makers and policy makers. Although economists have historically contributed to the study of labor, along with wages that should be paid, and how to possibly address labor shortages, they have not considered supply chains in a general context. In this paper, we construct, for the first time, a multiperiod supply chain network optimization model, using the firm’s NPV, and a supernetwork formalism, to identify both optimal product flows and labor productivity enhancements. Another novelty of the model is that we include not only operational links but also managerial links since management needs to be involved in supporting the supply chain activities over the time periods. The model is formulated and solved as a variational inequality problem. The solution of the model yields optimal product path and link flows, optimal labor values associated with the various supply chain activities, wages that should be paid, as well as the optimal productivity enhancements.

The new model in this paper extends recently published models on labor and supply chain networks in several directions, including to multiple time periods, elastic wages, and the incorporation of investments in link labor productivity enhancements. The numerical examples provide insights that include that both the firm and the laborers can gain when a firm invests in labor productivity enhancement, which here we consider to enhance health and safety.

Future research can include the construction of multiperiod game theory models with the inclusion of labor that can also address the competition for labor, another vivid feature of the pandemic. It would also be very interesting to formulate service supply chains, as in the healthcare sector, in particular, to include labor specifically. It would be interesting to consider extensions of the model in this paper (as well as related ones) to allow for the number of demand markets to change from time period to time period, perhaps, with appropriate additional investments in different time periods. It would also be interesting to allow the number of production sites and distribution centers to change from time period to time period, with an increase in number having an associated investment (but note that the type of investment would be different from the investments in the model in this paper, which are for productivity of labor). Finally, allowing for investments in each time period in terms of labor productivity would also be a valuable extension to the model in this paper.