## Abstract

This paper develops a simulation-based mission mobility reliability (MMR) analysis framework to account for uncertainty in mobility prediction of off-road ground vehicles in mission planning. A concept of MMR is first proposed to quantify reliability of a mission path which passes through different types of soils. A single-loop Kriging surrogate modeling method is then employed to overcome the computational challenge in MMR assessment caused by expensive mobility simulations. Built upon the surrogate model-based MMR analysis, a dynamic updating scheme is proposed to update the MMR estimation using online mobility data, during the course of a specific mission and for a particular vehicle. The online dynamic updating of MMR allows us for effective and dynamic decision-making during the mission phase, thus proactively avoid rare events of immobility during the mission. A case study demonstrates the efficacy of the proposed MMR analysis and updating framework.

## 1 Introduction

Off-road vehicle mobility analysis and its modeling and simulation (M&S) are important steps in developing and operating off-road vehicles [1]. In battlefield especially, using M&S tools to predict the capability of a vehicle to go through a certain area can provide important guidelines in designing off-road ground vehicles and planning operations [2]. In order to make good use of off-road ground vehicle M&S, the U.S. Army Tank Automotive Research, Development, and Engineering Center (TARDEC) which is now called Ground Vehicle Systems Center (GVSC) developed the North Atlantic Treaty Organization (NATO) Reference Mobility Model (NRMM) in the 1960s [26]. Other than the developed NRMM and its upgrades, a next-generation NRMM (NG-NRMM) using multibody dynamic modeling and new simulation techniques is being developed to further improve the prediction capability and accuracy of the mobility models [2,7,8].

In addition to the NRMM and NG-NRMM, many other mobility analysis models have also been developed in the past decades to predict ground vehicle performance under different off-road conditions. For instance, using Bekker's equations, an analytical tool called Bekker's Derived Terramechanics Model (BDTM) has been developed to evaluate vehicle off-road mobility [9]. Recuero et al. [10] developed a mobility simulation model using nonlinear finite element tires on granular material. Serban et al. [11] proposed a co-simulation framework for vehicle–terrain interaction simulation. By developing a computer program called WinMaku, Korlath et al. performed a simulation of mobility performance for off-road vehicles. Also, Hetherington [12] suggested a model using mean maximum pressure (MMP) to specify off-road performance of vehicles. The problem is that all models above treat road conditions as deterministic values. However, in reality, there are many sources of uncertainty such as the uncertainty in vehicle parameters, uncertainty in soil property, and uncertainty in the soil map [13], which lead to uncertainty in the predicted off-road vehicle mobility. In order to account for uncertainty in mobility prediction, Lessem et al. [14,15] and Priddy [16] converted deterministic NRMM into a stochastic model so that the variability of vehicle mobility can be included in the prediction.

Similarly, for NG-NRMM, the stochastic mobility modeling is still one of the most important topics in mobility prediction [17]. Several approaches have been proposed to meet the need of performing stochastic mobility prediction in NG-NRMM. For example, using a hardware-in-the-loop real-time simulator, Bradbury et al. [7] proposed a Monte Carlo simulation method (MCS) to predict the vehicle performance under uncertainty; Gonzalez et al. used the Kriging surrogate modeling approach to model terrain elevation in predicting off-road mobility over large spatial regions [18]; they also used Kriging surrogate model to develop a stochastic mobility map [19]; a dynamic Kriging-based method was developed in Ref. [17] to improve the effectiveness of Kriging surrogate modeling in generating stochastic off-road mobility map; Hu et al. [13] proposed a mobility testing design method to effectively reduce the uncertainty in mobility prediction and select the best mobility model from competing mobility models; Detweiler et al. [20] developed interpolation approaches for terrain maps to get a high-fidelity 3D terrain surface.

Despite the above progress in stochastic vehicle mobility prediction, there are several challenging issues that need to be addressed to meet the demand of NG-NRMM. For instance, (1) current approaches only considered the uncertainty in mobility prediction at specific locations and have not investigated effects of mobility prediction uncertainty on mission planning; (2) even though Kriging surrogate model has been used to reduce the computational time in vehicle mobility prediction, the required computational effort is still high; and (3) the integration of the off-road ground vehicle M&S and vehicle online operation data have been rarely studied.

In order to solve the aforementioned challenging issues and fill the void in simulation-based mobility prediction under uncertainty, this paper develops a new simulation-based mission mobility reliability (MMR) analysis framework for off-road ground vehicles. In contrast to the current approaches which focus only on mobility prediction uncertainty at specific locations, the proposed framework not only investigates the effects of mobility dependence over space on mission planning, but also overcomes the computational challenges in mobility prediction under uncertainty without sacrificing the accuracy. The main contributions can be summarized as (1) definition of MMR for the first time to quantify the effects of mobility prediction uncertainty on mission planning; (2) employment of a single-loop Kriging surrogate modeling method [21] to efficiently evaluate MMR by refining the Kriging surrogate model in regions which are critical for the MMR analysis; and (3) development of a Bayesian updating scheme which allows for dynamical updating of MMR over the course of a mission based on vehicle operational data. The synthesis of the aforementioned three contributions enables the off-road ground vehicle to proactively avoid rare events of immobility during a mission, and thus guarantee the safety of the vehicle on the battlefield.

The remainder of this paper is organized as follows: Sec. 2 introduces the background of off-road vehicle mobility model and MMR. Section 3 describes the proposed simulation-based off-road ground vehicle MMR analysis framework. Section 4 uses a case study to demonstrate the proposed framework. Section 5 concludes the results, contributions, and future works.

## 2 Background

### 2.1 Off-Road Vehicle Mobility Model.

Off-road vehicle mobility model is an analytical or simulation model that is used to predict the vehicle mobility characteristics for given off-road conditions and vehicle design [2,7,8]. At TARDEC (now called GVSC) or cooperating with TARDEC, many mobility models have been developed in the past, such as the Bekker's Derived Terrain-mechanics Model [9], the Army NATO Reference Mobility Model (NRMM), the absorbed power criteria, and the Dynamic Analysis and Design System (DADS) [6]. Different mobility models may have different quantities of interest. For instance, in the MMP model, the MMP is related to the vehicle parameters as follows [12]:
$MMP=1.26Wmwbpwdw$
(1)
where W is the weight of the vehicle, mw is the total number of road-wheels, b is the track width, pw is the track plate length and dw is the road-wheel diameter.
In the GO/NO-GO condition model [22], the decision-making parameter CLLwheels is given by
$CLL−wheels=1.85W⋅2naba0.8da0.8δ0.4$
(2)
in which na is the number of axles, ba and da are the overall width and diameter of the tire inflated but unloaded, and δ is the tire deflection when loaded.
In Bekker's Terramechanics model [9], the key parameter displacement x is given by
$x=A1e(−bx+(bx2−1))ωt+A2e(−bx+(bx2−1))ωt$
(3)
where bx is the coefficient of damping, ωt is the natural time frequency of an aperiodic vibration, and A1 and A2 are coefficients related to soil stress and soil deformation.

Among the developed models, NRMM is one of the most commonly used models. Its development started in the 1950s [23]. There are also other upgrades such as NRMM II, in order to increase the performance of the mobility prediction model [16,24]. Since NRMMs were developed several decades ago, the performance of these models is limited by the analysis methods and computing capabilities. With the development of powerful computational tools such as finite element analysis (FEA), high-performance computation, and simulation techniques, and with the need of developing advanced vehicle systems such as autonomous vehicles, there is an urgent need to develop the NG-NRMM [25]. One of the important aspects of the NG-NRMM is to account for the uncertainty in the vehicle mobility prediction.

No matter what kind of off-road mobility model is used, the inputs to the models (NRMM or NG-NRMM) can always be classified into two groups, namely, vehicle-related parameters and terrain-related parameters. The vehicle-related parameters include vehicle dimension, weight, power, movement type (wheel or track), and so on. Terrain-related parameters include slope, cohesive strength, friction coefficient, and bulk density. In addition, terrain-related parameters are usually space dependent (i.e., vary with spatial coordinates). Let X ∈ ℝn be a vector of vehicle parameters, Y ∈ ℝm be a vector of terrain parameters which are functions of spatial coordinates d, and a vehicle mobility model can be generalized as
$V=G(X,Y(d))$
(4)
where V ∈ ℝ is the vehicle mobility variable (e.g., the maximum attainable speed in this paper).

### 2.2 Mission Mobility Reliability.

