1 Introduction↩︎

Robotic exploration missions increasingly rely on autonomous systems to gather information where human access is limited or impossible, such as planetary surface exploration [1], environmental monitoring [2], and precision agriculture [3]. In these settings, robots must continually decide where to move and what data to collect under strict resource limits. Informative Path Planning (IPP) addresses this by finding paths that maximize information gain under travel-cost constraints.

In classical IPP(C-IPP) formulations, sensing is modeled as a digital operation: measurements such as images, temperature readings, or signal density are acquired without altering the robot’s physical state. Under this assumption, the amount of data collected is decoupled from the budget constraint, since traversing an edge while carrying many samples costs the same as traversing it while carrying none. As a result, energy and distance are not distinguished, and minimizing path length often serves as a reasonable surrogate for minimizing energy consumption. However, this abstraction breaks down in missions involving physical sample acquisition. Consider a lunar rover collecting regolith samples for laboratory analysis. Each additional sample increases the rover’s payload, raising the energy cost of all subsequent motion. The same effect arises for surface vehicles accumulating water samples or drones retrieving agricultural specimens for field inspection. In each of these settings, information is acquired in the form of physical samples rather than purely digital observations, and each collected sample modifies the robot’s physical state in a way that persists over the remainder of the mission.

In addition, many real-world missions collect multiple samples at scientifically significant locations, adding a further decision dimension: how much to collect at each site. In C-IPPthis choice would be inconsequential, since sampling does not affect traversal cost; under load-dependent energy, however, it fundamentally alters the planning problem, as heavier collection raises the cost of every subsequent edge. Routing, visitation order, and sampling amount thus become tightly intertwined and must be optimized jointly. To this end, we propose Load-aware IPP(LIPP), a strict generalization of C-IPPthat explicitly models the impact of physical sample mass on energy consumption.

The main contributions of this paper are as follows:

• We propose LIPP, a new class of IPP problem in which the amount of physical sampling affects estimation uncertainty and induces load-dependent, order-sensitive traversal costs.

• We derive an exact MIQP reformulation of LIPP.

• We derive theoretical bounds quantifying the path-length loss of LIPPrelative to distance-optimal solutions.

• We validate LIPPacross 2,000 simulated scenarios, showing up to \(3\times\) higher uncertainty reduction per unit energy than C-IPP. Our implementation code are publicly available at https://github.com/hojunieee/lipp-load-aware-ipp.

2 Related Work↩︎

IPPhas been widely applied to environmental monitoring, spatial field estimation, and target tracking, where paths are optimized under resource constraints to reduce estimation uncertainty or to accurately recover distributional statistics (e.g., quantiles) of the underlying field [4]–[10]. These methods differ chiefly in how the path is optimized. Discrete, optimization-based formulations operate on a graph and select sensing locations under a travel budget, including the sparse-regression mixed-integer program of Dutta et al.[11] that we build upon, while continuous-space methods optimize waypoints directly by gradient descent, using sparse Gaussian processes to scale to large environments under distance and kinematic budgets [12]. A complementary line addresses the adaptive, online setting, replanning as measurements arrive via deep reinforcement learning [6] or learned attention-based kernels [13]. Closest to our setting, Ott et al.[14] couple sensing and motion by letting the agent trade sensor accuracy against energy cost under a shared budget. Across all of these, however, measurements are digital: an observation incurs a fixed, action-local cost and leaves the robot’s physical state unchanged, so sensing and mobility remain weakly coupled and traversal cost is governed only by path geometry.

Separately, routing problems with load-dependent travel costs have been studied in the operations research literature. Representative examples include pickup-and-delivery problems, vehicle routing with load-dependent fuel consumption, and the Traveling Thief Problem (TTP), where item selection directly influences travel speed or cost [15]–[17]. These formulations explicitly capture the coupling between carried load and traversal effort, typically balancing transportation cost against profit or service objectives. However, in these settings, energy or fuel consumption is treated purely as an operational cost accumulated deterministically along a route, rather than as part of a sequential decision process that trades off information gain and future mobility. Moreover, they generally assume a fully known environment with deterministic rewards, and their objective is cost–profit optimization rather than uncertainty reduction. Consequently, they do not incorporate adaptive information gathering or the estimation-theoretic structure that underpins informative path planning.

