Instability-Aware Digital Twin Control for Grid-Constrained Low-Carbon Smart Ports

1. Introduction

Ports are becoming tightly coupled cyber-physical energy systems. Vessel traffic, berth operations, cargo-handling equipment, shore-power connections, terminal distribution networks, and digital control platforms interact under uncertain and time-varying conditions. Port authorities must therefore pursue service efficiency, electrification, energy security, and emission reduction within the same operational framework rather than as separate planning objectives (Acciaro et al. 2014; Davarzani et al. 2016; Bjerkan and Seter 2019; Lam 2019; Sdoukopoulos et al. 2019; Zhang et al. 2024).

Shore-side electricity can reduce local ship emissions during berthing, especially when supplied by lower-carbon electricity or coordinated with renewable generation (Winkel et al. 2016; Gutierrez-Romero et al. 2019; Dai et al. 2020; Stolz et al. 2021). Its operational value depends on vessel compatibility, connection timing, berth allocation, equipment demand, and grid headroom (Wu and Wang 2020; Zhen et al. 2022; Zhen et al. 2024; Wang et al. 2024a). The same coordination problem applies to photovoltaic, offshore-wind, wave, tidal, and OTEC resources because their benefit depends on how generation, curtailment, storage, shore power, and terminal loads are synchronized (Systems 2024; Wu et al. 2024).

The operational literature includes berth allocation, quay-crane assignment, yard operations, and integrated terminal scheduling under uncertainty, energy, and emission objectives (He, Huang, Yan and Wang, 2015; He, Huang and Yan, 2015; He, 2016; Wang, Wang and Meng, 2018; Wang et al., 2020; Rodrigues and Agra, 2022; Li et al., 2023; Chargui et al., 2023; Yu et al., 2023; Guo et al., 2024; Wang, Hu and Zhen, 2024; Cai et al., 2024). Port-energy studies have likewise advanced smart-grid operation, flexible berth allocation, coordinated energy-logistics scheduling, and uncertainty-aware multi-energy dispatch (Lam 2021; Mao et al. 2022; Zhang et al. 2022; Mo et al. 2024; Gao et al. 2025; Zhou et al. 2025).

Digital-twin research has established synchronized virtual representations for monitoring, prediction, disruption management, and decision support (Negri et al. 2017; Tao et al. 2019; Fuller et al. 2020; Jones et al. 2020; Ivanov and Dolgui 2021). Port applications extend these ideas to resilience, sustainability assessment, and vessel-energy evaluation (Zhou et al. 2021; Wang et al. 2024b; Zhou et al. 2024; Zhu et al. 2026; Zhang et al. 2025). Recent architecture and adoption studies clarify how port authorities can structure digital-twin use cases (Gil-Pereira et al. 2026; Saragani et al. 2026), while current reviews show a shift from visualization toward operational decision support (Gazzaneo et al. 2025). More broadly, recent research across marine, offshore, building-energy, and renewable-energy systems has demonstrated the growing convergence of digital twins, intelligent monitoring, uncertainty-aware decision support, adaptive control, renewable generation, and energy-storage technologies for improving operational awareness, resilience, and low-carbon performance across complex infrastructure systems (Samaei and Riffat 2025; Riffat et al. 2025a; Riffat et al. 2025b; Riffat and Samaei 2025c; Samaei and Riffat 2026d; Riffat et al. 2026b; Riffat et al. 2026b; Samaei and Riffat 2026b). The remaining gap is not logistics-energy integration by itself. It is the use of a continuously updated, interpretable operational-stress indicator inside a closed-loop joint controller under simultaneous congestion and grid-limit stress.

At the equipment level, Gao, Chang and Chen (2023) used a digital-twin environment and Q-learning to reduce automated stacking-crane energy consumption, while Li et al. (2024) developed a digital-twin-driven proactive-reactive framework for uncertain multi-equipment scheduling in container terminals. These studies establish digital-twin-enabled operational and energy optimization at equipment and terminal-subsystem levels, but they do not combine vessel flow, berth and crane control, shore power, flexible demand, time-varying grid limits, renewable storage, and an online port-wide operational-stress constraint within one synchronized controller.

More recently, Li et al. (2026) developed a digital-twin-based dynamic co-scheduling framework that integrates multi-equipment coordination with AGV charging and battery-energy management in sea–rail intermodal automated container terminals. That study extends operational-energy scheduling to dynamic intermodal equipment coordination, but it does not jointly address vessel queues, berth and crane control, shore power, flexible terminal demand, time-varying grid limits, renewable storage, and a port-wide online operational-stress constraint. Related port-energy digital-twin research has also progressed at system level. Pang et al. (2023) proposed a digital-twin architecture for integrated port energy systems, while Li, Fan and Qi (2025) developed a dynamic digital-twin-based energy-router optimization framework for zero-carbon ports. Neugebauer, Heilig and Voß (2024) showed that port digital-twin implementations differ substantially in synchronization, decision authority, physical–virtual coupling, and implementation maturity. These studies strengthen the port-energy and digital-twin context, but they do not combine vessel-level queue and service dynamics, berth and crane allocation, shore-power preservation, flexible terminal demand, time-varying grid limits, renewable storage, and an online operational-stress constraint within one rolling-horizon controller.

The classifications in Table 1 reflect the specific model formulation, synchronization scope, and implemented decision authority reported in each cited study. The table positions the present contribution and is not a systematic evidence-ranking exercise.

An arrival surge increases berth occupancy, crane deployment, yard activity, and shore-power demand. When requested power approaches a contractual grid limit, the controller may need to redistribute service resources and energy-intensive tasks. Poor coordination can lengthen berth occupation and amplify the anchorage queue. A schedule can therefore be physically connected and still perform poorly under coupled congestion and grid-limit exposure. The relevant requirement is coordinated operational feasibility and grid-operating compliance, not a claim of formal power-system dynamic stability.

The study makes three contributions. First, it formulates short-term operation of an electrified container terminal as a synchronized feedback system in which service activity shapes electrical demand and grid-operating constraints shape feasible scheduling. Second, it introduces a policy-weighted IPI as an operational regime indicator and embeds it as a joint constraint and event trigger in rolling-horizon control. Third, it evaluates policy-level performance within an externally grounded synthetic benchmark whose principal digital-twin architecture, quay-crane scale, and shore-power demand envelope are cross-checked against published real-terminal evidence, and examines the incremental contribution of the IPI, ramp-rate regulation, and coupled shore-power/equipment scheduling through paired experiments, sensitivity analyses, ablation, and a renewable-storage extension.

2. Integrated Instability-Aware Digital-Twin Architecture

2.1. System Architecture and Synchronization

The implementation is a simulation-based digital-twin control prototype for virtual commissioning and closed-loop policy testing (Figure 1). It contains physical and scenario inputs, a synchronized hybrid state, discrete-event vessel-flow and terminal-operation engines, an energy-balance engine, IPI assessment, and a rolling-horizon controller. Operational and energy states are updated every 60 s. The controller re-optimizes every 15 min over a 6 h prediction and optimized decision horizon. All 24 quarter-hour decision intervals are represented in the chromosome, but only the first 15 min segment is implemented before state synchronization and re-optimization.

