Add constrained accessibility wrapper

[soukhova]

This PR introduces an exported wrapper function constrained_accessibility() and three internal helpers implementing the family of accessibility measures introduced in the "Family of accessibility measures derived from spatial interaction principles" by Anastasia Soukhov, Rafael H. M. Pereira, Christopher D. Higgins and Antonio Páez (2025):

• total_constrained() - total-constrained balancing factor, can be interpreted as associated with Wilson (1971)'s 'unconstrained' model
• singly_constrained() - uses Wilson (1971)'s single constraint (demand or supply side)
• doubly_constrained() - uses Wilson (1971)'s doubly-constrained model using iterative proportional fitting

These functions follow the patterns of existing functions like spatial_availability() and support grouping, flexible decay functions, and optional detailed OD-level (flow) or aggregate (by origin or destination) outputs.

---

Mathematical Basis

Total-Constrained

kappa_tot_ij = (W_j * f(c_ij)) / sum_k(W_k * f(c_ik))
Accessibility_i = sum_j(kappa_tot_ij * Total Opportunities)

Sum of accessibility equals total opportunities (supply) in the region. Note, return_demand_side = TRUE (market potential variant) is also available.

Singly-Constrained

Supply-constrained (destination totals fixed) when return\_demand\_side = FALSE:
kappa_singly_ij = (Demand_i * f(c_ij)) / sum_k(Demand_k * f(c_kj))
Accessibility_ij = kappa_singly_ij * Supply_j

Demand-constrained (origin totals fixed) when return\\\_demand\\\_side = TRUE:
hatkappa_singly_ij = (Supply_j * f(c_ij)) / sum_k(Supply_k * f(c_ik))
Accessibility_ij = hatkappa_singly_ij * Demand_i

In either case, totals match either the demand at each origin or supply at each destination, depending on variant.

Doubly-Constrained

Accessibility_ij = A_i * B_j * O_i * D_j * f(c_ij)
A_i = 1 / sum_j(B_j * D_j * f(c_ij))
B_j = 1 / sum_i(A_i * O_i * f(c_ij))

Iterative proportional fitting updates (A_i) and (B_j) until convergence. This ensures that row sums equal (O_i) (demand) and column sums equal (D_j) (supply).

Note, only OD-level outputs are available (as aggregate outputs just match inputs).

---

Usage Examples

# Total-constrained
constrained_accessibility("total", travel_matrix, land_use_data,
  travel_cost = "travel_time",
  decay_function = decay_exponential(0.1),
  demand = NULL, supply = "jobs",
  return_demand_side = FALSE
)

# Singly-constrained
constrained_accessibility("singly", travel_matrix, land_use_data,
  travel_cost = "travel_time",
  decay_function = decay_exponential(0.1),
  demand = "population", supply = "jobs",
  return_demand_side = TRUE
)

# Doubly-constrained (toy dataset with matching totals)
constrained_accessibility("doubly", tm_small, lu_small,
  travel_cost = "travel_time",
  decay_function = decay_exponential(0.1),
  demand = "population", supply = "jobs",
  return_demand_side = NULL
)