Existing IPP methods do not jointly model repeated physical sampling and cumulative load-dependent traversal energy, while routing formulations with load-dependent travel costs do not incorporate probabilistic field estimation or information-theoretic objectives. Our work bridges these two directions by integrating load-dependent energy consumption into informative path planning for spatial field estimation, yielding a tightly coupled optimization problem that simultaneously determines sampling amount, traversal order, and uncertainty reduction.

3 Problem Formulation↩︎

3.1 Problem Setup↩︎

Consider a robot capable of collecting physical samples operating within a closed and bounded region \(\mathcal{X} \subset \mathbb{R}^2\). We discretize \(\mathcal{X}\) into a weighted directed graph \(G := (V, \mathcal{E}, d)\), where \(V\) is a set of \(n\) admissible sampling locations, \(\mathcal{E} \subseteq V \times V\) is a set of directed edges, and \(d : \mathcal{E} \rightarrow \mathbb{R}^{+}\) denotes the traversal cost per unit carried load along each directed edge. The cost function \(d\) may incorporate environment-dependent factors such as terrain slope or surface roughness for planetary rovers, as well as wind or water currents for aerial and aquatic robots.

In addition to the sampling vertices \(V\), we define a set of \(m\) test locations \(T \subset \mathcal{X}\) at which we aim to estimate an unknown static scalar field \(f : \mathcal{X} \rightarrow \mathbb{R}\). We model the field \(f\) as a Gaussian Process (GP), \(f \sim \mathcal{GP}(0, k),\) where \(k : \mathcal{X} \times \mathcal{X} \rightarrow \mathbb{R}\) is a known positive definite covariance kernel, similar to many GP-based IPP formulations [18].

The robot is assigned a start vertex \(s \in V\) and a goal vertex \(g \in V\). The planner produces a path for the robot, \(P = \bigl\langle (v_{1},\, l_{1}), (v_{2},\, l_{2}), \dots, (v_{p},\, l_{p}) \bigr\rangle,\) where \(v_{1} = s\), \(v_{p} = g\), \((v_{j}, v_{j+1}) \in \mathcal{E}\) for all \(j = 1,\dots,p-1\), and \(l_{j} \in \{0,1,\dots,S_{\max}\}\) denotes the number of unit-mass physical samples collected at vertex \(v_{j}\).
Problem 1 Load-Aware IPP(LIPP)
Given a set of test locations \(T\), determine a path \(P\) consisting of \((\text{vertex},\, \text{sample amount})\) tuples that minimizes the posterior uncertainty of the test locations subject to an accumulated energy budget: \[\begin{align} \min_{P} \quad & \operatorname{PostVar}_{T}(P) \\ \text{s.t.} \quad & \operatorname{Energy}(P) \le B. \end{align}\] Here, \(B\) denotes the available energy budget, and \(\operatorname{PostVar}_{T}(P)\) represents the trace of the posterior covariance matrix of the GP evaluated at \(T\) after collecting samples along path \(P\).

3.2 Measurement Noise and Energy Model↩︎

Conventional IPP formulations assume that the robot acquires a single digital noisy measurement \(y_j\) at each visited vertex \(v_j\), with constant measurement noise variance: \[y_j = f(v_j) + \eta_j, \qquad \eta_j \sim \mathcal{N}\!\left(0, \sigma^2\right), \label{eq:measurement95model951}\tag{1}\] where \(\sigma^2\) denotes the variance of the measurement noise associated with each sample.

In realistic sampling settings, where repeated sampling is allowed, collecting multiple independent samples at the same vertex can reduce the effective measurement noise variance due to the variance-reduction effect of averaging independent observations [19]. To capture this effect, we adopt the following measurement model: if \(l_j \in \mathbb{Z}_{\ge 0}\) samples are collected at vertex \(v_j\), then the resulting measurement noise is modeled as \(\eta_j \sim \mathcal{N}\!\left(0, \frac{\sigma^2}{l_j}\right).\) However, acquiring additional samples also increases the load carried by the robot. Consider a path consisting of an ordered sequence of vertices \(\{v_1, \dots, v_p\}\), where \(l_j \in \{0,1,\dots,S_{\max}\}\) denotes the number of unit-mass samples collected at vertex \(v_j\). Let \(\lambda\) denote the unit mass of a sample and \(R_0\) denote the mass of the robot without any samples. Then, the cumulative mass of the robot at the \(i\)-th vertex is defined as \(R_i := \sum_{j=1}^{i} (\lambda l_j) + R_0.\) Using this, we define the total energy expenditure along the path as \[E := \sum_{j=1}^{p-1} d(v_j, v_{j+1}) \, R_j,\] where \(d(v_j, v_{j+1})\) denotes the travel cost between consecutive vertices. This formulation provides a first-order approximation of load-dependent energy consumption based on the relation \(E = \text{force} \times \text{distance}\), assuming that the dominant traversal forces (e.g., friction and gravity) scale linearly with the robot mass.