2.2. Coupled Operational Stress and Adaptive Control

Vessel traffic and electrical demand are represented as one feedback system. Arrival pulses raise queue length, berth occupancy, equipment activity, and shore-power demand. The resulting grid-loading pressure can restrict the timing of energized service resources and flexible tasks, which can in turn lengthen berth occupation and reinforce congestion. Renewable generation and battery storage provide additional flexibility by reducing net import, absorbing surplus energy, and supporting ramp-rate compliance.

The IPI combines normalized queue length, grid loading, and the modelled emission consequence of the same state. It is a transparent policy-weighted composite, not a decomposition of statistically independent latent factors. An upward crossing of the event threshold initiates additional optimization of berth priorities, vessel and crane sequences, shore-power timing, flexible loads, renewable curtailment, and battery dispatch. The mechanism is illustrated in Figure 2. The feedback concept is qualitatively related to cascading interactions in interdependent systems, although the present index remains a port-specific operational composite rather than a network-cascade model (Buldyrev et al. 2010; Gao et al. 2012).

3. Mathematical Formulation of Coupled Operational-Stress Dynamics and IPI-Constrained Optimization

The model combines discrete vessel events with continuously valued energy states. Operational efficiency, service preservation, grid-operating compliance, renewable-energy use, storage, modelled emissions, and IPI-based regime control are evaluated in one synchronized structure. The term instability is used only for the defined operational-stress regime; it does not denote formal power-system or nonlinear-system stability.

3.1. State-Space Representation of the Coupled Smart-Port System

The complete synchronized state is

Here, $Q_k$ is the anchorage queue, $B_k$ is berth occupancy, $U_k$ is aggregate equipment utilization, $P_{SP,k}$ is shore-power demand, $P_{RE,k}$ is renewable power delivered to the terminal bus, $SO{C}_k$ is battery state of charge, $P_{g,k}$ is scheduled grid import in MW, $C_{e,k}$ is the modelled grid-plus-queue emission rate, and $I_k$ is the IPI. Power and energy are kept distinct: grid energy is obtained by accumulating $P_{g,k}\mathrm{\Delta }{t}_h$ and is reported in MWh.

The hybrid state transition is

where ${𝐝}_k$ contains stochastic arrivals, service requirements, productivity variation, renewable forecast error, flexible-task demand, and changes in the admissible grid operating limit. The control vector is

The first four blocks represent berth, crane, shore-power, and vessel-sequencing decisions. Battery and curtailment variables are active only in the renewable-storage extension. Throughout the formulation, $t_k=k\mathrm{\Delta }{t}_h$; notation written as $P\left(t_k\right)$ and $P_k$ refers to the same synchronized state, while vessel arrival, service-start, and departure events retain continuous timestamps.

3.2. Vessel Arrival and Service Dynamics

The vessel-traffic subsystem is an integer-valued discrete-event queue. Its time-varying arrival intensity is

where $\mathrm{\Delta }\lambda (t)$ contains daily, weekly, and scenario-specific variation. The nominal intensity is

At synchronization step $k$,

Each vessel enters service with residual normalized service work

where $T_{s,i}$ is the nominal service requirement under the standard two-crane reference allocation. Let $y_{i,j,k}\in \{\mathrm{0,1}\}$ indicate that quay crane $j$ is actively assigned to vessel $i$ at step $k$. Then

where ${\xi }_{i,k}$ is the vessel-specific productivity factor and $a_{i,k}^E\in \{\mathrm{0,1}\}$ indicates whether the scheduled energized service interval is accepted. The same assignment variables $y_{i,j,k}$ also enter the quay-crane power equation in Section 3.3; service progress therefore cannot be produced by an electrically inactive crane assignment. Residual work evolves as

The number of vessels admitted to berth service and the number completing service are

Queue length and the number of vessels in service are updated separately:

A vessel therefore leaves the anchorage queue when it is admitted to berth service, while berth capacity is released only at completion. Nominal service requirements follow

The value 6.8 h is the mean of the unbounded parent lognormal distribution. Applying the stated clipping rule gives a theoretical post-clipping mean of approximately 6.801 h. The coefficient of variation is 0.20 and the practical bounds are 4–14 h. Realized service duration is

The three vessels occupying berths at t = 0 are initialized at the start of service with their complete sampled residual requirements; no pre-horizon productive service is credited.

The state-based waiting indicator is

where μk is the rolling one-hour admission rate. Five berths and a 6.8 h parent mean imply a berth-limited nominal service ceiling of 0.735 vessel/h, equivalent to 17.65 vessels/day. Eight quay cranes under the standard two-crane allocation imply a lower crane-limited ceiling of 0.588 vessel/h, equivalent to 14.12 vessels/day. The crane subsystem is therefore the active nominal bottleneck. The baseline arrival intensity of 0.500 vessel/h corresponds to a traffic intensity of 0.85 and an approximately 15% mean crane-capacity margin. Strategy-specific completed-cohort averages are outcome-dependent; the complete realization-level service-work identity is reported in Supplementary Table S10.

3.3. Terminal-Equipment and Flexible-Load Dynamics

Quay-crane power is calculated from the active assignment variables used in Equation (8) and a non-optimizable electrical duty factor $d_{QC,j,k}\in [\mathrm{0,1}]$ generated by the common crane-cycle profile:

Yard-equipment power is

where $a_{Y,j,k}$ is the active-unit indicator produced by the terminal-operations engine and $d_{Y,j,k}\in [\mathrm{0,1}]$ is its common duty profile. Yard activity is driven by the same accepted cargo-service events as quay-crane progress and is not an independent energy-only decision. The aggregate handling-equipment demand and utilization ratio are

Duty-profile innovations are sampled once per synchronized minute and are common across policies. For quay cranes, zQC, k follows Beta (2.5, 6.0) and dQC, k = clip [0.65dQC, k−1 + 0.35zQC, k, 0.10, 0.85], with dQC,0 = 0.2941. For yard units, zY, k follows Beta (2.8, 5.2) and dY, k = clip [0.75dY, k−1 + 0.25zY, k, 0.10, 0.80], with dY,0 = 0.3500. The recursion is evaluated before clipping, and the duty factors scale electrical power without independently creating cargo progress.

The uncontrolled flexible-demand request is

where ${\epsilon }_f(t)\sim 𝒩(0,{0.20}^2)$ MW. Aggregate scheduled flexible demand is linked to task-level power and standby states by

For every task $𝓁$ with admissible window ${𝒲}_{𝓁}$,

The piecewise conversion efficiency is

Standby consumption contributes to grid-side energy but not to useful service. The 420 tasks comprise 240 reefer-service tasks released every 3 h with 24 h windows, 120 vehicle/equipment-charging tasks released four times per day with 8 h windows, and 60 workshop tasks released twice per day with 6 h windows. Their common useful-energy requirement is 2,214 MWh over 720 h: 1,512 MWh for reefer service, 486 MWh for charging, and 216 MWh for workshop and controllable services. Power is zero outside each admissible window; released and unfinished tasks draw standby power only when not active. Reefer and charging tasks may be pre-empted in 15 min blocks, whereas workshop tasks are non-pre-emptive after start with a 60 min minimum run. Every task must meet its useful-energy requirement before its deadline. Complete release, power, efficiency, standby, and interruption rules are reported in Supplementary Table S12.