As mentioned earlier, accounting for the uncertainty in vehicle mobility prediction is one of the most important aspects of NG-NRMM. The inherent uncertainty in the terrain and soil properties and variability of the vehicle performances lead to the uncertainty of the vehicle mobility. As shown in Fig. 1, if a vehicle design and a soil region are selected, a soil map can be used in conjunction with the vehicle mobility model to generate a mobility map. Due to the uncertainty in both the vehicle systems and the soil map, there is uncertainty in the vehicle mobility for any given spatial location d.

Fig. 1
Fig. 1
Close modal
In order to quantify the effects of uncertainty sources on the vehicle mobility, at a specific spatial location d, the off-road vehicle mobility reliability $RM(d)$ is defined as follows:
$RM(d)=∫∫G(X,Y(d))≥efX(x)fY(d)(y)dxdy$
(5)
where e is the threshold of decision-making parameter, if $V=G(X,Y(d))≥e$, it means that the vehicle can achieve the minimum required maximum attainable speed at this location; otherwise, the vehicle will be immobile (i.e., failed); $fX(x)$ and $fY(d)(y)$ are, respectively, the probability density function (PDF) of vehicle parameters and terrain parameters at d.
The above mobility reliability, which is also referred as the “probability of mobility,” has been employed in the literature to generate the probability of mobility maps [13]. Even though the mobility reliability is able to effectively quantify the uncertainty of mobility at any given spatial location, it cannot be directly applied to reliability-based mission planning, due to the fact that the stochastic mobility map does not provide information about the correlation of mobility over space. The correlation of vehicle mobility stems from correlations of soil properties over space. For a particular soil type and a certain soil property, the soil property needs to be modeled as a random field over space instead of a random variable. In this paper, the soil properties are assumed to be Gaussian random fields for the sake of illustration. The developed framework can also be extended to account for the situation that the random fields are non-Gaussian. Similar to many other random field problems, one way to model/simulate the variability of soil property over space is to employ the Karhunen–Loève (KL) expansion method [26]. For a given terrain parameter or soil property over a spatial domain $d∈Ω$, Y(d) is modeled using the KL expansion method as follows:
$Y(d)=μY(d)+σY(d)∑j=1NKLλjξjηj(d)$
(6)
in which $μY(d)$ and $σY(d)$ are, respectively, the mean and standard deviation of $Y(d)$, λj and $ηj(d)$ are the eigenvalues and eigenvectors of $ρ≜{ρ(di,dj),∀i,j=1,2,⋯,Nd}$, Nd is the number of discretized spatial points on a mission path Ω (as shown in Fig. 2), ξj, j = 1, 2, · · ·, NKL are standard normal random variables, NKL is the number of the largest eigenvalues that makes 95% of the summation of all eigenvalues, and the correlation $ρ(di,dj)$ between two locations is given by
$ρ(di,dj)=exp[−β1(d1i−d1j)2−β2(d2i−d2j)2]$
(7)
where $di=[d1i,d2i]$ are the two spatial coordinates at location i, β1 and β2 are respectively the correlation length parameters of coordinates 1 and 2. Note that the developed framework is not limited to above squared exponential function. The other types of correlation functions such as constant, linear, and Ornstein–Uhlenbeck [27] can also be employed to describe the correlation between two spatial points. The best correlation function needs to be determined according to the geostatistical data of the soil.
Fig. 2
Fig. 2
Close modal
Furthermore, this paper defines the mission mobility reliability (MMR) to relate the variability and correlation of vehicle mobility over space with mission planning. For a given mission plan Ω as shown in Fig. 2, the MMR, R(Ω), is defined as
$R(Ω)=Pr{V=G(X,Y(d))≥e,∀d∈Ω}$
(8)
where “$∀$” stands for “for all” and $Pr{⋅}$ is a probability operator. It implies that only if the speed of a vehicle is always greater than or equal to the threshold e on the mission path, the mission is considered as successful, otherwise is considered as failed.
The corresponding mission mobility failure probability is given by
$pf(Ω)=Pr{V=G(X,Y(d))
(9)
in which “$∃$” means “there exists.”

The aforementioned MMR quantifies the probability that the vehicle completes the mission path without immobility. Once an event of immobility happens, it implies that the mission is failed. Evaluating this MMR plays a vital role in guaranteeing the safety of the off-road ground vehicle in the battlefield and avoiding events of immobility, especially for autonomous ground vehicles.

Due to the involvement of computer simulation models of vehicle–terrain interaction and soil maps, the evaluation of Eqs. (8) or (9) is not straightforward. There are two main challenges that need to be addressed to perform the MMR analysis.

• First, the mobility prediction model given in Eqs. (8) and (9) is usually high-fidelity M&S. The evaluation of MMR needs to execute the high-fidelity M&S numerous times, which is computationally prohibitive.

• Second, there are many sources of uncertainty presented in the mobility analysis. Some of the uncertainty sources, such as the variability in the soil properties, can be reduced when a particular vehicle starts to proceed on a pre-specified mission path. As the vehicle proceeds over the path, more and more field information about the vehicle mobility will be collected. The probability of accomplishing the mission (i.e., MMR) will then change over time and needs to be updated. How to leverage the field information in the MMR analysis is the second challenging issue that needs to be addressed.

In order to tackle the above-mentioned challenges, in the subsequent section (Sec. 3), a simulation-based off-road ground vehicle MMR analysis approach is proposed.

## 3 Simulation-Based Off-Road Ground Vehicle Mission Mobility Reliability Analysis

In this section, we first provide a brief overview of the proposed simulation-based MMR analysis framework. Following that, we explain the proposed framework in details.

### 3.1 Overview.

Figure 3 shows an overview of the proposed framework. As shown in this figure, there are two phases of off-road ground vehicle MMR analysis, namely (1) mission-planning phase and (2) during-mission phase. The two phases are connected through the mobility M&S as indicated in the figure. The main objectives of the two phases are explained as below:

1. Phase 1-mission-planning phase: In the mission-planning phase, the MMR of a given mission path is evaluated with the consideration of the natural variability of soil properties over the mission path. The main objective is to address the challenges in solving Eq. (8) to provide the probability of completing a mission in the early planning stage.

2. Phase 2-during-mission phase: In the mission-planning phase, the MMR analysis is purely based on the M&S and soil maps. Even though there is uncertainty in the soil properties and vehicle system, the soil properties and vehicle system parameters are deterministic values in reality for a particular vehicle and a specific mission path. As a result, a mission can only have two states: success or fail. In the during-mission phase, vehicle mobility measurement data collected through sensors on the vehicle are used to update the MMR estimate dynamically based on the MMR analysis approach developed in Phase 1. This dynamic updating of MMR allows us to predict the failure (i.e., immobility) with high confidence before it happens and thus proactively avoid the situation of immobility. As illustrated in Fig. 3, the MMR gets close to the ground truth (i.e., success or fail) as the vehicle proceeds on the mission path.

Fig. 3
Fig. 3
Close modal

Next, we provide more details about the two phases of the simulation-based MMR analysis.

### 3.2 Simulation-Based Mission Mobility Reliability Analysis in Mission-Planning Phase.

In the mission-planning phase, the simulation-based MMR analysis is very similar to the time-dependent reliability analysis problems [28,29]. In recent years, various time-dependent reliability analysis approaches, such as first-order reliability method-based methods [30,31], sampling-based approaches [32,33], and adaptive surrogate-modeling-based methods [3438], have been proposed to overcome the computational issue in reliability analysis. For the MMR analysis, the reliability is space-dependent instead of time-dependent. Moreover, the MMR needs to be updated in the during-mission phase using the Bayesian method. Considering the fact that the M&S is shared by both the mission-planning and during-mission phases, the surrogate-modeling-based method is employed. More specifically, the single-loop Kriging (SILK) surrogate modeling approach [21] is employed in this paper, since SILK can easily handle the presence of both random variables and stochastic processes in reliability analysis.

The employment of SILK in the mission-planning phase, however, is not very straightforward. Two issues need to be addressed in order to perform MMR analysis with SILK. The first issue is the simulation of soil properties over the mission path, and the second one is the generation of training points for surrogate modeling. In what follows, we discuss how to tackle these two issues in detail.

#### 3.2.1 Simulation of Soil Properties Over the Mission Path.

Let $Y(d)≜[Y1(d),Y2(d),⋯,Ym(d)]$, $∀d∈Ω$, be the space-dependent soil properties or terrain parameters, where m is the number of (slope or soil) property variables, for a given mission path Ω as indicated in Fig. 2, $Y(d)$ may consist of multiple slope types or soil types which are described as random fields. It implies that even for a single soil property variable, $Yi(d),∀d∈Ω$, it consists of multiple random fields with different statistical properties. Taking the terrain slope as an example, to generate random realizations of $Yi(d),∀d∈Ω,i=1,2,⋯m$ for the MMR analysis, as illustrated in Fig. 4, the coordinates $d(p)Slope$ corresponding to the pth slope ID are first identified from Ω as
$d(p)Slope={d:whereS(d)=p,∀d∈Ω},p=1,2,⋯,Nslope$
(10)
where $S(d)$ is the slope ID at location d obtained from the terrain map (see Sec. 4.1 for a terrain map example), Nslope is the number of slope IDs, and we have
${d(1)Slope∪⋯∪d(Nslope)Slope}=Ω$
(11)
Fig. 4
Fig. 4
Close modal
Similarly, the coordinates $d(h)Soil$ of the hth soil ID are identified as
$d(h)Soil={d:whereSL(d)=h,∀d∈Ω},h=1,2,…,Nsoil$
(12)
in which $SL(d)$ is the soil ID obtained from the soil map (see Sec. 4.1), Nsoil is the number of soil types, and ${d(1)Soil∪⋯$$∪d(Nsoil)Soil}=Ω$.
After partitioning the coordinates Ω on the mission path into $d(1)Slope,d(2)Slope,…,d(Nslope)Slope$ for slope IDs and $d(1)Soil,d(2)Soil,…,d(Nsoil)Soil$ for soil types, the slope or soil property corresponding to a particular slope ID or soil ID is described using KL expansion given in Eq. (6). Taking $Yi(d),∀d∈Ω$ as an example, according to the partition of the coordinates on the mission path, $Yi(d),∀d∈Ω$ is expressed as a union of multiple random fields as follows:
$Yi(d)={Yi1∪⋯∪YiNT},∀d∈Ω$
(13)
in which $Yij≜{Yij(1),Yij(2),…,Yij(Ns(j))}$ where Yij(k), k = 1, 2, …, Ns(j) represents $Yi(d)$ associated with the jth slope ID or soil ID at the kth location, Ns(j) is the number of spatial points associated with the jth slope or soil ID (i.e., the number of elements in $d(j)Slope$ or $d(j)Soil$), and NT is given by
$NT={Nslope,ifYi(d)isaslope−relatedvariableNsoil,ifYi(d)isasoiltype−relatedvariable$
(14)
$Yij(k),∀j=1,2,…,NT;k=1,2,…,Ns(j)$ are modeled with KL expansion as follows:
$Yij(k)=μij+σij∑q=1Nijλqξqηq(k)$
(15)
where μij and σij are, respectively, the mean and standard deviation of Yij(k), Nij is the number of expansion terms, and λq and ηq are the eigenvalues and eigenvectors of the following correlation matrix:
$ρij=[1ρij(d(1)S,d(2)S)⋯ρij(d(1)S,d(NS(j))S)ρij(d(2)S,d(1)S)1⋯ρij(d(2)S,d(NS(j))S)⋮⋮⋱⋮ρij(d(NS(j))S,d(1)S)ρij(d(NS(j))S,d(2)S)⋯1]NS(j)×NS(j)$
(16)
in which ρij(·, ·) is the correlation function (i.e., Eq. (7)) of Yij(k), $d(k)S=d(j)Slope(k),∀k=1,2,…,Ns(j)$ if Yij(k) is a slope-related variable; Otherwise, $d(k)S=d(j)Soil(k),∀k=1,2,…,Ns(j)$ (i.e., if Yij(k) is a soil-type-related variable).

Using Eqs. (10) through (16), random realizations can be generated for $Yij,∀j=1,2,…,NT$. As illustrated in Fig. 4, the random realizations of $Yij$ are then assembled together according to the coordinates on the mission path to obtain random realizations of $Yi(d),∀d∈Ω$. The random realizations of $Yi(d),∀d∈Ω$ are used as the inputs in the SILK method [21] for (1) the selection of new training points in refining the surrogate model of the mobility model and (2) MMR analysis. Next, we will discuss how to generate initial training points for the surrogate modeling of the mobility model $V=G(X,Y(d))$.

#### 3.2.2 Generation of Training Points for Surrogate Modeling.

Generation initial training points. For any specific coordinates $d$, Fig. 5 shows the relationships between the map coordinates, soil type (i.e., ID), slope type/ID, and space-dependent parameters $Y(d)$.

Fig. 5
Fig. 5
Close modal

As shown in the previous figure, the soil property parameters share the same soil ID for any given spatial point, $d$. The soil ID represented as integer values governs the statistical parameters of soil properties (e.g., bulk density, friction coefficient, and cohesive strength). It indicates that there are strong dependences among the soil-type-related parameters. In order to maintain the strong dependences in generating initial training points for $Y(d)$, in this paper, two sets of initial training points are first generated using Latin hypercube sampling method [39].

The two sets of training points are, respectively, training samples of soil and slope IDs (i.e., integer) and training samples in the uniform domain of [0, 1]. Let the integer training samples of soil/slope IDs be $α=[α(1),α(2),⋯,α(Nin)]$ where Nin is the number of initial training points and the uniform domain training points be $u=[u(1),u(2),⋯,u(Nin)]$, $u$ are then transformed into training points of $Y(d)$ as follows:
$yt(l)=fs(u(l),θ(α(l))),∀l=1,2,…,Nin$
(17)
where fs(·, ·) is a function mapping $u(l)$ to training point $yt(l)$ of $Y(d)$ and $θ(α(l))$ are the statistical parameters (e.g., mean and standard deviation) related to soil/slope ID α(l). Note that $yt(l)$ here represents multiple slope/soil property parameters that share a same soil/slope ID α(l).
For instance, if $θ(k)=[θL(k),θU(k)],k=1,2,⋯,Nsoil$ are the lower and upper bounds of the training interval for a given slope/soil ID k, we have $yt(l)$ as
$yt(l)=θL(α(l))+u(l)(θU(α(l))−θL(α(l))),∀l=1,2,⋯,Nin.$
(18)
After defining the initial training points of $X$ and $Y(d)$ as $xt$ and $yt$, we then obtain the corresponding mobility responses as $vt$. Using $[xt,yt]$ and $vt$, a surrogate model is constructed for the computationally expensive vehicle mobility model $V=G(X,Y(d))$ using the Kriging surrogate modeling method as below [21,40]:
$V=G^(X,Y(d))$
(19)
in which $V=G^(X,Y(d))$ stands for the Kriging surrogate model which is an approximation of the original high-fidelity mobility model.
According to the property of the Kriging surrogate model, we have
$V∼N(μV(X,Y(d)),σV2(X,Y(d)))$
(20)
where $μV(X,Y(d))$ and $σV(X,Y(d))$ are, respectively, the mean and standard deviation of the Kriging surrogate model prediction, and N(·, ·) stands for normal distribution.

Since $V=G^(X,Y(d))$ may not accurately represent the original high-fidelity mobility model, it could lead to large error in the MMR analysis. To improve the accuracy of reliability analysis, many adaptive surrogate modeling methods have been proposed. As mentioned earlier, the SILK method is employed in this paper for the MMR analysis because it is capable to handle time-dependent problem which is similar to the space-dependent problem. Next, we discuss how to adaptively add new training points to refine the surrogate model.

Generation of new training points through adaptive training. In SILK [21], the surrogate model is adaptively refined in important regions to guarantee the accuracy of the reliability analysis. Let the random samples of $X$ be $xMCS≜{xiMCS,i=1,2,…,nMCS}$ and the random realizations generated for $Y(d)$ using Eqs. (10) through (16) be $yMCS≜{yijMCS,i=1,2,…,nMCS;j=1,2,…,Nd}$ where nMCS is the number of realizations in MCS, Nd is the number of discretized spatial points on the mission path Ω, and $yijMCS$ indicates the ith realization of $Y(d)$ at jth spatial location; the new training point for $X$ and $Y(d)$ is identified in SILK as follows [21]:
$xnew,ynew=[xi*MCS,yi*(j*|i*)MCS]$
(21)
where $i*$ and $j*|i*$ are identified by
$i*=argmini=1,2,⋯,nMCS{Qmin(i)}$
(22)
$Qmin(i)={Qe,ifμV(xiMCS,yijMCS)
(23)
in which Qe is a value greater than 2 and $Q(xiMCS,yijMCS)$ is given by [21]
$Q(xiMCS,yijMCS)=|μV(xiMCS,yijMCS)−e|σV(xiMCS,yijMCS)$
(24)
and
$j*|i*=argminj=1,2,⋯,Nd{Q(xi*MCS,yi*jMCS)}$
(25)
After the new training point $xnew,ynew$ is identified, the corresponding mobility response is obtained. The new training point is then added to the training point data set and the surrogate model $V=G^(X,Y(d))$ is re-trained. The process continues until the accuracy criterion is satisfied. In SILK, the accuracy criterion is defined as below [21]
$εrmax<0.05$
(26)
where
$εrmax=maxNf2*∈[0,N2]{|Nf2−Nf2*|Nf1+Nf2*}$
(27)
in which $N2=∑i=1nMCSId{Qmin(i)<2},Nf1=$$∑i=1nMCSId{(μV(i,j), $Nf2=Id{(μV(i,j), where $μV(i,j)≜μV(xiMCS,yijMCS)$ and $Id{E}$ is an indicator function. $Id{E}=1$, if event E is true. Otherwise, $Id{E}=0$.
When the accuracy criterion given in Eq. (26) is satisfied, we terminate the adaptive training of $V=G^(X,Y(d))$. After that, the MMR of a given mission path Ω is estimated as follows:
$R(Ω)≈1nMCS∑i=1nMCSId{μV(i,j)≥e,∀j=1,2,…,Nd}$
(28)
where $μV(i,j)≜μV(xiMCS,yijMCS)$ is the mean prediction obtained by plugging $xMCS$ and $yMCS$ into the surrogate model $V=G^(X,Y(d))$.

Until now, we have discussed all the details for simulation-based MMR analysis in the mission-planning phase. By overcoming the computational challenge in this phase, we are able to efficiently perform the MMR analysis and ensure the MMR in mission planning. While the planned mission path usually has a high MMR, there is still a risk of having immobility (i.e., failure) events when the vehicle proceeds on the path. The dynamic updating of MMR in the during-mission phase (as discussed in Sec. 3.3) allows us to sense the immobility events before it happens and thus proactively avoid rare events of failures (i.e., immobility).

### 3.3 Dynamic Updating of Mission Mobility Reliability in the During-Mission Phase.

The goal of dynamic updating of MMR in the during-mission phase is to proactively forecast the failure based on vehicle measurement data given that a vehicle has already accomplished part of the mission path (as illustrated in Fig. 6).

Fig. 6
Fig. 6
Close modal
Denoting the spatial points of the finished part and unfinished part of the mission path as, respectively, Ωf and Ωuf, we have $Ωf∪Ωuf=Ω$. Let the measured vehicle mobility data of the finished part Ωf be $v1:c≜{v1,v2,⋯,vc},∀c=1,2,⋯,Nd$, where $vi=v(di)$ are the mobility data collected at spatial point $di$, and $dc$ is the current location: the task is to dynamically estimate the following MMR conditioned on $v1:c$
$R(Ω)|v1:c=Pr{V=G^(X|v1:c,Y(d)|v1:c)≥e,∀d∈Ωuf}$
(29)
where $V=G^(X,Y(d))$ is the refined surrogate model from Phase 1, $X|v1:c$ and $Y(d)|v1:c$ are, respectively, $X$ and $Y(d),d∈Ωuf$ conditioned on $v1:c$.

A critical step in evaluating Eq. (29) is to update $X$ and $Y(d),d∈Ωuf$ based on $v1:c$. The challenge of performing this updating is that there is no direct relationship between the random variables $X$ and $Y(d),d∈Ωuf$ and $v1:c$. Fortunately, the vehicle mobility M&S provides such a bridge between these variables and the measurement data. In this section, a sequential updating scheme is proposed to dynamically update MMR (as shown in Fig. 7). It consists of three main steps, namely, (1) Bayesian inference of uncertainty parameters, (2) uncertainty sources updating, and (3) MMR updating and remaining mobilable distance estimation.

Fig. 7
Fig. 7
Close modal

#### 3.3.1 Bayesian Inference of Uncertain Parameters in the Mobility Model.

The proposed framework starts with the online mobility measurement data $v1:c$ and sequentially estimates $x$ and $y(d)$ based on $v1:c$. The estimated $x$ and $y(d)$ are then used to update $X$ and $Y(d),d∈Ωuf$ using conditional Gaussian random field. For given $v1:c$, the posterior PDF of $X$ and $Y(d)$ at $dc$ is estimated using the recursive Bayesian method as below
$f(x,y(dc)|v1:c)=f(vc|x,y(dc))f(x,y(dc)|v1:c−1)∫∫f(vc|x,y(dc))f(x,y(dc)|v1:c−1)dxdy(dc)∝f(vc|x,y(dc))f(x,y(dc)|v1:c−1)$
(30)
where “∝” stands for “proportional to,” $f(vc|x,y(dc))$ is the likelihood function of observing $vc$ for given $x$ and $y(dc)$ at $dc$, and $f(x,y(dc)|v1:c−1)$ is given by
$f(x,y(dc)|v1:c−1)=∫f(y(dc)|y(dc−1))f(x,y(dc−1)|v1:c−1)dy(dc−1)$
(31)
in which $f(y(dc)|y(dc−1))$ is the PDF of $y(dc)$ conditioned on $y(dc−1)$.
As illustrated in Fig. 7, the mobility prediction model $V=G^(X,Y(d)),$ is employed to connect mobility measurement data $vc$ with $x$ and $y(dc)$. $f(vc|x,y(dc))$ is thus computed by
$f(vc|x,y(dc))=∏j=1Ncϕ(vcj−μV(x,y(dc))σε2+σV2(x,y(dc)))$
(32)
where ϕ(·) is the PDF of a standard normal variable, $vc=[vc1,vc2,…,vcNc]$, σɛ is the standard deviation of the observation noise, Nc is the number of observations at location $dc$, $μV(x,y(dc))$ and $σV(x,y(dc))$ are respectively the mean and standard deviation of the vehicle mobility prediction obtained from the surrogate model $V=G^(X,Y(d)),$ (i.e., Eq. (20)).
The equation given in Eqs. (30) and (31) are analytically intractable. In this paper, the particle filtering (PF) method is employed to solve Eqs. (30) and (31) recursively [41]. Since the goal of during-mission phase is to perform online updating of MMR, the PF method can also be replaced with unscented Kalman filter method [42] to further improve the efficiency of online updating. Let the prior samples of $X$ and $Y(d)$ at $dc$ be $xcprior=[xc(1),xc(2),⋯,xc(NPF)]$ and where NPF is the number of particles in PF, the posterior samples $xcpost=[xcpost(1),xcpost(2),…,xcpost(NPF)]$ and $ycpost=[ycpost(1),ycpost(2),⋯,ycpost(NPF)]$ are obtained by resampling $xcprior$ and $ycprior$ based on the following weights
$wj=f(vc|xc(j),yc(j))∑i=1NPFf(vc|xc(i),yc(i)),∀j=1,2,…,NPF$
(33)
where $f(vc|xc(i),yc(i)),i=1,2,⋯,NPF$ is the likelihood function computed using Eq. (32).

The posterior samples $xcpost$ will be directly used as the prior samples for $X$ at the next spatial point $dc+1$ since random variables are constant over space. However, $ycpost$ cannot be directly used as prior samples of $Y(d)$ at $dc+1$ due to the spatially varying characteristics of $Y(d)$. In addition, the discontinuity of soil type and slope ID also complicates the generation of prior samples for $Y(d)$ at $dc+1$. This paper develops a numerical procedure (as depicted in Fig. 8) to tackle these challenges.

Fig. 8
Fig. 8
Close modal

#### 3.3.2 Uncertainty Sources Updating.

Let the posterior samples of the finished part of the mission path (as indicated in Fig. 6) be $y1:cf≜{y1post,y2post,…,ycpost}$, where $yipost=[yipost(1),yipost(2),…,$$yipost(NPF)],∀i=1,2,…,c$, the soil properties and slope, $Y(d),∀d∈Ωf$ of the finished part of the mission path be $Y1:cf$, and $Y(d),∀d∈Ωuf$ of the unfinished part of the mission path be $Yc+1:Nduf$, we first partition $y1:cf$ into $y(1)f,y(2)f,⋯,y(Nsoil)f$ if they are soil-type-related or $y(1)f,y(2)f,…,y(Nslope)f$ if they are slope-ID-related. Similarly, we partition $Y1:cf$ and $Yc+1:Nduf$ according the slope ID or soil type as follows:
$Y1:cf⇒Partition{Y(1)f,Y(2)f,…,Y(Nsoil)f}or{Y(1)f,Y(2)f,…,Y(Nslope)f}Yc+1:Nduf⇒Partition{Y(1)uf,Y(2)uf,…,Y(Nsoil)uf}or{Y(1)uf,Y(2)uf,…,Y(Nslope)uf}$
(34)
After the partition of $y1:cf$, $Y1:cf$, and $Yc+1:Nduf$, for the ith slope ID or soil type and given sample of $y1:cf$ as below
$y1:cf(j)≜{y1post(j),y2post(j),…,ycpost(j)}={y(1)f(j)∪⋯∪y(NT)f(j)},j=1,2,…,NPF$
(35)
the conditional mean $μ(i)uf|y(i)f(j)$ and covariance matrix $Σ(i)uf|y(i)f(j)$ of $Y(i)uf|y(i)f(j)$ are obtained as follows:
$μ(i)uf|y(i)f(j)=μuf(i)+Σfu(i)(Σff(i))−1(y(i)f(j)−μf(i))$
(36)
$Σ(i)uf|y(i)f(j)=Σuu(i)−Σfu(i)(Σff(i))−1(Σfu(i))T$
(37)
where i = 1, 2, …, NT, NT = Nsoil if the variable is a soil-type-related variable and NT = Nslope if it is a slope-related variable, $μuf(i)$ is the unconditional mean values of $Y(i)uf$, $μf(i)$ is the unconditional mean values of $Y(i)f$, $Σfu(i)$ is the covariance matrix between $Y(i)f$ and $Y(i)uf$, $Σff(i)$ is the auto-covariance matrix of $Y(i)f$, and $Σuu(i)$ is the auto-covariance matrix of $Y(i)uf$. The covariance matrices are obtained similarly to Eq. (16), according to the coordinates of the spatial points.
Based on $μ(i)uf|y(i)f(j)$ and $Σ(i)uf|y(i)f(j)$, we can then generate a group of samples for $Y(i)uf|y(i)f(j)$ using multivariate Gaussian distribution as follows:
$y(i)uf(j)=[y(i)uf(1,j),y(i)uf(2,j),…,y(i)uf(Nui,j)]$
(38)
where Nui is the number of variables in $Y(i)uf$.

Implementing Eqs. (35) through (38) for all i = 1, 2, …, NT and j = 1, 2, …, NPF, we obtain NPF samples for $Y(1)uf,Y(2)uf,…,Y(Nsoil)uf$ or $Y(1)uf,Y(2)uf,…,Y(Nslope)uf$. The samples are then assembled together according to their coordinates (as indicated in Fig. 8) to form random samples of $Yc+1:Nduf|y1:cf$. Defining the obtained samples as $yijuf,i=1,2,…,NPF;j=1,2,…,(Nd−c)$, $y⋅1uf=[y11uf,y21uf,…,yNPF1uf]$ is then the prior samples of $Y(d)$ at spatial point $dc+1$. By implementing the above procedure recursively, we are able to sequentially estimate the soil properties and vehicle parameters based on the vehicle mobility measurement data.

#### 3.3.3 Mission Mobility Reliability Updating and Remaining Mobilable Distance Estimation.

Since $y1:cf$ are estimated posterior samples of $Y1:cf$ conditioned on $v1:c$, $yijuf,i=1,2,⋯,NPF;j=1,2,⋯,(Nd−c)$ are also samples of $Yc+1:Nduf|v1:c$. As indicated in Fig. 7, the samples of $Yc+1:Nduf|v1:c$ and $xcpost$ are then used as inputs to Eq. (29) to estimate the conditional MMR $R(Ω)|v1:c$ as below
$R(Ω)|v1:c≈1NPF∑i=1NPFId{μV(i,j)≥e,∀j=1,2,…,Nd−c}$
(39)
where μV(i, j) is the mean mobility prediction $μV(xcpost(i),yijuf)$ obtained by plugging $xcpost(i)$ and $yijuf$ into the mobility surrogate model $v=G^(xcpost(i),yijuf)$.
Alternatively, the mission mobility failure probability is computed by
$pf(Ω)|v1:c≈1NPF∑i=1NPFId{μV(i,j)
(40)
In addition to $R(Ω)|v1:c$ and $pf(Ω)|v1:c$, we can also obtain the cumulative density function (CDF) of the remaining mobilable distance (RMD) at current location $dc$ as
$Pr{RMD≤Ne|v1:c}=1NPF∑i=1NPFId{μV(i,j)
(41)
in which NeNdc is a specific threshold of RMD and Id{ · } is an indicator function defined in Eq. (27). The above CDF value presents the probability that the vehicle will be immobilable in the future Ne spatial points conditioned on current vehicle measurements $v1:c$.

As the vehicle proceeds on the mission pass, $R(Ω)|v1:c$ and $Pr{RMD≤Ne|v1:c}$ will be updated. If the updated MMR shows that the $R(Ω)|v1:c$ is too low, the decision-maker needs to decide whether the vehicle should continue this path or a new path needs to be identified based on $Pr{RMD≤Ne|v1:c}$. Next, a case study is used to demonstrate the two phases of the proposed simulation-based MMR analysis.

## 4 Case Study

In this section, an off-road map is used to demonstrate the efficacy of the proposed framework. The application problem is first introduced in Sec. 4.1. Following that, Sec. 4.2 implements the MMR analysis in mission-planning phase, and Sec. 4.3 demonstrates the dynamic updating of MMR during the mission.

### 4.1 Description of the Mission Mobility Reliability Analysis Problem.

A map taken from ArcGIS/ENVI database [43] and US Geological Survey database [44] is selected to demonstrate the proposed framework. The slope and soil maps of the selected area are shown in Fig. 9 accordingly, which are discretized by 1046 points in the coordinate d1 and 707 points in the coordinate d2. As shown in Fig. 9, there are 10 different slope intervals and soil types represented by different IDs, with slope ID zero and soil ID zero represents water. Additionally, each slope ID is associated with a slope parameter, and each soil ID is associated with three different soil parameters, i.e., cohesive strength, friction coefficient, and bulk density. In this paper, all these parameters are considered as stationary Gaussian random field. Table 1 presents statistical information of different parameters corresponding to different slope/soil IDs, which are assumed based on the data obtained from the US Geological Survey database [44].

Fig. 9
Fig. 9
Close modal
Table 1

Distribution and correlation parameters of slope and soil properties

Slope/soil IDSlopeSoil
Cohesive strengthFriction coefficientBulk density
μ**β1*2*μ/σβ12μ/σβ12μ/σβ12
0WaterWater
15/0.23.3/3.00.2/0.0140/390.01/0.00133/320.05/0.00127/28
210/0.53.4/3.62/0.242/430.4/0.0528/331.2/0.123/27
312/0.54.3/2.25/0.441/410.7/0.0228/272.1/0.120/25
414/13.6/4.46/0.241/390.76/0.0231/292.3/0.126/24
516/23.5/5.12/0.140/420.56/0.0133/321.6/0.0427/21
618/22.8/5.011/0.238/370.8/0.0227/272.1/0.0526/28
722/23.2/2.25/0.0239/420.7/0.0328/301.2/0.0224/28
824/14.7/4.04/0.137/390.45/0.00529/281.45/0.0124/27
928/14.5/6.78/0.442/350.78/0.0332/322.35/0.0526/24
Slope/soil IDSlopeSoil
Cohesive strengthFriction coefficientBulk density
μ**β1*2*μ/σβ12μ/σβ12μ/σβ12
0WaterWater
15/0.23.3/3.00.2/0.0140/390.01/0.00133/320.05/0.00127/28
210/0.53.4/3.62/0.242/430.4/0.0528/331.2/0.123/27
312/0.54.3/2.25/0.441/410.7/0.0228/272.1/0.120/25
414/13.6/4.46/0.241/390.76/0.0231/292.3/0.126/24
516/23.5/5.12/0.140/420.56/0.0133/321.6/0.0427/21
618/22.8/5.011/0.238/370.8/0.0227/272.1/0.0526/28
722/23.2/2.25/0.0239/420.7/0.0328/301.2/0.0224/28
824/14.7/4.04/0.137/390.45/0.00529/281.45/0.0124/27
928/14.5/6.78/0.442/350.78/0.0332/322.35/0.0526/24

*Note: μ is the mean, σ is the standard deviation, β1 is the correlation length in coordinate 1, and β2 is the correlation length in coordinate 2.

In this paper, the failure event is defined as the maximum attainable speed is less than 2 m/s (i.e., e = 2 in Eq. (9)). The off-road mobility model is assumed to be
$V(d)=G(X,Y(d))=exp(0.013Y22(d)+0.8cos(π⋅(0.5Y1(d)−1.5Y3(d))/180)+sin(πY4(d)/45))+Y3(d)Y4(d)−0.207$
(42)
where $V(d)$ is the vehicle speed at a spatial point $d$, $Y(d)=[Y1(d),Y2(d),Y3(d),Y4(d)]$ is a vector of terrain parameters which has fur parameters, i.e., slope $(Y1(d))$, cohesive strength $(Y2(d))$, friction coefficient $(Y3(d))$, bulk density $(Y4(d))$. For the sake of illustration, the uncertainty of vehicle parameters $X$ is not considered in this case study.

For the mission planning of the ground vehicle at the selected map, a path is selected from coordinates (250,150) to (900,450) as shown in Fig. 10. As a result, there are 400 spatial points (coordinates) on the selected route. Figure 11 gives the slope IDs and soil IDs of each point.

Fig. 10
Fig. 10
Close modal
Fig. 11
Fig. 11
Close modal

Next, we will demonstrate how to perform surrogate-modeling-based MMR analysis in the mission-planning phase and how to update the MMR in the during-mission phase.

### 4.2 Mission-Planning Phase: Surrogate-Modeling-Based Mission Mobility Reliability Analysis.

As discussed earlier, the parameters $Yi(d),i=1,2,3,4$ are space-dependent random fields. Before estimating MMR, random realizations of these parameters are generated using the KL expansion method using Eqs. (13) through (16). In addition, 15 initial training points are generated for $Yi(d),i=1,2,3,4$ using the approach discussed in Sec. 3.2. After that, an initial surrogate $V=G^(X,Y(d))$ is constructed based on the training points.

The surrogate model is then adaptively refined using the SILK method as discussed in Sec. 3.2. Figure 12 shows the predicted mission mobility failure probability (pf(Ω)) with respect to the added number of function evaluations. The result shows that the SILK-based MMR analysis method only needs additional four function evaluations to get a very accurate estimate of the mission mobility failure probability. In this figure, the true mission mobility failure probability is estimated using the brutal-force MCS method with Ns = 104 samples. This demonstrates the high efficiency of the adaptive surrogate-modeling-based MMR analysis method.

Fig. 12
Fig. 12
Close modal

### 4.3 During-Mission Phase: Dynamic Updating of Mission Mobility Reliability.

From the mission-planning phase, as depicted in Fig. 12, it shows that the mission mobility failure probability of the selected mission path (i.e., Fig. 10) is 0.1022 from both the MCS and adaptive Kriging surrogate-modeling-based method. Even though the failure probability is quite low, rare events of immobility are still possible to happen. To proactively avoid the rare events, MMR is dynamically updated based on the online mobility data as discussed in Sec. 3.3.

In order to demonstrate the dynamic updating process, we first select a realization of the terrain properties that leads to failure of the mission. The selected terrain properties are then assumed to be unknown during the updating process. Figure 13 depicts the vehicle speed corresponds to the selected realization of terrain properties at each point. It shows that if the vehicle proceeds on the mission path (i.e., Fig. 10), an immobility event (i.e., maximum attainable speed <2 m/s) will happen at the 232th point. If the proposed framework can predict the immobility event before the vehicle reaches to the 232-th point, it could save the vehicle from the immobility event.

Fig. 13
Fig. 13
Close modal

In this paper, the slope at each location is assumed to be directly measurable when a vehicle passes a certain point, while the other soil properties are not directly measurable. Using the Bayesian method discussed in Sec. 3.2, we can estimate the posterior distributions of the soil properties as mobility data of the vehicle are collected. Figure 14 shows the estimated posterior distributions of soil parameters at the first point, when the vehicle starts to proceed on the selected mission path. Comparing with the prior distribution, it is observed that the uncertainty in the soil property parameters is reduced. Based on the Bayesian inference at the first point, we can then update the distribution of the remaining points on the mission path by following the procedure as shown in Fig. 8. Based on the updated distributions, we can further calculate the mission mobility failure probability of the remaining part of the mission path. As discussed in Sec. 3.2, the updated distributions of the parameters are also used as the prior distributions for the next spatial point. The aforementioned process continues as the vehicle proceeds on the mission path.

Fig. 14
Fig. 14
Close modal

Figure 15 presents the 95% confidence interval and mean value of the posterior distributions of the soil property parameters after the vehicle proceeds 50 points on the mission path. The results show that the mean values of the posterior distribution of soil parameters are close to the underlying true values which are assumed to be unknown. Following that, Fig. 16 presents the updating history of the mission mobility failure probability as the vehicle proceeds over the mission path. It shows that $pf(Ω)|v1:c$ goes up close to 1 at a point close to the 30th point which is way ahead the 232th point where the failure (i.e., immobility) actually happens. Figure 17 presents the updating history of the RMD over the mission path as the vehicle proceeds on the route. It can see that the proposed dynamic updating framework is able to accurately estimate the RMD as the vehicle proceeds on the route and more observations of the vehicle mobility are collected. This demonstrates the effectiveness of proposed framework in proactively identifying the rare events of immobility before it happens.

Fig. 15
Fig. 15
Close modal
Fig. 16
Fig. 16
Close modal
Fig. 17
Fig. 17
Close modal

## 5 Conclusions and Future Work

This paper proposes a simulation-based MMR analysis framework to account for the uncertainty and improve the prediction accuracy of off-road ground vehicle mobility. The framework has two phases, including a mission-planning phase and a during-mission phase. For the mission-planning phase, a single-loop Kriging surrogate modeling approach is employed to reduce the computational cost by selecting the most critical training points as well as maintaining the prediction accuracy of MMR. For the during-mission phase, a Bayesian updating framework is proposed to reduce the uncertainty of terrain parameters with observed vehicle mobility information. Using the conditional Gaussian process, we can also reduce the uncertainty of the path that has not been accomplished and thus dynamically update the MMR. The dynamic updating of MMR allows us to proactively avoid the rare occurrence of failures before it happens. A case study using the data from ArcGIS/ENVI database and US Geological Survey database demonstrates the efficacy of the proposed framework.

Even though the developed framework is applied to off-ground vehicle MMR analysis, it is also applicable to other type of engineering systems where time-dependent or space-dependent variables are presented. In this paper, only several sources of uncertainty (e.g., slope and soil properties) are considered in the MMR analysis and dynamical updating of MMR. In reality, however, heterogenous sources of uncertainty may affect the mobility prediction. One of the most important uncertainty sources is the model discrepancy between the simulation model and the underlying true physics. The model discrepancy will affect the effectiveness of Bayesian inference and thus the MMR prediction, since the M&S is the basis of MMR analysis and dynamical updating in the proposed framework. In a recent work, a sequential calibration and validation (SeCAV) framework has been proposed to tackle model uncertainty in Bayesian inference [45]. Integrating the SeCAV framework with the proposed framework will make the proposed framework robust to the inadequacy in the M&S and is worth studying in the future. Moreover, this paper focuses on the MMR analysis and MMR updating for any given mission plan. Performing mission-planning optimization based on the developed framework will allow us to identify the optimal mission path that is reliable for the off-road ground vehicle. Reliability-based mission planning is one of our future work.

In addition, synthetic data are used to demonstrate the effectiveness of the proposed framework. The application of the proposed framework to the actual off-road ground vehicles using real experiments is also an important topic that needs to be studied in the future. In this paper, the soil and slope maps are assumed to be correct. Considering the uncertainty in the soil and slope maps is another topic that worth pursuing.

## Acknowledgment

This research was supported by the Automotive Research Center (ARC) in accordance with Cooperative Agreement W56HZV-04-2-0001 U.S. Army CCDC Ground Vehicle Systems Center (GVSC), Warren, MI. The support is gratefully acknowledged.

## Nomenclature

• b=

track width

•
• e=

minimum required maximum attainable speed

•
• m=

number of random field variables related to soil property or slope

•
• $d$=

spatial coordinates with $d=[d1,d2]$ with d1 and d2 represent the coordinate in direction 1 and 2 respectively.

•
• V=

vehicle mobility variable (e.g., the maximum attainable speed in this paper)

•
• W=

weight of vehicle

•
• $X$=

a vector of vehicle parameters

•
• ba=

overall width of the tire inflated but unloaded

•
• da=

overall diameter of the tire inflated but unloaded

•
• dw=

•
• mw=

•
• na=

number of axles

•
• nMCS=

number of MCS samples

•
• pw=

track plate length

•
• $dc$=

current location

•
• $dc−1$=

previous location

•
• $dc+1$=

next location

•
• $di$=

the spatial coordinates at location i, $di=[d1i,d2i]$

•
• $vi$=

mobility data collected at spatial point $di$

•
• $v1:c$=

measured vehicle mobility data of the finished part Ωf

•
• Nc=

number of observations at location $dc$

•
• Nd=

number of discretized spatial points on a mission path

•
• Nij=

number of expansion terms in KL expansion for Yij(·)

•
• Nin=

number of initial training points

•
• NKL=

number of expansion terms in KL expansion

•
• NPF=

number of particles in particle filtering

•
• Nslope=

number of slope IDs

•
• Nsoil=

number of soil types

•
• NT=

NT = Nslope if $Yi(d)$ is a slope-related variable

•
• NT=

Nsoil if $Yi(d)$ is a soil-type-related variable

•
• $y1:cf$=

posterior samples of soil/slope properties of the finished part of the mission path

•
• $Y1:cf$=

$Y(d)$ of the finished part of the mission path

•
• $Y(i)f$=

$Y(d)$ associated with the ith slope ID or soil type in the finished part of the mission path

•
• $Yc+1:Nduf$=

$Y(d)$ of the unfinished part of the mission path

•
• $Y(i)uf$=

$Y(d)$ associated with the ith slope ID or soil type in the unfinished part of the mission path

•
• fs(·, ·)=

a function mapping $u(l)$ to training point $yt(l)$ of $Y(d)$

•
• pf(Ω)=

mission failure probability of a mission path Ω

•
• $RM(d)$=

vehicle mobility reliability at location $d$

•
• G(·, ·)=

vehicle mobility simulation model

•
• Ns(j), j = 1, …, NT=

number of spatial points associated with the jth slope or soil ID on a mission path Ω

•
• R(Ω)=

mission mobility reliability of a mission path Ω

•
• $R(Ω)|v1:c$=

MMR conditioned on $v1:c$

•
• $Yi(d),i=1,…,m$=

the ith space-dependent soil property or terrain parameter

•
• Yij(k), k = 1, 2, …, Ns(j)=

$Yi(d)$ associated with the jth slope ID or soil ID at the kth location

•
• $d(p)Slope,p=1,…,Nslope$=

coordinates corresponding to the pth slope ID

•
• $d(h)Soil,h=1,…,Nsoil$=

coordinates corresponding to the hth soil type

•
• $yt(l),l=1,…,Nin$=

the lth training point of $Y(d)$

•
• $Y(d)$=

a vector of terrain parameters

•
• α(l), l = 1, …, Nin=

the lth training point of soil/slope ID

•
• $α=[α(1),…,α(Nin)]$=

integer training samples of soil/slope IDs

•
• βi, i = 1, 2=

correlation length parameter of coordinate i

•
• $Σfu(i)$=

covariance matrix between $Y(i)f$ and $Y(i)uf$

•
• $Σff(i)$=

auto-covariance matrix of $Y(i)f$

•
• $Σuu(i)$=

auto-covariance matrix of $Y(i)uf$

•
• Ω=

domain of a mission path

•
• Ωf=

spatial points of the finished part of the mission path

•
• Ωuf=

spatial points of the unfinished part of the mission path

## References

1.
Korlath
,
G.
,
2007
, “
Mobility Analysis of Off-Road Vehicles: Benefits for Development, Procurement and Operation
,”
J. Terramech.
,
44
(
5
), pp.
383
393
. 10.1016/j.jterra.2007.10.007
2.
McCullough
,
M.
,
Jayakumar
,
P.
,
Dasch
,
J.
, and
Gorsich
,
D.
,
2017
, “
The Next Generation NATO Reference Mobility Model Development
,”
J. Terramech.
,
73
(
1
), pp.
49
60
. 10.1016/j.jterra.2017.06.002
3.
Ciobotaru
,
T.
,
2009
, “
Semi-empiric Algorithm for Assessment of the Vehicle Mobility
,”
Leonardo Electron. J. Pract. Technol.
,
8
(
15
), pp.
19
30
.
4.
McCullough
,
M.
, and
Haug
,
E.
,
1986
, “
Dynamics of High Mobility Track Vehicles
,”
ASME. J. Mech., Trans., and Automation
,
108
(
2
), pp.
189
196
.
5.
Petrick
,
E.
,
Janosi
,
Z.
, and
Haley
,
P.
,
1981
, “
The Use of the Nato Reference Mobility Model in Military Vehicle Procurement
,” No. 0148-7191, SAE Technical Paper.
6.
Singh
,
A
.,
2012
, “
Mobility of Military Vehicles at TARDEC
,”
Army Tank Automotive Research Development and Engineering Center
,
Warren, MI
.
7.
,
M.
,
Dasch
,
J.
,
Gonzalez
,
R.
,
Hodges
,
H.
,
Jain
,
A.
,
Iagnemma
,
K.
,
Letherwood
,
M.
,
Mccullough
,
M.
,
Priddy
,
J.
, and
Wojtysiak
,
B.
,
2016
,
Next-Generation Nato Reference Mobility Model (Ng-Nrmm)
,
Tank Automotive Research, Development and Engineering Center (TARDEC)
,
Warren, MI
.
8.
McCullough
,
M.
,
Jayakumar
,
P.
,
Dasch
,
J.
, and
Gorsich
,
D.
,
2016
,
Developing the Next Generation NATO Reference Mobility Model
,
US ARMY TARDEC
,
Warren, MI
.
9.
Laughery
,
S.
,
Gerhart
,
G.
, and
Goetz
,
R.
,
1990
,
Bekker’s Terramechanics Model for Off-Road Vehicle Research
,
Tacom Research Development and Engineering Center
,
Warren, MI
.
10.
Recuero
,
A.
,
Serban
,
R.
,
Peterson
,
B.
,
Sugiyama
,
H.
,
Jayakumar
,
P.
, and
Negrut
,
D.
,
2017
, “
A High-Fidelity Approach for Vehicle Mobility Simulation: Nonlinear Finite Element Tires Operating on Granular Material
,”
J. Terramech.
,
72
(
1
), pp.
39
54
. 10.1016/j.jterra.2017.04.002
11.
Serban
,
R.
,
Olsen
,
N.
,
Negrut
,
D.
,
Recuero
,
A.
, and
Jayakumar
,
P.
,
2017
, “
A Co-simulation Framework for High-Performance, High-Fidelity Simulation of Ground Vehicle-Terrain Interaction
,”
North Atlantic Treaty Organization
,
Brussels, Belgium
, Report No. Avt-265, https://www.researchgate.net/publication/317598264_A_co-simulation_framework_for_high-performance_high-fidelity_simulation_of_ground_vehicle-terrain_interaction
12.
Hetherington
,
J.
,
2001
, “
The Applicability of the MMP Concept in Specifying Off-Road Mobility for Wheeled and Tracked Vehicles
,”
J. Terramech.
,
38
(
2
), pp.
63
70
. 10.1016/S0022-4898(00)00010-0
13.
Hu
,
Z.
,
Mourelatos
,
Z. P.
,
Gorsich
,
D.
,
Jayakumar
,
P.
, and
Majcher
,
M.
,
2020
, “
Testing Design Optimization for Uncertainty Reduction in Generating Off-Road Mobility Map Using a Bayesian Approach
,”
ASME J. Mech. Des.
,
142
(
2
), p.
021402
. 10.1115/1.4044111
14.
Lessem
,
A.
,
Ahlvin
,
R.
,
Mason
,
G.
, and
Mlakar
,
P.
,
1992
, “
Stochastic Vehicle Mobility Forecasts Using the NATO Reference Mobility Model. Report 1. Basic Concepts and Procedures
,”
Army Engineer Waterways Experiment Station, Geotechnical Lab
,
Vicksburg, MS
.
15.
Lessem
,
A.
,
Mason
,
G.
, and
Ahlvin
,
R.
,
1996
, “
Stochastic Vehicle Mobility Forecasts Using the NATO Reference Mobility Model
,”
J. Terramech.
,
33
(
6
), pp.
273
280
. 10.1016/S0022-4898(97)00010-4
16.
Priddy
,
J. D.
,
1995
, “
Stochastic Vehicle Mobility Forecasts Using the NATO Reference Mobility Model. Report 3. Database Development for Statistical Analysis of the NRMM II Cross-Country Traction Empirical Relationships
,”
Army Engineer Waterways Experiment Station, Geotechnical Lab
,
Vicksburg, MS
.
17.
Choi
,
K.
,
Jayakumar
,
P.
,
Funk
,
M.
,
Gaul
,
N.
, and
Wasfy
,
T. M.
,
2019
, “
Framework of Reliability-Based Stochastic Mobility Map for Next Generation Nato Reference Mobility Model
,”
ASME J. Comput. Nonlinear Dyn.
,
14
(
2
), p.
021012
. 10.1115/1.4041350
18.
González
,
R.
,
Jayakumar
,
P.
, and
Iagnemma
,
K.
,
2017
, “
Stochastic Mobility Prediction of Ground Vehicles Over Large Spatial Regions: A Geostatistical Approach
,”
Auton. Robots
,
41
(
2
), pp.
311
331
. 10.1007/s10514-015-9527-z
19.
Gonzalez
,
R.
,
Jayakumar
,
P.
, and
Iagnemma
,
K.
,
2017
, “
Generation of Stochastic Mobility Maps for Large-Scale Route Planning of Ground Vehicles: A Case Study
,”
J. Terramech.
,
69
(
1
), pp.
1
11
. 10.1016/j.jterra.2016.10.001
20.
Detweiler
,
Z. R.
, and
Ferris
,
J. B.
,
2010
, “
Interpolation Methods for High-Fidelity Three-Dimensional Terrain Surfaces
,”
J. Terramech.
,
47
(
4
), pp.
209
217
. 10.1016/j.jterra.2010.01.002
21.
Hu
,
Z.
, and
,
S.
,
2016
, “
A Single-Loop Kriging Surrogate Modeling for Time-Dependent Reliability Analysis
,”
ASME J. Mech. Des.
,
138
(
6
), p.
061406
. 10.1115/1.4033428
22.
Maclaurin
,
E.
,
1997
, “
Proposed Revisions to MMP Based on the Results of Tractive Performance Trials With Single Pneumatic Tyres and a Modular Track System
,”
Defence Evaluation and Research Agency
, DERA/LS4/TR970122/1.0.
23.
Priddy
,
J. D.
,
1999
,
Improving the Traction Prediction Capabilities in the NATO Reference Mobility Model (NRMM)
,
Army Engineer Waterways Experiment Station, Geotechnical Lab
,
Vicksburg, MS
.
24.
Jones
,
R. A.
,
Price
,
S. J.
, and
Ahlvin
,
R. B
,
2004
,
Mission Level Mobility Analysis of the US Marine Corps HIMARS Vehicles
,
Engineer Research And Development Center, Geotechnical And Structures Lab
,
Vicksburg, MS
.
25.
Shoop
,
S.
,
Kestler
,
K.
, and
Haehnel
,
R.
,
2006
, “
Finite Element Modeling of Tires on Snow
,”
Tire Sci. Technol.
,
34
(
1
), pp.
2
37
. 10.2346/1.2169827
26.
Huang
,
S.
,
,
S.
, and
Rebba
,
R.
,
2007
, “
Collocation-Based Stochastic Finite Element Analysis for Random Field Problems
,”
Probab. Eng. Mech.
,
22
(
2
), pp.
194
205
. 10.1016/j.probengmech.2006.11.004
27.
Gneiting
,
T.
,
2002
, “
Compactly Supported Correlation Functions
,”
J. Multivar. Anal.
,
83
(
2
), pp.
493
508
. 10.1006/jmva.2001.2056
28.
Hu
,
Z.
, and
Du
,
X.
,
2013
, “
Time-dependent Reliability Analysis With Joint Upcrossing Rates
,”
Struct. Multidiscip. Optim.
,
48
(
5
), pp.
893
907
. 10.1007/s00158-013-0937-2
29.
Zhang
,
D.
, and
Han
,
X.
,
2020
, “
Kinematic Reliability Analysis of Robotic Manipulator
,”
ASME J. Mech. Des.
,
142
(
4
), p.
044502
. 10.1115/1.4044436
30.
Du
,
X.
, and
Hu
,
Z.
,
2012
, “
First Order Reliability Method With Truncated Random Variables
,”
ASME J. Mech. Des.
,
134
(
9
), p.
091005
.10.1115/1.4007150
31.
Hu
,
Z.
, and
Du
,
X.
,
2014
, “
First Order Reliability Method for Time-Variant Problems Using Series Expansions
,”
Struct. Multidiscip. Optim.
,
51
(
1
), pp.
1
21
. 10.1007/s00158-014-1132-9
32.
Mourelatos
,
Z. P.
,
Majcher
,
M.
,
Pandey
,
V.
, and
Baseski
,
I.
,
2015
, “
Time-dependent Reliability Analysis Using the Total Probability Theorem
,”
ASME J. Mech. Des.
,
137
(
3
), p.
031405
. 10.1115/1.4029326
33.
Hu
,
Z.
, and
Du
,
X.
,
2013
, “
A Sampling Approach to Extreme Value Distribution for Time-Dependent Reliability Analysis
,”
ASME J. Mech. Des.
,
135
(
7
), p.
071003
. 10.1115/1.4023925
34.
Liu
,
Y.
,
Zhao
,
Y.
,
Hu
,
Z.
,
Mourelatos
,
Z. P.
, and
,
D.
,
2019
, “
Collision-Avoidance Reliability Analysis of Automated Vehicle Based on Adaptive Surrogate Modeling
,”
ASCE-ASME J. Risk Uncertainty Eng. Syst. Part B: Mech. Eng.
,
5
(
2
), p.
020906
.
35.
Jiang
,
C.
,
Qiu
,
H.
,
Gao
,
L.
,
Wang
,
D.
,
Yang
,
Z.
, and
Chen
,
L.
,
2020
, “
Real-time Estimation Error-Guided Active Learning Kriging Method for Time-Dependent Reliability Analysis
,”
Appl. Math. Model.
,
77
(
1
), pp.
82
98
. 10.1016/j.apm.2019.06.035
36.
Hu
,
Z.
, and
Du
,
X.
,
2015
, “
Mixed Efficient Global Optimization for Time-Dependent Reliability Analysis
,”
ASME J. Mech. Des.
,
137
(
5
), p.
051401
. 10.1115/1.4029520
37.
Wang
,
Z.
, and
Wang
,
P.
,
2013
, “
A new Approach for Reliability Analysis With Time-Variant Performance Characteristics
,”
Reliab. Eng. Syst. Saf.
,
115
(
1
), pp.
70
81
. 10.1016/j.ress.2013.02.017
38.
Zhang
,
D.
,
Han
,
X.
,
Jiang
,
C.
,
Liu
,
J.
, and
Li
,
Q.
,
2017
, “
Time-dependent Reliability Analysis Through Response Surface Method
,”
ASME J. Mech. Des.
,
139
(
4
), p.
041404
. 10.1115/1.4035860
39.
Stein
,
M.
,
1987
, “
Large Sample Properties of Simulations Using Latin Hypercube Sampling
,”
Technometrics
,
29
(
2
), pp.
143
151
. 10.1080/00401706.1987.10488205
40.
Simpson
,
T. W.
,
Mauery
,
T. M.
,
Korte
,
J. J.
, and
Mistree
,
F.
,
2001
, “
Kriging Models for Global Approximation in Simulation-Based Multidisciplinary Design Optimization
,”
AIAA J.
,
39
(
12
), pp.
2233
2241
. 10.2514/2.1234
41.
Arulampalam
,
M. S.
,
,
S.
,
Gordon
,
N.
, and
Clapp
,
T.
,
2002
, “
A Tutorial on Particle Filters for Online Nonlinear/non-Gaussian Bayesian Tracking
,”
IEEE Trans. Signal Process.
,
50
(
2
), pp.
174
188
. 10.1109/78.978374
42.
Wan
,
E. A.
, and
Van Der Merwe
,
R.
,
2000
, “
The Unscented Kalman Filter for Nonlinear Estimation
,”
Proceedings of the IEEE 2000 Adaptive Systems for Signal Processing, Communications, and Control Symposium (Cat. No. 00EX373)
,
,
Oct. 4
, IEEE, pp.
153
158
.
43.
Clausen
,
A.
,
Wang
,
F.
,
Jensen
,
J. S.
,
Sigmund
,
O.
, and
Lewis
,
J. A.
,
2015
, “
Topology Optimized Architectures With Programmable Poisson’s Ratio Over Large Deformations
,”
,
27
(
37
), pp.
5523
5527
44.
Larsen
,
U. D.
,
Signund
,
O.
, and
Bouwsta
,
S.
,
1997
, “
Design and Fabrication of Compliant Micromechanisms and Structures With Negative Poisson’s Ratio
,”
J. Microelectromech. Syst.
,
6
(
2
), pp.
99
106
. 10.1109/84.585787
45.
Jiang
,
C.
,
Hu
,
Z.
,
Liu
,
Y.
,
Mourelatos
,
Z.
,
Gorsich
,
D.
, and
Jayakumar
,
P.
,
2020
, “
A Sequential Calibration and Validation Framework for Model Uncertainty Quantification and Reduction
,”
Comput. Methods Appl. Mech. Eng.
,
368
, p.
113172
. 10.1016/j.cma.2020.113172