Note that the proposed load-aware formulation strictly generalizes C-IPP. In particular, when the robot base mass is normalized to \(R_0 = 1\) and the sample unit mass approaches zero \(\lambda \to 0\), the cumulative mass satisfies \(R_j \to R_0 = 1\) for all \(j\). In this limit, the energy model reduces to \(E(P) = \sum_{j=1}^{p-1} d_j\), and the energy budget constraint becomes equivalent to the standard distance budget constraint used in C-IPP: \(\sum_{e} d(e)\chi_e \le b\). Therefore, C-IPP can be viewed as a special case of the proposed LIPP formulation in which sampling does not affect mobility.

3.3 Objective Function↩︎

We quantify information using the posterior variance [11] over the test set, \[\operatorname{PostVar}_{T}^M(P) = \operatorname{trace}\!\left( M \, \bar{k}_{TT} \right), \label{eq:postvar951}\tag{2}\] where \(M \in \mathbb{R}^{m \times m}\) is a diagonal weight matrix representing the importance of each test vertex. The posterior covariance matrix given by kernel is \[\bar{k}_{TT} = k_{TT} - k_{TV} \left( k_{VV} + N \right)^{-1} k_{VT},\] where \(k_{VV} \in \mathbb{R}^{n \times n}, \quad k_{TV} \in \mathbb{R}^{m \times n},\) and \(N \in \mathbb{R}^{n \times n}\) is a diagonal noise matrix whose entries depend on the sampling amount decision \(l_j\), and its diagonal element is defined as \(N_{jj}=\frac{\sigma^2}{l_j}\) using our new measurement noise model.

However, the objective in 2 involves a matrix inverse whose entries depend nonlinearly on the discrete sampling decisions through \(N\), resulting in a highly nonlinear and nonconvex problem not directly amenable to standard optimization solvers.

To address this, we leverage the equivalence between the GP posterior mean and the linear least-squares estimator (LLSE) under Gaussian assumptions (Theorem 5 in [11]). Crucially, this is not a linear approximation but an exact algebraic identity: under jointly Gaussian assumptions, the LLSE coincides with the conditional expectation, so the reformulation preserves the original objective exactly. The key insight is that the optimal GP predictor at each test location is a linear combination of the observations at visited vertices; collecting this weight matrix into \(A\in\mathbb{R}^{m\times n}\), the posterior covariance can be written as a polynomial in \(A\) and \(N\) without the matrix inverse.

Specifically, we introduce a mixed-integer formulation of Problem 1 as follows.