3.4. Shore-Power Allocation Model

Shore-power allocation is constrained by vessel compatibility, berth presence, shore-point availability, the vessel-specific hoteling requirement, and terminal electrical capacity. Total demand is

where $z_i\in \{\mathrm{0,1}\}$ is the compatibility flag, ${\delta }_{i,k}\in \{\mathrm{0,1}\}$ is the connected state, and $P_i^{hotel}\sim U(\mathrm{2,8})$ MW. The common compatibility fraction is 45.3%. The connection state is

subject to

Setup and disconnection occupancy is inserted deterministically by the shore-power decoder. Ten consecutive setup steps precede the first energized step, and five consecutive disconnection steps follow the last energized step. Both transition states occupy a shore-power point and deliver zero hoteling energy. The three initially occupied vessels begin unconnected and require the complete setup interval. An incomplete end-of-horizon transition is retained in the no-arrival continuation, and disconnection must finish before berth release and recorded departure.

For each compatible vessel, ${𝒯}_i^{SP}$ denotes the required energized hoteling window specified by its common manifest demand and constrained to lie within its realized berth-occupancy set ${𝒯}_i^B$:

Required and delivered hoteling energy are

with the service-preservation condition

The vessel-specific energized-window duration is generated as H_iSP = clip (0.70T, _i + ν_i, 2, 10) h, where ν_i is normally distributed with zero mean and 0.75 h standard deviation, and the resulting window is constrained to the realized berth stay. Connection setup and disconnection occupy one shore-power point for 10 and 5 min, respectively; no hoteling energy is credited during either transition. At every minute, the connected, setup, and disconnection states together cannot exceed five points, while connected hoteling demand cannot exceed 35 MW. The required energized window is not the vessel’s complete cargo-service duration. The same compatibility, hoteling power, and required-energy attributes are supplied to all strategies, and unmet shore-power energy is reported explicitly.

3.5. Renewable-Generation and Battery-Storage Model

The renewable-energy layer represents photovoltaic generation, offshore wind, wave-energy conversion, tidal-stream generation, and ocean thermal energy conversion. Total available renewable power is

The deterministic 168 h base profiles are generated from the technology-specific equations and scale factors reported in Supplementary Table S37, repeated four times, and followed by the first 48 h of a fifth template. Hourly base values are held constant within the 60 synchronized one-minute steps. Forecast errors follow technology-specific AR(1) recursions with contemporaneously correlated Gaussian-copula innovations, common across policies within each replication.

The correlation matrix defines contemporaneous correlation among the pre-clipping Gaussian innovations. AR filtering and technology-specific clipping can alter the realized correlations of the bounded power-error trajectories; the input matrix is therefore not interpreted as the achieved output-correlation matrix.

For each technology $r$,

Renewable curtailment and delivered renewable power are

Renewable power accepted by the terminal is allocated explicitly between direct supply to terminal loads and battery charging. Battery charging is restricted to accepted renewable energy; grid-to-battery charging and power export are disabled. Renewable availability that is neither supplied directly to the terminal nor accepted by the battery is recorded as curtailment. These source-allocation variables prevent implicit grid charging or unreported disposal of renewable surplus.

The battery state is updated as

with

Simultaneous charging and discharging are prohibited:

A terminal condition prevents the initial battery inventory from becoming unaccounted energy:

The verified initial and terminal SOC are both 60%, so the net change in stored energy is zero. The maximum replication-level terminal-SOC residual is 2.45 × 10^-4, and no replication exceeds the 10^-3 tolerance.

3.6. Operational Power Balance and Grid-Operating Constraints

The gross operational power requirement is

Auxiliary demand covers lighting, communication, administration, refrigeration support, and other non-handling services; it is centred at 5 MW and bounded within 4.75–5.25 MW. Because export is outside the benchmark scope, scheduled grid import is

The admissible grid operating limit is time dependent:

with

The Scenario 4 limits are 54, 48, 45, and 51 MW. These values are contractual and operational scheduling limits. For the controlled synthetic benchmark, the upstream physical import rating is set to

which is equal to the maximum connected-load envelope. Scheduled grid import is physically bounded by

Exceedance of the time-dependent operating limit is measured as

The proposed controller enforces v₍g,k₎ = 0, whereas unconstrained benchmarks may produce v₍g,k₎ > 0, which is retained as an operating-limit violation. Because every scheduled peak remains below the 66.85 MW physical rating, the experiments contain no involuntary load shedding or unmet electrical energy.

Minute-level ramp feasibility is enforced by

The connected-load envelope is $14.4+7.2+35+5.25+5.0=66.85$ MW. Subsystem energy is

and gross operational energy is

Net grid-import energy, mean demand, and observed peak are

For time-varying limits, average utilization is evaluated from minute-level ratios:

Load factor and peak-to-average ratio are

Gross operational energy and net grid-import energy are identical in grid-only experiments. In the renewable-storage extension, gross operational energy remains the subsystem requirement, while renewable generation and battery operation reduce net grid import.

3.7. Carbon-Emission Formulation

The reported emission boundary comprises electricity-related grid emissions and anchorage-queue congestion emissions. The instantaneous modelled rate is

where

The factor 1000 converts MW to kW. Cumulative modelled grid-plus-queue emissions are

Operational renewable generation is assigned zero point-of-use emissions. The boundary excludes auxiliary-engine emissions from vessels already at berth, emissions during shore-power connection transitions, and lifecycle emissions of energy infrastructure. The reported values therefore represent modelled grid-plus-queue operational emissions, not total port emissions. The coefficients ${\alpha }_g=0.52$ kgCO₂e/kWh and ${\beta }_q=180$ kgCO₂e/(vessel h) are scenario-specific comparative parameters rather than site-calibrated inventory factors; absolute emission values require recalibration for a deployment port.

3.8. Policy-Weighted Operational-Stress Index

In this study, instability denotes an operational-stress regime produced by simultaneous queue growth, high grid loading, and their modelled emission consequence. It is not a small-signal, transient, voltage, frequency, or Lyapunov stability measure. The IPI is used only for operational regime classification and control triggering within the defined benchmark.

The principal modelling contribution of the framework is the explicit representation of coupled congestion and electrical stress with emission consequences. In a highly utilized smart port, an arrival disturbance may increase berth occupancy, crane deployment, yard-equipment activity, and shore-power demand. The resulting grid stress may restrict service resources, extend berth occupation, and further increase the queue. The Instability Propagation Index quantifies the policy-weighted intensity of these mutually reinforcing mechanisms and their modelled emission consequences. The index is defined as

where

where Q_crit = 18 vessels, P_g,lim(t) is the time-dependent admissible grid capacity, and C_target is the fixed emission-normalization target. The normalized weights satisfy

Under the baseline configuration,

The emission-normalization target is calculated from the nominal 60 MW grid limit and the critical queue threshold:

Equal nominal weights do not imply equal marginal influence of queue length and grid loading because the emission-consequence term reweights both variables. Under the equal-weight configuration, the reduced-form effective coefficient of normalized queue length is approximately 0.365. The corresponding effective coefficient of normalized grid loading varies from approximately 0.560 to 0.635 as the admissible Scenario-4 grid limit varies from 45 to 60 MW. The equal-weight specification is therefore interpreted as equal weighting of the three policy terms, not equal sensitivity to the two underlying state variables.

A fixed C_target is retained across all scenarios to preserve comparability. Substitution of Equation (62) into Equation (66) shows that the emission-consequence term is algebraically dependent on queue length and grid import. It is retained for policy interpretability and carbon-consequence weighting, not as a third statistically independent state dimension. Because the emission rate is calculated from grid import and anchorage queue length, the three terms are not statistically independent latent factors. The IPI is therefore interpreted as a transparent policy-weighted operational composite: the first two terms describe congestion and electrical loading, while the third assigns the environmental consequence of the same coupled state. The thresholds are operational decision boundaries rather than universal physical constants and are tested through weight and threshold robustness analyses. The operating regime is classified as

The short-term instability-growth rate is evaluated discretely as

A positive value of $g_I(k)$ indicates instability amplification, whereas a negative value indicates recovery toward a lower-stress operating condition. This explicit expression replaces an unspecified nonlinear propagation operator and allows the instability trend to be reproduced directly from the synchronized state trajectory. Unless otherwise stated, the equal-weight configuration is used in the principal simulations. Weight robustness is evaluated by varying each coefficient between 0.20 and 0.50 while preserving the weight-normalization condition in Equation (68). In the term Instability Propagation Index, propagation refers only to temporal amplification of the defined composite across synchronized steps; it does not denote identification of a causal network-cascade process or a formal propagation operator.

3.9. IPI-Constrained Multi-Objective Optimization Formulation

The control problem minimizes composite vessel delay, gross operational energy, and cumulative modelled emissions while enforcing operational, electrical, storage, service-preservation, and instability constraints:

For each vessel arriving during the 720-h evaluation horizon, anchorage waiting and observed service exposure are right-censored at T. The operational metric is defined by

For a vessel not admitted by T, t_i^start is treated as infinite and anchorage exposure ends at T. For a vessel still in service, t_i^dep is treated as infinite and observed service exposure ends at T; no unobserved post-horizon duration is imputed. The two vessels initially present in the anchorage queue are included in state evolution, queue-emission accounting, conservation, and terminal-horizon checks but are excluded from the horizon-arrival delay numerator and denominator. The three initially occupied vessels are treated analogously for state and workload conservation. This convention is applied identically to all policies.

The remaining objectives are

The optimization is subject to berth-capacity and non-overlap constraints,

crane availability,

shore-power limits and vessel-level service preservation,

grid capacity and ramp-rate limits,

battery feasibility,

renewable curtailment,

and avoidance of the instability-dominated regime,

Queue non-negativity is imposed as

while $Q_{crit}=18$ vessels is used as an instability-normalization and early-warning threshold rather than a physical queue-capacity limit. Vessel sequencing satisfies

The formulation also enforces flexible-task useful-energy requirements, shore-power coverage, terminal SOC, vessel-berth compatibility, service precedence, and deterministic repair of infeasible mixed-variable chromosomes.

For performance reporting, berth-utilization efficiency is

where setup, idle connection, and post-service occupancy are excluded from productive berth time.

Each non-dominated solution $m$ is normalized as

with distance from the utopia point

The implemented policy is

where P is the NSGA-II non-dominated set. The implemented policy is the feasible knee solution nearest the normalized utopia point and must satisfy the numerical IPI margin at every synchronized step.

4. Computational Framework, Scenario Configuration, and Implementation

The framework was implemented as a high-resolution co-simulation prototype that couple’s stochastic vessel flow, berth and quay-crane scheduling, yard-equipment activity, shore-power allocation, flexible-load management, renewable generation, battery dispatch, grid-operating-limit checks, IPI assessment, and multi-objective optimization. It is a virtual commissioning environment rather than a deployed operational twin: no live terminal asset, online state-estimation interface, or site-calibrated bidirectional data connection is used.

4.1. Digital Twin Computational Environment

The implementation combines five synchronized engines: vessel-flow simulation, terminal operations, energy balance, IPI assessment, and rolling-horizon control. Physical and energy states are updated every 60 s, producing 43,200 state records in each 720-h replication. The control problem is reconsidered every 15 min over a 6 h prediction and optimized decision horizon, represented by 24 quarter-hour control intervals. Only the first 15 min segment is implemented before the state is updated and the rolling horizon advances.

The first interval uses a 120-individual, 300-generation solve. Scheduled updates are warm-started from the preceding non-dominated population and stop under the convergence rule in Section 4.4. Each replication contains one initial solve and 2,879 scheduled refinements averaging 11.6 generations. Additional event-driven refinement occurs only on an upward crossing of IPI = 0.50 and after a 15 min refractory interval. The proposed controller averages 14.2 event-triggered refinements with 38.4 generations per triggered solve. This additional computation is part of the proposed method and is reported separately rather than treated as an equal-budget property of the benchmark comparison.

Random seeds 202601–202630 define the paired common-random-number experiment. Within each replication, all strategies receive the same vessel arrivals, service requirements, shore-power compatibility and hoteling demand, crane-productivity factors, flexible tasks, disturbance pulses, renewable errors, and grid-limit trajectory. The optimizer is initialized from a deterministic substream associated with the same replication seed and applied consistently within each paired comparison. The primary 30-replication variance therefore reflects the combined stochastic experiment under this fixed seed mapping. A separate 10 stochastic-seed × 5 optimizer-substream audit quantifies optimizer-initialization variance and is reported in Section 4.6 and Supplementary Table S36.

The optimization is non-anticipative. At decision time $k$, the controller has access only to synchronized states observed up to $k$, manifest attributes of vessels already available to the controller, and forecasts over the common 6 h prediction horizon. Future realized arrivals outside the horizon, future productivity realizations, and future renewable forecast errors are not exposed. Supplementary Figure S1 provides the complete workflow.

4.2. Operational Definition and Benchmark Configuration of the Smart Port

For the benchmark, a smart port is represented as a container-terminal cyber-physical system with modelled observability of vessel, berth, equipment, energy, and emission states at the 60 s synchronization interval; controllable electrified handling, shore-power, renewable, and storage assets; explicit monitoring of grid capacity, net-load ramping, emissions, and operational stress; and closed-loop authority over berth, crane, vessel, and energy decisions.

The benchmark terminal has five deep-water berths, five shore-power points, eight electrified quay cranes, twelve electrified yard-handling units, and a 60 MW nominal grid operating limit (Table 2). Each occupied vessel nominally receives two quay cranes when equipment is available. Based on the 6.8 h parent mean, the five berths provide a nominal ceiling of 0.735 vessel/h, whereas the eight-crane fleet under the standard two-crane allocation provides a lower ceiling of 0.588 vessel/h. The quay-crane fleet is therefore the active nominal service bottleneck. Shore-power demand is sampled from 2–8 MW for a common compatible-vessel fraction of 45.3%. The initial state contains two queued vessels and three occupied berths.