\[\begin{align} \min_{A,N} \quad & \operatorname{tr}\!\left( M \bigl(A(k_{VV}+N)A^\top - 2k_{TV}A^\top + k_{TT}\bigr) \right) \tag{3} \\[4pt] \text{s.t.}\quad & \sum_{e\in\mathcal{E}_v^{\mathrm{in}}} \chi_e = \sum_{e\in\mathcal{E}_v^{\mathrm{out}}} \chi_e \le 1 \tag{4} \\ & \sum_{e\in\mathcal{E}_s^{\mathrm{out}}} \chi_e = 1, \quad \sum_{e\in\mathcal{E}_g^{\mathrm{in}}} \chi_e = 1 \tag{5} \\ & \sum_{e\in\mathcal{E}_s^{\mathrm{in}}} \chi_e = 0, \quad \sum_{e\in\mathcal{E}_g^{\mathrm{out}}} \chi_e = 0 \tag{6} \\ & y_v = \sum_{e\in\mathcal{E}_v^{\mathrm{in}}} \chi_e \tag{7} \\ & -A_{\max} y_v \le A_{t,v} \le A_{\max} y_v \tag{8} \\ & o_s = 0 \tag{9} \\ & 0 \le o_v \le |V|-1 \tag{10} \\ & o_v \ge o_u + 1 - \mathcal{M}^{(o)}(1-\chi_{uv}) \tag{11} \\ & L_s = \lambda \sum_{c=1}^{S_{\max}} c\, z_{s,c} \tag{12} \\ & L_v \ge L_u + \lambda \sum_{c=1}^{S_{\max}} c\,z_{u,c} - \mathcal{M}^{(L)}(1-\chi_{uv}) \tag{13} \\ & 0 \le L_v \le L_{\max} \tag{14} \\ & R_v = R_0 + L_v \tag{15} \\ & l_v = \sum_{c=1}^{S_{\max}} c\,z_{v,c} \tag{16} \\ & \sum_{c=1}^{S_{\max}} z_{v,c} = y_v \tag{17} \\ & \sum_{(u,v)\in\mathcal{E}} d_{uv} \,R_u\, \chi_{uv} \le B \tag{18} \\ & \chi_e \in \{0,1\}, \quad e\in\mathcal{E} \tag{19} \\ & y_v \in \{0,1\}, \quad v\in V \tag{20} \\ & z_{v,c} \in \{0,1\}, \quad v\in V,\; c=1,\dots,S_{\max} \tag{21} \\ & A \in \mathbb{R}^{m\times n},\quad l_v\in\mathbb{Z}_{\ge0},\quad L_v\in\mathbb{R}_{\ge0},\quad R_v\in\mathbb{R}_{\ge0}. \tag{22} \end{align}\] where constraints 4 and 7 hold for all \(v\in V\setminus\{s,g\}\), constraint 8 holds for all \(t\in T\) and \(v\in V\setminus\{s,g\}\), and constraints 11 and 13 hold for all \((u,v)\in\mathcal{E}\).

In this formulation, \(A \in \mathbb{R}^{m\times n}\) is a continuous decision variable representing the linear estimator matrix obtained from the LLSE reformulation. The discrete decision variables are: (i) binary edge variables \(\chi_e\) indicating whether edge \(e\) is selected, (ii) binary vertex variables \(y_v\) indicating whether vertex \(v\) is visited, (iii) binary sampling variables \(z_{v,c}\) indicating whether \(c\) unit samples are collected at vertex \(v\), (iv) integer order variables \(o_v\) representing the position of vertex \(v\) along the path, (v) integer sample count variable \(l_v\) representing the number of samples taken at vertex \(v\), (vi) continuous load variables \(L_v\) representing the mass of the cumulative sample load carried upon departing vertex \(v\), and (vii) continuous mass variables \(R_v\) representing the mass of the robot upon departing vertex \(v\).

Constraints 4 –7 encode the flow and path feasibility requirements. Specifically, 4 enforces flow conservation at intermediate vertices, 5 and 6 ensure a single outgoing edge from the start vertex and a single incoming edge to the terminal vertex, and 7 maintains consistency between vertex and edge activation variables. Constraint 8 links the estimator variables to the visitation decisions via a linear big-\(M\) formulation, enforcing \(A_{t,v}=0\) whenever vertex \(v\) is not selected (i.e., \(y_v=0\)). Thus, estimator coefficients are active only at visited vertices.

We introduce ordering constraints to eliminate subtours and capture the order-dependent nature of the load–energy coupling induced by repeated physical sampling. Constraint 9 fixes the order of the start vertex to zero, while 10 bounds the order variable within \(0\) and \(|V|-1\). Constraint 11 enforces precedence along selected edges: if edge \((u,v)\) is active, then the order of vertex \(v\) must be at least one greater than that of \(u\). Here, \(\mathcal{M}^{(o)}\) denotes a sufficiently large big-\(M\) constant, e.g., \(|V|\). Together, constraints 9 –11 form a Miller–Tucker–Zemlin (MTZ) subtour elimination scheme [20], preventing disconnected cycles and ensuring a single directed path from \(s\) to \(g\).

We next introduce the load and energy constraints. Constraint 12 initializes the sample load at the start vertex to the load collected at the first vertex. Constraint 13 ensures that the cumulative load at each visited vertex equals the load carried upon arrival plus the samples collected at that vertex. The Constraint 14 bounds the load variable within its feasible range, and Constraint 15 defines the mass of the robot using the cumulative sample load and its initial mass. Constraint 16 defines the integer sampling amount \(l_v\) at vertex \(v\) in terms of the binary variables \(z_{v,c}\), which is used in the noise matrix \(N\) in the objective. Constraint 17 ensures consistency between vertex activation and sampling decisions by enforcing that exactly one sampling level is selected if and only if the vertex is visited. Finally, constraint 18 imposes the total energy budget \(B\), where the energy expenditure depends on both traversal distance and cumulative robot mass along the path.