The benchmark represents a controlled, reproducible smart-port environment rather than a calibrated replica of one commercial terminal. Site deployment would require recalibration using terminal-operating-system records, AIS histories, vessel calls, equipment telemetry, smart-meter data, renewable-resource measurements, and grid-connection agreements. The present implementation is therefore described as a simulation-based digital-twin control prototype; the reported results establish computational performance only within the defined synthetic benchmark.

Selected principal operating scales are externally cross-checked against published real-terminal evidence in Section 4.3 and Table 3. The full parameter set and input-generation rules are reported in Supplementary Tables S1 and S2. The connected-load envelope is 66.85 MW, exceeding the nominal 60 MW contractual limit. Values above the applicable time-dependent operating limit are recorded as violations; every policy remains below the physical envelope. Although the benchmark is synthetic in realization, its principal operational and electrical scales were not selected without external reference. To assess external plausibility, the vessel–terminal digital-twin structure, quay-crane configuration, shore-power demand envelope, and energy-management context were cross-checked against independently published studies using operational data from real container terminals, as described in Section 4.3.

4.3. External Empirical Validation of the Benchmark Operating Envelope

The synthetic benchmark was subjected to an external empirical validation of its principal operating envelope using independently published evidence from real container terminals. This assessment examines whether the major scales governing the coupled logistics–energy problem—digital-twin synchronization, vessel–berth coordination, quay-crane configuration, shore-power demand, and terminal-energy management—are consistent with documented operating conditions. The objective is not to reproduce or retrospectively calibrate a single commercial terminal, but to determine whether the benchmark occupies an empirically defensible region of real smart-port operation. At system level, Eom et al. (2023) developed a port digital twin using operational data from the Pusan Newport International Terminal at Busan New Port and evaluated the proposed digital-twin scheduling approach against actual terminal operations recorded in September 2022. Their framework explicitly links vessel information, terminal operations, scheduling, and carbon-emission consequences and uses real-time data to support recurrent operational decision making. The reported comparison with actual operation provides an empirical precedent for the central modelling premise adopted here: that synchronized vessel–terminal information and recurrent digital-twin-based scheduling can influence waiting, operational coordination, and emission outcomes in a real container-port environment. Independent support for the digital-twin energy layer is provided by Zhang et al. (2025), who constructed and evaluated a digital-twin energy-management platform at an operating container terminal in the Wuhan Port Area. That study integrated multi-source operational information with electric quay cranes, electric rail-mounted gantry cranes, automated guided vehicles, and vessel-service energy assessment at minute-level temporal resolution. This evidence supports the practical relevance of coupling synchronized terminal operations with energy monitoring and management, which forms the second major layer of the present framework. The quay-crane fleet scale also has direct empirical support. Geerlings, Heij and van Duin (2018) developed an energy-demand model for an operational container-terminal context with six to eight ship-to-shore cranes and used real terminal information to assess crane electricity demand and peak-shaving opportunities. Their results show that an eight-crane operating configuration represents a realistic intermediate container-terminal scale and that crane operating rules materially affect electrical peak demand. The eight electrified quay cranes used in the present benchmark are therefore consistent with an independently studied real-terminal operating configuration rather than being an arbitrary equipment count. The vessel-specific shore-power range provides the clearest direct parameter-level external validation. In an empirical cold-ironing assessment based on container-vessel and Port of Genoa information, large container vessels exceeding 140 m were reported to have an average at-berth electrical demand of approximately 2 MW and peak demand of up to 8 MW. The 2–8 MW shore-power envelope adopted in the present benchmark therefore directly spans a published real-world operating range for large container vessels. This agreement is particularly important because shore-power demand is one of the primary drivers of the coupled berth–grid interaction examined in the proposed controller.

Taken together, the external evidence indicates that the benchmark is synthetic in stochastic realization but empirically grounded in its principal physical and operational scales. The shore-power envelope directly matches published container-vessel operating data; the eight-crane equipment scale is consistent with independently studied terminal configurations; and the synchronized operation–energy architecture is supported by digital-twin implementations evaluated in operating container terminals. The benchmark should therefore be interpreted as a controlled, reproducible integration of empirically supported operating conditions rather than an unconstrained hypothetical test case. The remaining distinction is deployment-specific calibration: application to a particular port would require local arrival histories, service distributions, grid agreements, emission factors, equipment telemetry, and operating rules.

4.4. Operational Scenarios and Disturbance Conditions

The 168-h input template is repeated over the 720-h horizon. Scenario 1 is a nominal finite-horizon reference with the baseline nonhomogeneous Poisson process and a 60 MW limit; its nominal arrival intensity of 0.500 vessel/h remains below the 0.588 vessel/h crane-limited mean-rate ceiling. Scenario 2 introduces arrival pulses of 18 vessels/day and peaks of 20 vessels/day. Scenario 3 retains nominal arrivals and reduces the grid operating limit to 48 MW. Scenario 4 combines four 20-vessel/day arrival pulses with temporary limits of 54, 48, 45, and 51 MW. Exact intervals are given in Supplementary Table S3.

The 18-test confirmatory comparison is restricted to Scenario 4, the coupled congestion-and-grid-limit disturbance benchmark. Scenarios 1–3 define boundary and stress-reference inputs and are not pooled into the confirmatory statistical family. This separation avoids treating heterogeneous scenario definitions as repeated observations of one estimand.

4.5. IPI-Constrained Optimization Strategy

NSGA-II uses binary tournament selection with tournament size two, simulated binary crossover with probability 0.90 and distribution index 15, and polynomial mutation with distribution index 20. The grid-only joint controller has 96 active variables: 12 vessel-sequence keys, 12 berth assignments, and 24 control values each for crane allocation, shore-power priority, and flexible demand. The renewable-storage controller adds 24 battery-power and 24 renewable-curtailment variables, giving 144 active variables. The per-variable mutation probabilities are therefore $1/96=0.0104167$ and $1/144=0.0069444$.

Mixed decision blocks are repaired before objective evaluation. Berth incompatibility and duplicate assignments are corrected deterministically; crane totals are restricted to eight; shore-power assignments are restricted to five points and 35 MW; flexible tasks are repaired to satisfy useful-energy requirements; continuous energy variables are clipped to physical bounds; and total normalized constraint violation is used in Deb’s feasibility comparison.

Population-generation combinations of 80/200, 120/300, and 160/400 were tested under Scenario 4. Using fixed normalization bounds of 6–20 h for composite delay, 14,000–19,000 MWh for energy, and 7,500–10,500 tCO₂e for modelled emissions, with the reference vector (1.10, 1.10, 1.10), the selected 120/300 configuration achieved a normalized hypervolume of 0.6414, 4.8% above 80/200. Increasing the search to 160/400 raised hypervolume only to 0.6459 and reduced mean maximum IPI by 0.01 while runtime increased from 17.4 to 26.9 min. The 120/300 choice was retained in a held-out check using seeds 202631–202640. Convergence is monitored using

where εHV = 10^-6 is a denominator safeguard used only in the relative hypervolume-change calculation. The convergence criterion must hold for 50 consecutive generations, or until the generation limit is reached. The normalization bounds and reference vector were fixed before comparing the three tuning configurations and were not recalculated by configuration. Infeasible vectors are ranked through Deb’s total normalized-violation rule. Complete encoding, tuning, convergence, and execution details are reported in Supplementary Section S4 and Tables S4–S6B and S24; optimizer-initialization variance components are reported separately in Supplementary Table S36.

4.6. Benchmark Models and Comparison Scope

Three benchmarks represent distinct operating policies (Table 4). The operation-focused policy minimizes composite vessel delay and berth idle time and uses earliest-feasible energy dispatch without active grid coordination. The energy-focused policy minimizes energy and modelled emissions using grid and shore-power repair and penalization, with vessel-level service preservation taking precedence; any residual operating-limit or ramp-rate violation after repair is retained and reported. The decoupled policy fixes the operational schedule before energy dispatch and retains residual violations when stage-two dispatch cannot repair the schedule. The proposed policy jointly minimizes composite vessel delay, gross energy, and modelled emissions while enforcing service preservation, time-dependent grid limits, ramping, and the IPI condition.

All strategies receive the same observations, forecasts, 6 h prediction horizon, 15 min scheduled update interval, and stochastic realizations. The scheduled NSGA-II budget is 120 × 300 for the operation-focused, energy-focused, and proposed policies; the decoupled policy uses 120 × 150 in each of two stages. The proposed policy additionally uses the event-triggered refinements quantified in Section 4.1. The four-policy results are therefore policy-level engineering comparisons rather than a pure equal-computation algorithm contest. A delay-only diagnostic using the complete proposed operational decision space and the same effective evaluation budget achieved a mean composite vessel delay of 7.4 h. Because NSGA-II is heuristic, this value is an optimization-quality reference rather than a certified lower bound or guaranteed global optimum. The without-IPI-constraint configuration uses the same nominal decision space as the proposed controller after removal of the hard IPI constraint. It is treated as a controller-configuration diagnostic rather than an exact causal isolation of the constraint because NSGA-II is heuristic and the closed-loop state and optimization trajectories diverge after the first differing decision.

Supplementary Table S6A gives decision blocks and information controls, while Table S6B reports the proposed controller’s additional event-triggered computation.

4.7. Numerical and Statistical Analysis

Capacity, ramp-rate, service-preservation, storage, and IPI conditions are checked at the 60 s synchronization resolution. The proposed implementation requires 17.4 min on average for one 720 h replication. Computation used Python 3.12.4, NumPy 2.0.1, SciPy 1.14.0, pymoo 0.6.1.3, NumPy PCG64 streams, and OpenBLAS 0.3.27 on Ubuntu 22.04.4 LTS with an AMD Ryzen 9 7950X processor, 64 GB RAM, and 30 replication-level workers. Runtime comprises 8.6 s for the initial solve, 964.5 s for scheduled warm-start refinements, 14.9 s for event-triggered refinements, and 54.9 s for state updating and bookkeeping.

For each primary metric, the mean and two-sided 95% confidence interval are

The confirmatory family contains 18 paired comparisons: six primary metrics against three benchmarks. Shapiro-Wilk screening is supplemented by a two-sided paired sign-flip permutation test with 100,000 permutations and BCa paired bootstrap intervals based on 10,000 resamples. Holm correction is applied across the same 18-test family, and Cohen’s dz is reported. All inferential routes preserve the direction and significance of the primary comparisons. Across the 30 independent Scenario 4 seeds, total arrivals have a mean of 454.6 vessels, a sample standard deviation of 21.20, and a t29 95% confidence interval of 446.684–462.516 vessels. In the 10 stochastic-seed × 5 optimizer-substream audit, optimizer variance is estimated from the mean within-seed variance, while simulation variance is the between-seed variance of substream means after subtracting the within-seed contribution divided by five. Optimizer initialization accounts for at most 2.211% of total variance for every primary metric, and policy rankings are unchanged. Relative confidence-interval half-widths remain below the prespecified 7.5% tolerance.

5. Results and Discussion

The primary comparison is the 720-h Scenario 4 coupled congestion-and-grid-limit disturbance experiment. The benchmark used for these comparisons has been externally cross-checked against published operating-terminal evidence as described in Section 4.3; the reported effects nevertheless remain estimates for the controlled benchmark rather than direct measurements from a deployment port. All means, confidence intervals, paired tests, energy quantities, emissions, grid-limit indicators, and IPI summaries are evaluated from the underlying 60 s records. Figures 3 and 4 visualize only quantities that can be reconstructed directly from the reported replication-level summaries; they are not synthetic time histories.

5.1. Comparative Performance under the Coupled Congestion and Grid-Limit Disturbance

Figure 3 expresses each policy’s change relative to the operation-focused benchmark. Positive values denote improvement: reductions in composite vessel delay, energy, emissions, peak grid import, and maximum IPI, and an increase in berth-utilization efficiency. The proposed policy is the only strategy that improves all six reported metrics relative to the operation-focused baseline. The energy-focused policy reduces energy-related indicators but increases composite delay and lowers berth utilization; the decoupled policy provides intermediate performance.

5.2. Operational-Stress Regime and Grid-Limit Exposure

Figure 4 relates the mean replication-level maximum IPI to mean time above the applicable time-dependent grid operating limit. Operation-focused, energy-focused, and decoupled scheduling have mean maximum IPI values of 1.42, 1.16, and 1.09. The proposed policy has a mean maximum of 0.82 and 0.00 h above the operating limit. Across its 30 replications, the empirical 95th percentile and largest observed minute-level IPI are 0.97 and 0.99; no replication reaches IPI = 1.00. The largest uncontrolled pre-constraint import request is 63.7 MW. The scheduled proposed-policy import remains below the applicable time-dependent operating limit at every synchronized step, and no physical-limit violation occurs. Marker size represents the final queue and shows that lower electrical exposure is accompanied by lower terminal congestion.

The IPI values in Table 5 are means of replication-level maxima, not maxima of pointwise mean trajectories. The implemented feasibility condition is I_k ≤ 1 − 10^-6. Supplementary Table S26 reports empirical percentiles, largest observed values, exceedance counts, recovery times, and physical-limit checks for all policies.

5.3. Comparative Operational, Energy, Modelled-Emission, and IPI Performance

Relative to operation-focused, energy-focused, and decoupled scheduling, the proposed framework reduced energy by 21.2%, 8.9%, and 14.0%; modelled grid-plus-queue emissions by 21.0%, 8.5%, and 13.6%; peak loading by 15.1%, 5.7%, and 9.1%; and mean maximum IPI by 42.3%, 29.3%, and 24.8%. Mean composite vessel delay fell by 9.2%, 30.7%, and 19.4%, while berth-utilization efficiency increased by 3.2, 10.9, and 6.4 percentage points.

The complete 18-test family remains significant after Holm correction, and the permutation and BCa-bootstrap checks give the same conclusions. The policy tests do not by themselves isolate the IPI from objective and decision-authority differences. In the paired same-space ablation, removing the IPI constraint increases composite delay by 0.60 h (95% CI 0.32–0.88; dz = 0.80), energy by 470 MWh (260–680; dz = 0.84), emissions by 270 tCO₂e (145–395; dz = 0.81), peak import by 3.9 MW (2.8–5.0; dz = 1.32), and maximum IPI by 0.26 (0.20–0.32; dz = 1.62); all Holm-adjusted p-values are below 0.001.