With this formulation, we explicitly capture the coupling between uncertainty reduction from repeated sampling and the order-dependent energy cost induced by transporting the accumulated physical load. Nevertheless, the resulting optimization problem remains computationally challenging. The objective function is nonconvex due to the cubic term \(A (k_{VV}+N) A^\top\), and the energy budget constraint 18 contains bilinear terms in \(R_u\) and \(\chi_{uv}\). Consequently, the resulting formulation is a mixed-integer nonconvex program and is NP-hard. To enable an efficient solution using a commercial solver such as Gurobi, we reformulate the problem into an equivalent MIQP: the cubic terms in the objective are reduced to quadratic form via variable disaggregation (auxiliary variables \(A_{t,v,c}\)), and the bilinear energy constraint is exactly linearized via a McCormick envelope.

4 Mixed Integer Quadratic Programming Reformulation↩︎

The objective in 3 contains the cubic interaction \(z_{v,c}\,A_{t,v}^{2}\) arising from the product \(A\,N\,A^{\top}\), where the diagonal noise matrix \(N\) depends on the discrete sampling decisions through \(z_{v,c}\). To reduce this to a quadratic program, we disaggregate each estimator coefficient into per-sampling-level components by introducing auxiliary continuous variables \[A_{t,v,c}\in\mathbb{R}, \qquad t\in T,\;v\in\mathcal{V},\;c=1,\dots,S_{\max},\] with the aggregation and big-\(M\) linking constraints \[A_{t,v} = \sum_{c=1}^{S_{\max}} A_{t,v,c}, \qquad -A_{\max}\,z_{v,c} \;\le\; A_{t,v,c} \;\le\; A_{\max}\,z_{v,c}. \label{eq:A95aggregation95bigM}\tag{23}\] These constraints ensure \(A_{t,v,c}=0\) whenever \(z_{v,c}=0\), linking the estimator coefficients to the discrete sampling choices.

Because both \(M\) and \(N\) are diagonal, the trace objective \(\operatorname{tr}\!\bigl(M(A(k_{VV}+N)A^{\!\top} -2\,k_{TV}A^{\!\top}+k_{TT})\bigr)\) decomposes into a sum over test locations \(t\). Expanding and substituting \(N_{vv}=\sum_{c=1}^{S_{\max}}\frac{\sigma^{2}}{c}\,z_{v,c}\) yields cubic terms of the form \(\frac{\sigma^{2}}{c}\,z_{v,c}\,A_{t,v}^{2}\). By 23 , each such term is equivalently expressed as the quadratic \(\frac{\sigma^{2}}{c}\,A_{t,v,c}^{2}\), since \(A_{t,v,c}\) is nonzero only when \(z_{v,c}=1\). The resulting MIQP objective is \[\begin{align} \min_{A_{\cdot,\cdot,\cdot},\,z}\; \sum_{t=1}^{m} M_{tt} \Big(\, &\sum_{v_1,v_2\in\mathcal{V}} k_{VV}(v_1,v_2)\, A_{t,v_1}\,A_{t,v_2} \nonumber\\ &+\sum_{v\in\mathcal{V}} \sum_{c=1}^{S_{\max}} \frac{\sigma^2}{c}\,A_{t,v,c}^{2} \nonumber\\ &-2\!\sum_{v\in\mathcal{V}} k_{TV}(t,v)\,A_{t,v} \;+\; k_{TT}(t,t) \Big). \label{eq:MIQP95objective} \end{align}\tag{24}\]

All constraints 4 –21 carry over unchanged, except that the bilinear energy budget constraint 18 is replaced by its exact McCormick linearization [21]. Introducing auxiliary variables \(w_{uv}\in\mathbb{R}_{\ge0}\) for every \((u,v)\in\mathcal{E}\), we write \[\begin{align} \sum_{(u,v)\in\mathcal{E}} d_{uv}\, w_{uv} &\;\le\; B, \tag{25} \\[1pt] w_{uv} &\le R_u, \qquad w_{uv} \le R_{\max}\,\chi_{uv}, \nonumber\\ w_{uv} &\ge R_u - R_{\max}\,(1-\chi_{uv}), \qquad w_{uv} \ge 0, \tag{26} \end{align}\] where \(R_{\max}=R_0+L_{\max}\). When \(\chi_{uv}=1\) the envelope forces \(w_{uv}=R_u\); when \(\chi_{uv}=0\) it forces \(w_{uv}=0\). Since \(R_u\in[R_0,R_{\max}]\) is bounded and \(\chi_{uv}\) is binary, the McCormick relaxation is exact, and the overall formulation remains an MIQP solvable by an off-the-shelf solver such as Gurobi.