Vessel population and nominal service work are conserved across admission, service, completion, residual work, and the terminal queue. Within each replication, the common cohort is formed from the actual vessel IDs completed under all four policies. The mean common-cohort size is 392.0 vessels (SD 5.206; 95% CI 390.056–393.944), the mean compatible-vessel count is 177.6 (SD 13.124; 95% CI 172.699–182.501), the mean nominal service-work sum is 2,665.6 work-h (95% CI 2,648.488–2,682.712), and the mean required and delivered hoteling energy is 4,420.0 MWh (95% CI 4,276.643–4,563.357). Population expectations are not substituted for observed cohort attributes. A no-arrival continuation completes the four systems at 828.4, 781.6, 760.9, and 744.8 h, respectively. Detailed balances are reported in Supplementary Tables S10, S11, and S27.

Compatible-vessel shore-power comparisons use the required energized hoteling window defined in Section 3.4, not the entire cargo-service duration. The proposed policy admits 199.8 compatible vessels on average compared with 183.0 under operation-focused scheduling, while delivering every admitted vessel’s required shore-power energy within the 0.1% tolerance. The reported connection-hours therefore reflect strategy-specific admitted vessel sets and hoteling requirements rather than deletion of required service.

5.4. Trade-Off Analysis Using the Composite-Delay-Plus-Nominal-Service Indicator

The horizontal metric combines mean composite vessel delay with the 6.8 h nominal service requirement. It is a policy-level turnaround indicator and is not an event-level vessel turnaround measurement.

Figure 5 uses the four 30-replication policy means reported in Table 5. It replaces the seed-specific Pareto visualization so that every displayed coordinate is directly traceable to the confirmatory summary dataset.

The proposed controller occupies the lowest-delay and lowest-emission position among the four policy means. The figure is used only to visualize the policy-level trade-off; multi-objective selection and feasibility remain defined by Equations (92–94).

5.5. Operational-Energy Requirement, Workload Fairness, and Service Preservation

Over 720 h, operation-focused, energy-focused, decoupled, and proposed scheduling consume 18,420, 15,930, 16,870, and 14,510 MWh. The largest aggregate difference relative to operation-focused scheduling occurs in shore-power energy, followed by quay-crane, yard-equipment, and flexible-load energy. All strategies receive the same 454.6 mean realized arrivals, two initially queued vessels, and three vessels initially in service.

Population accounting gives mean completions of 399.7, 417.3, 427.5, and 436.7 vessels and final queues of 56.1, 38.2, 28.1, and 18.9 vessels for operation-focused, energy-focused, decoupled, and proposed scheduling, respectively. The expected workload obtained from the 6.8 h parent mean is 3,125.28 work-h; conservation is evaluated instead from the directly sampled vessel requirements. Across the 30 replications, the mean realized total service workload is 3,126.942 work-h. Completed work is 2,877.8, 2,858.5, 2,855.7, and 2,825.4 work-h; residual in-service work is 7.6, 8.2, 8.0, and 36.0 work-h; and queued work is 241.542, 260.242, 263.242, and 265.542 work-h. The maximum absolute replication-level work-conservation residual is $3.06\times {10}^{-10}$ work-h. Energy per completed productive-service hour is 6.40, 5.57, 5.91, and 5.14 MWh/work-h, respectively.

The mean within-replication common-cohort size is 392.0 vessels, with 2,665.6 nominal work-h, 177.6 compatible vessels, and 4,420.0 MWh of required and delivered hoteling energy. Common-cohort crane-plus-yard energy is 4,370, 3,890, 4,130, and 3,620 MWh for operation-focused, energy-focused, decoupled, and proposed scheduling. The corresponding no-arrival drain-horizon completion times are 828.4, 781.6, 760.9, and 744.8 h. These checks confirm the policy ordering without assigning system-level auxiliary demand or queue emissions to individual vessels.

Required useful flexible service is 2,214 MWh under every policy. Mean compatible vessels entering berth service are 183.0, 191.0, 195.6, and 199.8, with zero unmet shore-power energy and complete coverage of the required energized window for every admitted compatible vessel. Detailed checks are reported in Supplementary Tables S9–S14.

5.6. High-Load Exposure and Grid-Operating-Limit Analysis

Grid stress was evaluated from all 43,200 one-minute observations in each replication. The operation-focused, energy-focused, decoupled, and proposed strategies remained above 54 MW for 64.60, 25.38, 38.92, and 5.29 h, respectively. Fixed 60 MW exceedance lasted 27.47 h under operation-focused scheduling and 2.35 h under decoupled optimization, while the other strategies recorded none.

Feasibility was also assessed against the actual time-dependent operating limit rather than only the fixed planning levels. Mean hours above the applicable limit are 52.84, 18.62, 26.47, and 0.00 h for operation-focused, energy-focused, decoupled, and proposed scheduling. Maximum operating-limit-utilization ratios are 1.438, 1.184, 1.286, and 0.986, and cumulative violation magnitudes are 218.4, 39.1, 74.2, and 0.0 MW·h. The corresponding ramp-rate violation counts are 81, 24, 43, and 0, with maximum ramps of 9.8, 8.4, 9.1, and 6.7 MW/min.

5.7. Renewable-Energy and Battery-Storage Extension

The renewable-storage extension uses 6 MW PV, 10 MW offshore wind, 3 MW wave, 2 MW tidal, and 2 MW OTEC, together with an 8 MW/16 MWh battery. Over 720 h, gross generation is 6,508.8 MWh. After 240.0 MWh of curtailment, 6,268.8 MWh is accepted: 4,947.8 MWh serves demand directly and 1,321.0 MWh charges the battery. Storage subsequently delivers 1,192.2 MWh, producing a 128.8 MWh round-trip loss and 6,140.0 MWh of total renewable service. Battery throughput is 2,513.3 MWh, equivalent to 78.54 full cycles over 30 d, or approximately 2.62 cycles/day. Observed SOC ranges from 21.3% to 88.7%, initial and terminal SOC are 60%, the mean absolute terminal energy residual is 0.0003 MWh, and the simultaneous charge/discharge residual is zero. Battery degradation, cycle-dependent efficiency loss, replacement cost, and lifetime constraints are not represented; the extension therefore supports short-horizon operational and energy-balance feasibility, not lifecycle-economic or investment viability. After 240.0 MWh of curtailment, 6,268.8 MWh of renewable energy is accepted by the terminal. Of this amount, 4,947.774 MWh is supplied directly to terminal loads and 1,321.026 MWh is used for battery charging. The battery subsequently delivers 1,192.226 MWh, producing a round-trip storage loss of 128.8 MWh. Grid-to-battery charging and exported energy are both zero. Direct renewable supply and battery discharge therefore provide $\mathrm{4,947.774}+\mathrm{1,192.226}=\mathrm{6,140.0}$ MWh of renewable service to terminal loads. The renewable-storage comparison is a secondary exploratory paired analysis outside the prespecified 18-test confirmatory benchmark family. Holm adjustment is applied separately across its six reported outcomes.