5 Analysis↩︎

Incorporating repeated physical sampling and load-aware energy modeling fundamentally changes the structure of the IPP problem, enabling improved energy efficiency while maintaining low posterior variance. These benefits, however, introduce trade-offs, including potentially longer geometric path lengths and increased computational complexity. Although path length no longer directly corresponds to total energy expenditure in the LIPP setting, it remains practically relevant as a proxy for mission execution time. In this section, we provide a theoretical analysis of both trade-offs.

5.1 Theoretical Bound on Execution Path Length↩︎

We bound the execution path length of LIPP relative to C-IPP, where C-IPP denotes 3 –8 with \(N_{jj} = \sigma^2\) and the distance constraint \(\sum_{e} d(e)\chi_e \le b\), as in [11].

Let \(P_D\) and \(P_E\) denote the paths generated by C-IPP and LIPP, respectively, with \(p_D := |P_D|\) and \(p_E := |P_E|\) vertices, and let \(S_{\max}\) be the maximum number of unit samples per vertex. We assume that the LIPP solution consumes no more energy than the C-IPP solution, i.e., \(E(P_E) \le E(P_D)\). For a path \(P = \langle v_1, \dots, v_p \rangle\), define \[E(P) := \sum_{j=1}^{p-1} R_j d_j, \qquad D(P) := \sum_{j=1}^{p-1} d_j,\] where \(d_j\) is the travel cost of edge \((v_j, v_{j+1})\) and \(R_j\) is the robot mass while traversing it. Since each visited vertex contributes between \(1\) and \(S_{\max}\) unit samples of mass \(\lambda > 0\) to the base mass \(R_0\), \[R_0 + \lambda j \;\le\; R_j \;\le\; R_0 + \lambda S_{\max} j,\] and hence, for any path \(P\), \[\sum_{j=1}^{p-1} (R_0 + \lambda j)\, d_j \;\le\; E(P) \;\le\; \sum_{j=1}^{p-1} (R_0 + \lambda S_{\max} j)\, d_j . \label{eq:energy95sandwich}\tag{27}\]

Since \(1 \le j \le p-1\), bounding \(j\) termwise in 27 and using \(E(P_E) \le E(P_D)\) gives \[\begin{align} (R_0 + \lambda)\, D(P_E) &\le E(P_E) \le E(P_D) \\ &\le \big(R_0 + \lambda S_{\max} (p_D - 1)\big) D(P_D). \end{align}\] Dividing by \(R_0 + \lambda > 0\) yields \[D(P_E) \;\le\; C\, D(P_D), \qquad C := \frac{R_0 + \lambda S_{\max} (p_D - 1)}{R_0 + \lambda}, \label{eq:general95bound}\tag{28}\] which holds for any arbitrary graph structure and vertex selection.

If the graph lives on an equidistant grid, i.e., \(d_j = D(P)/(p-1)\) along each path, then \(\sum_{j=1}^{p-1} j\, d_j = \tfrac{p}{2} D(P)\) exactly, and 27 sharpens 28 to \[\frac{D(P_E)}{D(P_D)} \;\le\; \frac{R_0 + \lambda S_{\max} \tfrac{p_D}{2}}{R_0 + \lambda \tfrac{p_E}{2}}. \label{eq:grid95bound}\tag{29}\]

Finally, although we have defined a theoretical bound, the path length can also be controlled explicitly by adding the linear constraint \(\sum_{(u,v)\in\mathcal{E}} d_{uv}\,\chi_{uv} \le b\) for a user-specified budget \(b\). This limits execution time while still optimizing load-dependent energy, without increasing the number of decision variables. We thus obtain theoretical execution-time bounds for arbitrary graphs, together with a simple mechanism for regulating execution time directly.

5.2 Computational Complexity↩︎

Both formulations result in mixed-integer quadratic programs (MIQPs), but they differ in relaxation strength and solver behavior. In terms of nominal problem size, the C-IPP formulation scales as \(O(m n^2)\), whereas the LIPP formulation introduces additional sampling-indexed variables and bilinear terms, increasing the size to \(O(m n^2 S_{\max}^2)\).