The extension uses the same 30 operational replication seeds and paired renewable-error streams (Table 6). Relative to the grid-only proposed policy, the mean paired reductions are 6,140 MWh for net grid import (95% CI 5,970.062–6,309.938), 8.7 MW for peak import (7.849–9.551), 3,193 tCO₂e for modelled emissions (3,091.526–3,294.474), 0.14 for maximum IPI (0.119–0.161), and 0.30 h for composite delay (0.149–0.451); berth-utilization efficiency increases by 0.70 percentage points (0.318–1.082). All six Holm-adjusted p-values are below 0.001. Technology-specific profiles, uncertainty models, correlations, the exact battery-energy balance, and the paired inference are reported in Supplementary Section S9 and Tables S17 and S38.

5.8. Sensitivity, Weight, and Threshold Robustness

The one-at-a-time sensitivity, weight, coupled threshold-pair, emission-factor, candidate-set, and scenario analyses are engineering robustness checks. Candidate-set sizes of 16 and 24 change every proposed-policy headline metric by less than 1.5% relative to the selected size of 12. Descriptive mean ordering is unchanged across Scenarios 1–4; no confirmatory inference is assigned to the S1–S3, sensitivity, weight, threshold-pair, emission-factor, or candidate-set analyses. Varying αg over 0.35–0.70 t/MWh and βq over 0.12–0.24 t/(vessel·h) preserves the policy ordering of total modelled emissions. Complete values are reported in Supplementary Tables S18–S20, S22–S23, and S30. The IPI-weight analysis is a post hoc recalculation on fixed policy trajectories; the controller is not re-optimized for each alternative weight vector. The emission-factor analysis is likewise a fixed-trajectory inventory sensitivity. Changes in ${\alpha }_g$ and ${\beta }_q$ alter the reported emission accounting and post-processed IPI values but do not represent closed-loop policy sensitivity after re-optimization.

The threshold analysis varies the critical-transition and instability boundaries as a coupled pair while preserving a 0.40 separation, with I_inst = I_crit + 0.40. It does not isolate the individual effect of either threshold and is not evidence for universal constants. Complete values are reported in Supplementary Table S20.

5.9. Ablation Analysis

Table 7 reports descriptive outcomes for independently re-optimized closed-loop controller configurations. Supplementary Table S28 provides paired confidence intervals, effect sizes, raw p-values, and Holm-adjusted p-values, while Table S39 reports replication-level IPI feasibility. Holm adjustment is applied jointly across the 25 exploratory ablation comparisons and separately from the 18 confirmatory benchmark tests.

Relaxing a mathematical constraint enlarges the feasible set, but the independently re-optimized NSGA-II populations and closed-loop state trajectories can diverge after the first altered decision. The observed differences therefore combine component removal, heuristic search-path effects, and state-path effects and are not interpreted as exact causal decompositions.

The without-IPI-constraint configuration removes Equation (87) while retaining IPI monitoring and the upward-crossing trigger. The without-ramp-rate configuration removes Equation (83) alone. The without-coupled-scheduling configuration optimizes shore-power and equipment decisions sequentially. The objective-only stripped controller removes the IPI constraint and event trigger, the ramp-rate constraint, and joint shore-power/equipment coordination while retaining the three objectives and standard physical and service-preservation bounds.

The without-event-triggered-refinement configuration also removes approximately 65,434 additional objective evaluations per replication. Its reported difference therefore represents the combined effect of refinement timing and additional computation rather than the isolated informational value of the IPI crossing signal.

5.10. Scientific and Practical Implications

The results distinguish contractual and operational grid compliance from the broader operational-stress state represented by the IPI. The index adds queue and modelled emission consequence to grid loading, but it should not be read as a power-flow or transient-stability measure. Its practical role is to identify a high-stress scheduling regime in which congestion and electrical exposure occur together.

From an operational perspective, the framework provides a single decision layer for berth priorities, crane deployment, shore-power allocation, flexible demand, renewable supply, storage, and grid constraints. The strongest evidence is provided jointly by the paired Scenario 4 comparison, the same-space IPI ablation, the robustness analyses, and the external empirical benchmark validation. The latter shows that the principal operating scales are consistent with evidence from real container terminals, while the paired experiments establish the comparative behaviour of the proposed controller under controlled conditions. Port-specific deployment would then require local calibration rather than reformulation of the underlying framework.

6. Conclusions

This study developed an instability-aware, simulation-based digital-twin control framework for coordinated vessel operations and terminal-energy management under congestion and time-varying grid constraints. Its central contribution is the integration of a policy-weighted Instability Propagation Index within a rolling-horizon joint logistics–energy controller that simultaneously coordinates berth allocation, crane deployment, vessel sequencing, shore-power scheduling, flexible electrical demand, renewable generation, and battery operation.

Within the 720-h Scenario 4 coupled congestion-and-grid-limit benchmark, the proposed controller achieved the lowest composite vessel delay, final queue, operational energy demand, modelled grid-plus-queue emissions, peak grid import, and mean replication-level maximum IPI, together with the highest berth-utilization efficiency. Relative to the operation-focused, energy-focused, and decoupled benchmark policies, operational energy was reduced by 8.9–21.2%, modelled emissions by 8.5–21.0%, peak grid import by 5.7–15.1%, and mean maximum IPI by 24.8–42.3%. These improvements were supported consistently by paired parametric tests, sign-flip permutation tests, and BCa bootstrap confidence intervals. The same-space without-IPI ablation further showed that the hard IPI constraint contributes materially to the observed operational and electrical improvements, while the delay-only diagnostic provided an independent optimization-quality reference within the same decision space.

The renewable-storage extension further demonstrated the ability of the framework to coordinate low-carbon supply and storage without compromising terminal service. With a 23 MW renewable portfolio and an 8 MW/16 MWh battery, net grid import decreased from 14,510 to 8,370 MWh, peak import from 54.9 to 46.2 MW, modelled emissions from 8,040 to 4,847 tCO₂e, and mean maximum IPI from 0.82 to 0.68, while preserving gross terminal energy requirements and terminal state of charge. These results support the operational feasibility of coordinated renewable and storage integration over the simulated horizon; lifecycle battery degradation, replacement economics, and long-term investment performance remain outside the present assessment.

Importantly, the benchmark is synthetic in stochastic realization but is not an unconstrained hypothetical test case. Its principal operating scales are externally grounded in published evidence from real container-terminal environments, including digital-twin-based vessel and terminal coordination, empirically studied quay-crane configurations, and reported container-vessel shore-power demand ranges. This external empirical grounding strengthens the physical and operational plausibility of the benchmark while preserving the controlled and reproducible conditions required for rigorous policy comparison.

The framework should therefore be interpreted as a validated benchmark-level control architecture rather than a calibrated digital replica of a specific commercial port. Port-specific deployment would require local calibration of vessel-arrival processes, service-time distributions, equipment duty factors, shore-power profiles, grid operating limits, emission coefficients, IPI weights and thresholds, and renewable-resource characteristics using terminal-operating-system records, AIS histories, equipment telemetry, smart-meter measurements, grid agreements, and site-specific operational data. Such calibration represents the next deployment stage rather than a reformulation of the underlying framework.

The present model does not resolve distribution-network power flow, voltage and frequency dynamics, protection behaviour, communication failures, cyberattack, equipment degradation, electricity-market participation, battery ageing cost, or container-level cargo flows. These extensions define the principal directions for future development toward a fully deployed operational digital twin.