In practice, however, the computational gap is driven less by polynomial growth and more by differences in relaxation quality. In the classical formulation, routing and sampling variables are largely separable, leading to comparatively tighter continuous relaxations. In contrast, the LIPP formulation—through the bilinear coupling introduced in Constraint 26 —links load accumulation, sampling decisions, and traversal order. When integrality constraints are relaxed, these interactions permit fractional routing and load propagation, which can weaken the resulting lower bound after linearization (e.g., via McCormick envelopes) [22]. A weaker root relaxation forces the branch-and-bound solver to explore a larger portion of the search tree to certify optimality. The empirical solve-time trends in Fig. 5 suggest this behavior is realized in practice as problem size grows.

6 Experimental Results↩︎

In this section, we empirically validate four properties of LIPPagainst C-IPPand a Greedy baseline over 2,000 randomly generated scenarios spanning various graph densities and sizes: (i) the added sampling-amount dimension changes the optimal vertex set, not merely visitation order; (ii) LIPPcoincides with C-IPPat \(\lambda=0\) and outperforms it as \(\lambda\) grows; (iii) it achieves comparable paths at lower energy under reasonable budgets; and (iv) its overhead stems mainly from weaker relaxations.

6.1 Experiment Setup↩︎

The experiment emulates a lunar rover mission collecting regolith samples for subsequent estimation of thorium (Th) surface concentration, which provides critical geochemical information such as crustal composition and evolution history [23]. Terrain elevation is incorporated by modeling traversal cost as a function of elevation difference between vertices, \(d_{uv} = d_{\mathrm{euclid}} \cdot (1 + \alpha(\mathrm{height}_v - \mathrm{height}_u))\), where \(\alpha\) is a constant scaling factor. The total energy to traverse an edge is the traversal cost multiplied by the robot’s accumulated mass at the start of that edge. The elevation profile and Th concentration field are held fixed across all trials. The C-IPPformulation follows Section V-A, and the Greedy heuristic selects at each step the vertex that maximizes posterior variance reduction normalized by travel distance. For C-IPPand Greedy, we evaluate uniform sample counts of 1, 2, and 3 per visited vertex, denoted by the suffixes “_S1”, “_S2”, and “_S3”, respectively. For LIPP, the suffix “_Bk” indicates an energy budget of \(k\) units.

6.2 Qualitative Results↩︎

To validate property (i), Fig. 2 compares the paths each method produces on the same instance. The Greedy algorithm recursively selects the vertex with the greatest distance-normalized posterior-variance reduction, yielding a myopic path that exhausts its budget before reaching vertex d. C-IPPplans over a longer horizon and, despite the tight distance budget, routes to include vertex d; its distance constraint, however, prevents it from also sampling at vertex e. LIPPdefers heavy sampling to the end of the path: it collects a single unit at vertices d, c, and g, far from the test locations, then increases the sample count in the important region around vertices f and e, where the test locations are clustered. It achieves this while expending roughly half the energy of C-IPP. Crucially, because sampling amount is an added decision dimension, LIPPdoes more than redistribute samples over a fixed vertex set—it selects a different set of vertices altogether, altering both sampling allocation and routing rather than merely permuting visitation order.

6.3 Quantitative Results↩︎

To validate property (ii), i.e., that LIPPrelaxes to C-IPPas discussed in Section III-B, we set the energy budget of LIPPto \(B = 2.0\), the distance budget of C-IPPto \(b = 2.0\), the robot mass to \(R_0 = 1.0\), and sweep the unit sample mass \(\lambda\) from \(0\) to \(1\). Fig. 3 (a) and 3 (c) show that as \(\lambda \to 0\), the posterior variance reduction and efficiency of LIPPconverge to those of C-IPPwith \(S = S_{\max}\). This matches the theoretical equivalence in Section III-B: when samples are weightless, LIPPnaturally collects the maximum number at every visited vertex. The posterior variance reduction of C-IPP, by contrast, stays flat as sample mass grows, showing that its path selection is entirely decoupled from the physical cost of carrying samples. This decoupling drives a steep rise in total energy consumption (Fig. 3 (b)). LIPPinstead explicitly accounts for the growing load through its joint optimization (Section IV), adapting routing order and sampling amount to the energy budget. At \(\lambda = 1.0\), Fig. 3 (c) shows LIPPattaining roughly three times the posterior variance reduction per unit energy of the C-IPPbaseline under uniform sampling with \(S_{\max} = 3\). This demonstrates that explicitly modeling load-dependent coupling yields substantial efficiency gains when the physical burden of carrying samples is non-negligible.

Regarding property (iii), like any informative path planning problem, LIPP’s behavior depends heavily on the budget it is given. To examine how much path length LIPPsacrifices to save energy, we test three energy budgets, each defined relative to \(B_{\mathrm{CIPP}}\), the energy C-IPPconsumes on the same graph: an excessive budget equal to the full \(B_{\mathrm{CIPP}}\), a moderate budget of \(0.6\,B_{\mathrm{CIPP}}\), and a tight budget of \(0.35\,B_{\mathrm{CIPP}}\). Given the full \(B_{\mathrm{CIPP}}\) (Fig. 4 (a)), LIPPhas enough energy to follow the same path as C-IPPand collect the maximum number of samples at each vertex; instead, it tends to explore additional vertices rather than gather more information at a few, though never exceeding the \(\frac{R_0 + \lambda S_{\max}(p_D - 1)}{R_0 + \lambda}\) distance bound of Section V-A. As the budget tightens to \(0.6\,B_{\mathrm{CIPP}}\) (Fig. 4 (b)), this slack disappears: LIPPmatches C-IPP’s travel distance while attaining comparable posterior variance at lower energy. Under the tight \(0.35\,B_{\mathrm{CIPP}}\) budget (Fig. 4 (c)), LIPPcan afford neither heavy sampling nor every useful vertex, which constrains it to shorter paths than C-IPPwhile still reaching the same posterior variance. Note that we retain only graphs on which all methods achieve similar posterior variance reduction; consequently, the \(0.35\,B_{\mathrm{CIPP}}\) cases tend to be graphs where a short path already suffices, and the full-budget cases those where it does not. These results indicate that LIPP’s longer paths are not an inherent limitation but a consequence of an overly generous budget; under a suitably tight budget, it matches or improves upon C-IPPin both distance and energy.

Lastly, regarding property (iv), the primary limitation of LIPPis computational cost. As shown in Fig. 5 (a), both C-IPPand LIPPscale steeply with graph size, as expected for NP-hard mixed-integer programs, but LIPPgrows markedly faster than C-IPP. This overhead has the two sources identified in Section V-B: the \(S_{\max}^2\)-fold increase in decision variables and, dominantly, a weaker continuous relaxation that enlarges the branch-and-bound search space. Fig. 5 (b) disentangles the two against the \(S_{\max}^2 \times\) C-IPPreference line, which captures the variable-count factor alone: when the coupling is negligible (small \(\lambda\)), LIPPtracks this reference, whereas at non-trivial sample mass it rises well above it, showing that the relaxation—not the variable count—drives the excess cost. Fig. 5 (b) further reveals a non-monotonic trend: once samples become heavy enough that the energy budget sharply restricts the feasible set, the search space shrinks and solve time falls again. Although planning is performed offline, scaling LIPPto the finer discretizations needed for real-world deployment will require more efficient exact and approximate solvers.

7 Conclusion and Future Work↩︎

In this work, we introduced LIPP(load-aware IPP), motivated by a fundamental limitation of classical C-IPP: when information acquisition physically alters the robot through cumulative sampling load, the traversal cost of future edges depends on past sampling decisions, making the problem inherently order-dependent — a property C-IPP cannot model. We showed that LIPP strictly generalizes C-IPP, recovering it exactly as sample unit mass \(\lambda \to 0\), and derived a MIQP formulation that jointly optimizes routing, visitation order, and sampling allocation under an energy budget. We further established theoretical bounds on the path-length increase of LIPP relative to C-IPP and validated these properties across 2,000 diverse mission scenarios, demonstrating that LIPP achieves comparable posterior variance at significantly lower energy cost as sample mass grows.

LIPP opens several directions for future work. First, the weaker LP relaxations introduced by load-dependent coupling motivate faster solution methods—such as tailored heuristics, tighter relaxations, or approximation schemes—to scale to larger environments. Second, LIPP lays the groundwork for collaborative, heterogeneous multi-robot informative path planning with physical sampling. For example, a mission might pair agile scouts capable of rapid sampling with slower carrier robots that transport large accumulated loads at lower relative energy cost. Such settings amplify the importance of explicitly modeling load-dependent mobility, and constitute a promising direction for future research.

References↩︎