Instability-Aware Digital Twins for Robust Co-Optimization of Operational Carbon Emissions, Thermal Comfort, and Structural Serviceability Reliability under Uncertainty in Smart Buildings

1. Introduction

The global transition towards low-carbon, healthy, and resilient infrastructure has fundamentally changed the way modern buildings are designed, operated, and evaluated. Buildings are no longer regarded as static energy-consuming assets exposed to fixed environmental and operational conditions. They increasingly function as dynamic cyber-physical systems in which energy demand, indoor environmental quality, occupant behavior, equipment performance, grid-carbon intensity, and structural response evolve simultaneously under uncertainty. This transformation has created a need for decision frameworks that can address operational sustainability, occupant comfort, indoor-air quality, and long-term structural performance within a unified and physically interpretable environment.

Building energy optimization, digital-twin technology, uncertainty quantification, indoor-environmental assessment, and structural reliability analysis have each advanced substantially. However, their integration remains limited, particularly when operational decisions must remain effective under climatic variability, stochastic occupancy, changing grid conditions, equipment uncertainty, sensor noise, and progressive structural degradation (Tian et al. 2018; Cespedes-Cubides and Jradi 2024; Torzoni et al. 2024; Xu et al. 2025; Samaei and Riffat 2026a). As a result, many apparently optimal building-control strategies remain dependent on nominal assumptions and may perform differently when implemented under realistic variability.

Digital-twin technology provides a promising foundation for addressing this challenge. A building digital twin can combine physics-based or data-driven models with observations from temperature, humidity, CO₂, airflow, energy, and structural sensors to maintain an updated numerical representation of building operation. Recent studies have demonstrated the potential of digital twins for energy prediction, operational control, occupancy monitoring, information integration, fault detection, and predictive building management (Cespedes-Cubides and Jradi 2024; Boukaf et al. 2024; Hwang et al. 2025; Sayed et al. 2025; Samaei and Riffat 2026a). White-box and hybrid digital-twin implementations have also shown that physics-based state estimation and real-time optimization can improve the transparency and adaptability of building-energy control (Guerrero Ramírez et al. 2025; Nuevo-Gallardo et al. 2025).

Despite this progress, most building digital twins remain concentrated on energy efficiency, equipment operation, and short-term thermal performance. Indoor-air quality may be monitored, but it is not always represented as a dynamic state with explicit ventilation constraints. Structural response is even less commonly incorporated as a time-dependent operational consideration. Long-term structural serviceability reliability and the instability of operational decisions under uncertainty are therefore rarely integrated as co-equal components of the same digital-twin decision process (Tahmasebinia et al. 2023; Arsecularatne et al. 2024; Yeom et al. 2024).

Operational carbon emissions, thermal comfort, indoor-air quality, and structural serviceability are interdependent. Multi-objective building studies have shown that energy use, carbon emissions, lifecycle cost, thermal comfort, and envelope performance frequently compete with one another, requiring Pareto-based analysis rather than single-objective optimization (Chen et al. 2022; Chen et al. 2023; Benaddi et al. 2024; Kang et al. 2024a; Yang and Jradi 2025). Recent reviews similarly identify multi-objective optimization as an important tool for energy-efficient building design, retrofit, and control (Alexakis et al. 2025; Li et al. 2025). However, the majority of these studies focus on energy, carbon, cost, and comfort. Structural reliability is commonly omitted or treated as a fixed design-stage requirement rather than as a time-dependent operational objective.

Indoor-air quality introduces an additional operational constraint. Reducing total supply airflow, outdoor-air intake, or HVAC operating time can lower fan and conditioning energy, but it can also increase indoor CO₂ concentration and reduce ventilation adequacy. Conversely, increasing outdoor airflow may improve IAQ while increasing heating, cooling, and fan electricity. The timing of ventilation is also important when grid carbon intensity varies over time. A carbon-aware controller may use periods of relatively low grid-carbon intensity for enhanced ventilation and use the resulting IAQ buffer during more carbon-intensive periods. Such control is acceptable only when minimum outdoor-air requirements and absolute indoor CO₂ limits remain active. Therefore, ventilation cannot be represented solely as an energy-control variable without an accompanying dynamic IAQ model.

Building operation may also influence structural serviceability over long time horizons. Wide indoor-temperature deadbands, intermittent HVAC operation, external temperature variation, and solar exposure can alter the thermal history of restrained structural and façade-related components. These effects may change serviceability demand through restrained thermal expansion, slab restraint, façade–frame interaction, and temperature-dependent deformation. At the same time, uncertain resistance and progressive material degradation can reduce the available serviceability margin.

The structural objective considered in this study does not represent sudden global collapse or replace conventional code-based structural design. Instead, it represents the probability that the time-dependent demand of an equivalent reinforced-concrete structural component exceeds its degraded serviceability resistance over a 20-year horizon. This distinction allows operational thermal histories to be linked to long-term structural performance without making unsupported claims regarding whole-building collapse.

Structural digital-twin research has demonstrated that model updating, condition monitoring, and probabilistic integrity assessment can support reliability-informed decisions in civil infrastructure (Torzoni et al. 2024; Sakr and Sadhu 2024; Li and Brennan 2024; Riffat et al. 2026). Probabilistic deterioration models for reinforced-concrete structures further show that uncertainty in material properties, resistance, loading, environmental exposure, and degradation must be represented explicitly in long-term reliability assessment (Val et al. 2025; Samaei and Riffat 2026a; Riffat et al. 2025a; Riffat et al. 2025b). Nevertheless, these structural advances are rarely linked directly to building-energy, IAQ, and comfort optimization.

Uncertainty further complicates the operational decision problem. Building performance depends on uncertain weather, solar irradiance, occupancy, internal equipment loads, infiltration, envelope conductance, thermal capacitance, heat-pump performance, fan efficiency, grid-carbon intensity, ventilation demand, material properties, structural restraint, resistance, degradation, and sensor accuracy. Uncertainty analysis has long been recognized as essential in building-performance assessment (Tian et al. 2018). More recent studies have extended uncertainty treatment to machine-learning-based energy prediction, occupant-centered control, robust model-predictive control, retrofit assessment, and urban climate variability (Kang et al. 2024b; Yamamoto et al. 2024; Zhang et al. 2025; Choi et al. 2025; Xu et al. 2025). Sensitivity analysis is also important because it identifies the uncertain inputs that dominate individual objectives and reveals interactions that cannot be captured through one-factor-at-a-time assessment (Pang et al. 2020).

However, uncertainty is frequently evaluated after optimization rather than embedded within the optimization process itself. In such formulations, candidate solutions are initially selected using nominal or mean conditions and are only subsequently examined under uncertainty. This sequence can preserve policies whose apparent performance is favorable only because uncertainty was excluded during decision generation.

A more fundamental limitation concerns decision instability. Even when expected values, confidence intervals, or sensitivity indices are available, the location of the selected Pareto solution may change considerably across stochastic realizations or independent optimization runs. This issue is particularly important in nonlinear multi-objective systems, where small changes in weather, occupancy, equipment performance, or material parameters can alter both objective values and the geometry of the Pareto front.

Previous studies have cautioned that optimization algorithms in building-energy applications may produce ineffective or misleading solutions when the problem formulation, uncertainty structure, and algorithmic settings are not treated carefully (Si et al. 2019). Robust optimization studies have begun to address uncertainty in operational control and retrofit decisions (Yamamoto et al. 2024; Zhang et al. 2025). Nevertheless, explicit and dimensionally consistent quantification of both objective-response instability and Pareto-solution migration remains limited, particularly when operational carbon, comfort, IAQ, and structural serviceability reliability are considered together.

Expected performance and operational robustness are not equivalent. A policy may be Pareto efficient in terms of expected carbon emissions, discomfort, and reliability while remaining highly sensitive to uncertain inputs. Similarly, a solution may exhibit acceptable mean performance but migrate substantially in the operational decision space when the optimization is repeated. A decision-grade digital twin should therefore quantify not only the expected performance of an operational policy but also the stability of its outcomes and decision variables.

To address these limitations, this study develops an instability-aware digital twin framework for the robust co-optimization of operational carbon emissions, occupied thermal comfort, and structural serviceability reliability under uncertainty. Indoor-air quality is represented explicitly through a dynamic CO₂ balance, minimum outdoor-air requirements, demand-controlled ventilation, and an absolute occupied CO₂ constraint. The framework combines stochastic uncertainty propagation, recursive digital-twin state updating, probabilistic reliability assessment, surrogate-assisted multi-objective optimization, and explicit objective- and decision-instability metrics.

The physical model is formulated as a multi-rate, sequentially coupled thermal–energy–IAQ–structural system. The term sequentially coupled is used because operational decisions first determine indoor thermal conditions, ventilation, IAQ, and electrical demand. The resulting indoor, outdoor, and solar thermal histories are then propagated through an equivalent structural serviceability model. The calculated reliability response is subsequently returned to the multi-objective optimization layer. The formulation does not imply a fully bidirectional physical coupling at every numerical time step.

The methodology is evaluated using a controlled numerical benchmark representing a six-story office building with a gross floor area of 10,200 m². The building is represented by six equivalent thermal and IAQ zones and is served by a recirculating all-air VAV system supplied by a reversible air-to-water heat pump. The benchmark includes a one-year operational horizon for carbon, comfort, and IAQ evaluation and a 20-year horizon for structural degradation and serviceability reliability assessment.

Only operational policy variables are optimized: occupied heating setpoint, occupied cooling setpoint, high-carbon CO₂ target, and carbon-aware ventilation gain. Building geometry, envelope properties, structural characteristics, internal-load densities, and installed HVAC capacities remain fixed. The study therefore addresses operational control rather than simultaneous building-design and operational-design optimization.

The proposed framework is compared with six matched ablation baselines. B1 and B2 evaluate nominal energy–comfort and nominal-operational-input, reliability-coupled optimization; B3 adds the full modeled uncertainty set without instability criteria; B4 adds operational objective-instability screening without structural reliability; B5 retains uncertainty, reliability, and objective instability but omits decision instability; and B6 retains uncertainty, reliability, and decision instability but omits objective instability. Objective- and decision-instability response engines are trained from an independent pilot ensemble and frozen before the 30 paired final runs. The controls estimate framework-level effects under matched optimizer-evaluation and paired-seed budgets. Because each formulation generates its own Pareto archive, the comparisons are interpreted as matched framework comparisons rather than an exact causal decomposition of independent terms.

The study is based on a controlled numerical benchmark rather than a field-deployed building. Synthetic observations are generated using a separate numerical truth model and are contaminated with prescribed measurement noise before state assimilation. The results therefore provide numerical evidence of methodological performance and should not be interpreted as field validation.

The principal contributions are sixfold.

• Operational carbon emissions, occupied zone-time thermal discomfort, and the 20-year structural serviceability limit-state exceedance probability are treated as co-equal objectives, while IAQ remains a hard feasibility condition.
• Total supply airflow and outdoor-air ventilation are represented separately, enabling carbon-aware ventilation shifting without relaxing minimum outdoor-air or occupied CO₂ requirements.
• Annual carbon and comfort are evaluated using Monte Carlo simulation, whereas structural reliability is evaluated using FORM and verified through design-point importance sampling.
• Objective instability is calculated from normalized stochastic dispersion in carbon, discomfort, and conditional serviceability exceedance probability; decision instability measures normalized migration of matched Pareto policies across 30 independent runs.
• A 21-component augmented EKF/prognostic vector is propagated and corrected within a reduced physics model that also contains six auxiliary humidity-ratio states. Synthetic observations are generated by a structurally distinct truth model containing 32 primary states and six auxiliary humidity-ratio states.
• ARD Matérn-5/2 Gaussian-process surrogates accelerate NSGA-II, followed by full-model re-evaluation of retained nondominated policies and instability-aware robust ranking.

Section 2 reviews the relevant research streams. Section 3 presents the methodology. Section 4 defines the benchmark and numerical implementation. Section 5 reports verification, Pareto, uncertainty, IAQ, instability, sensitivity, and ablation results. Section 6 presents the conclusions, limitations, and future research directions.

2. Literature Review

The development of decision-oriented frameworks for smart-building operation has evolved through several overlapping research areas, including building digital twins, multi-objective performance optimization, indoor-environmental control, uncertainty quantification, robust decision-making, structural digital twins, and time-dependent reliability assessment. Each area has generated important methodological advances; however, the literature remains fragmented.

Energy-oriented building digital twins rarely incorporate structural reliability or long-term degradation. Structural digital twins generally focus on condition monitoring and integrity management without addressing operational carbon emissions, occupant comfort, or indoor-air quality. Multi-objective building optimization commonly considers energy, carbon, cost, and comfort but seldom evaluates whether the resulting operational policies remain stable under stochastic variability. Robust optimization studies address uncertain performance but do not generally combine recursive digital-twin updating, IAQ constraints, structural serviceability reliability, and explicit decision-instability analysis.

This section reviews the principal research streams relevant to the present study and identifies the methodological gap addressed by the proposed instability-aware framework.

Digital twins have become an important component of the digital transformation of building operation. In the building sector, a digital twin is generally understood as an evolving virtual representation of a physical asset or process that combines computational models, operational data, sensing, state estimation, and feedback between predicted and observed system behavior.

Recent reviews indicate that building digital twins are increasingly used to improve operational energy efficiency, integrate heterogeneous sensor information, support predictive control, detect system faults, monitor occupancy, and assist building-management decisions (Cespedes-Cubides and Jradi 2024; Boukaf et al. 2024; Arsecularatne et al. 2024). Early implementations were primarily focused on energy modeling and HVAC control. These systems typically connected building-energy models with historical or real-time measurements to improve demand prediction, operational scheduling, fault detection, and demand-response performance.

More recent research has extended digital-twin architectures towards information fusion, occupancy monitoring, intelligent building-management systems, and real-time operational optimization. Sayed et al. (2025) and Hwang, Kim and Yoon (2025) demonstrated the growing role of integrated sensing and information-processing architectures in building energy management. Li et al. (2023) showed how operational building data can be linked with predictive analytics through an intelligent monitoring platform.

White-box and hybrid digital twins have also received increasing attention because they retain physical interpretability while permitting recursive calibration and operational control. Guerrero Ramírez et al. (2025) developed a white-box building-energy digital twin for blind control and passive-heating applications. Nuevo-Gallardo et al. (2025) examined real-time digital twins for energy optimization and predictive building control. Yeom et al. (2024) connected digital-twin and extended-reality technologies for building energy management, illustrating the potential for digital twins to become interactive decision-support environments rather than static simulation tools.

Despite this progress, the scope of most building digital twins remains focused on energy consumption, equipment operation, occupancy, and short-term thermal response. Reviews by Tahmasebinia et al. (2023) and Arsecularatne, Rodrigo and Chang (2024) indicate that structural degradation, serviceability reliability, and long-term integrity are rarely incorporated into the same operational twin.

Indoor-air quality is also frequently treated as a monitored output rather than as a dynamic state that constrains operational optimization. Temperature, humidity, and CO₂ sensors may be integrated within the digital-twin architecture, but ventilation energy, occupancy-generated CO₂, minimum outdoor-air requirements, and time-varying IAQ targets are not always represented explicitly within the optimization problem.

This omission is significant because energy, comfort, IAQ, and structural response evolve over different but interacting time scales. A digital twin that optimizes only energy consumption may select wider temperature deadbands, reduced HVAC operation, or reduced airflow without evaluating their consequences for occupant comfort, CO₂ concentration, structural temperature variation, or long-term serviceability demand.

A decision-oriented building digital twin should therefore move beyond energy prediction and operational visualization. It should integrate physical-state evolution, indoor-air quality, uncertainty propagation, recursive updating, long-term reliability, and stability of the resulting control decisions.

Multi-objective optimization has become a major methodology in building-performance research because building operation and design involve competing objectives. Energy consumption, operational carbon emissions, thermal comfort, indoor-air quality, lifecycle cost, daylighting, renewable-energy utilization, and retrofit performance cannot generally be optimized independently.

Chen, Tsay and Ni (2022) proposed an integrated multi-objective framework combining carbon emissions, thermal comfort, and global cost. Their findings demonstrated the importance of Pareto-based decision-making when environmental, economic, and occupant-related criteria conflict. Chen, Tsay and Zhang (2023) extended this perspective by considering carbon emissions from a whole-life-cycle viewpoint.

More recent studies have broadened the use of multi-objective optimization in building systems. Yang and Jradi (2025) examined the balance between CO₂ emissions and thermal comfort in a central heat-pump system. Kang et al. (2024a) developed an integrated passive-design method for zero-carbon buildings considering carbon emissions, economics, and thermal comfort. Benaddi, Boukhattem and Tabares-Velasco (2024) optimized building-envelope components using economic, environmental, and thermal-comfort criteria.

These studies confirm that carbon and thermal comfort are strongly interconnected. Reducing HVAC operation or widening temperature deadbands may lower energy use and emissions but increase discomfort. Conversely, tighter thermal regulation may improve comfort while increasing heating, cooling, and fan electricity.

Systematic reviews further demonstrate the rapid development of this field. Alexakis et al. (2025) reviewed genetic-algorithm-based multi-objective optimization for energy-efficient building retrofits. Li, Shari and Ab Kadir (2025) provided bibliometric evidence that multi-objective building-performance optimization has become a central research direction.

Nevertheless, most studies continue to focus on energy, carbon, cost, comfort, and envelope performance. Structural reliability is rarely treated as a dynamic operational objective. It is usually considered implicitly through code-compliant design or as a fixed constraint established before operation.

This separation may be inadequate when operational strategies alter the thermal history of restrained structural and façade-related components. Wide thermal deadbands, intermittent conditioning, solar exposure, and changing external temperatures can affect restrained thermal expansion, slab response, façade–frame interaction, and serviceability demand. If resistance also changes through long-term degradation, the structural margin can evolve throughout the operating life.

Indoor-air quality presents another multi-objective interaction. Reducing outdoor airflow can decrease fan power and conditioning demand but may increase indoor CO₂ concentration. Increasing ventilation can improve IAQ but may raise energy use and carbon emissions, particularly during periods of high grid-carbon intensity. Therefore, ventilation-related energy optimization requires an explicit dynamic IAQ constraint rather than an unconstrained airflow-reduction objective.

The present study extends conventional carbon–comfort optimization by treating operational carbon emissions, occupied thermal discomfort, and structural serviceability reliability as three co-equal objectives while imposing indoor-air quality as a hard feasibility condition.

Indoor-air quality is an essential component of operational building performance because ventilation decisions influence occupant exposure, fan energy, heating and cooling loads, and carbon emissions. Indoor CO₂ concentration is a useful indicator of occupancy-related ventilation adequacy in office buildings, particularly in demand-controlled ventilation applications. However, CO₂ is not a comprehensive measure of indoor-air quality because it does not directly characterize particulate matter, volatile organic compounds, ozone, biological contaminants, or other indoor and outdoor pollutant exposures. Accordingly, the present framework uses CO₂ as an explicitly modeled ventilation-control and feasibility indicator rather than as a complete representation of occupant exposure to all indoor-air pollutants.

A ventilation controller must distinguish between total supply airflow and outdoor airflow. In an all-air VAV system, total supply airflow primarily responds to thermal demand, whereas the outdoor-air component is governed by occupancy, indoor contaminant concentration, minimum ventilation requirements, and equipment limitations. Treating ventilation airflow as a single unconstrained optimization variable may therefore produce physically ambiguous or infeasible control strategies.

Time-varying electricity carbon intensity introduces an additional opportunity for operational control. Outdoor-air demand can be partially shifted towards periods of lower grid-carbon intensity by temporarily lowering the indoor CO₂ target and increasing ventilation. The resulting IAQ buffer may permit a higher, but still acceptable, CO₂ target during carbon-intensive periods. Such a strategy is valid only if minimum outdoor-air requirements and an absolute occupied CO₂ limit remain continuously enforced.

Existing energy and comfort optimization studies do not always represent these ventilation dynamics explicitly. Indoor-air quality may be excluded, treated as a fixed airflow assumption, or represented by a simplified constraint disconnected from occupancy and control actions. Consequently, some apparently low-carbon strategies may obtain energy savings by reducing ventilation rather than by improving the timing and efficiency of operation.

The present study addresses this limitation through a dynamic zone CO₂ balance, stochastic occupancy, explicit minimum outdoor-air calculation, demand-controlled ventilation, carbon-aware CO₂ target adjustment, and a hard occupied CO₂ limit. IAQ is not introduced as a fourth objective because the purpose is not to trade indoor-air quality against carbon, comfort, or reliability. Instead, it defines the admissible operating space within which the three primary objectives are optimized.

Recent studies have moved toward explicit integration of energy, thermal comfort, and indoor-air-quality objectives. Al Mindeel, Spentzou and Eftekhari (2024) reviewed multi-objective building formulations treating these domains as competing performance targets. Andersen and Jradi (2025) demonstrated ventilation-system optimization combining indoor-air quality, energy cost, and electricity-related environmental impact. Ageli and Chassiakos (2026) further identified real-time integration of indoor-environmental sensing, energy-performance modeling, and operational control as an important development pathway for building digital twins. Nevertheless, these studies do not incorporate time-dependent structural serviceability reliability together with separate objective- and decision-instability metrics.

Uncertainty is one of the most persistent challenges in building-performance modeling. Tian et al. (2018) reviewed uncertainty analysis in building-energy assessment and showed that uncertainty arises from weather, occupancy behavior, schedules, physical parameters, equipment performance, model assumptions, and operational conditions. Their review demonstrated that deterministic evaluation can misrepresent building performance when uncertain inputs are not propagated explicitly.

Recent studies have extended uncertainty analysis towards occupant-centered control, machine-learning-based prediction, robust optimization, and urban-scale climate modeling. Kang et al. (2024b) investigated occupant-related uncertainty in office-building energy simulation and comfort assessment, demonstrating that occupancy behavior can materially alter thermal and energy predictions. Xu et al. (2025) systematically reviewed uncertainty quantification in machine-learning-based building-energy models and emphasized the importance of representing uncertainty in both prediction and decision systems.

Choi, Bergés and Pozzi (2025) quantified uncertainty in building energy use using a probabilistic spatiotemporal model of urban temperature, highlighting the importance of climatic and spatial variability. These studies demonstrate that uncertain weather and occupancy conditions affect not only predicted energy demand but also the confidence with which operational decisions can be implemented.

Sensitivity analysis is an important companion to uncertainty propagation. Pang et al. (2020) critically reviewed sensitivity analysis in building-performance modeling and showed that identifying influential parameters is essential for credible simulation, calibration, uncertainty reduction, and decision-making. A digital twin may propagate many uncertain inputs, but without global sensitivity analysis it remains unclear which variables dominate carbon, comfort, IAQ, or structural reliability.

Robust optimization studies have begun to embed uncertainty more directly in operational and retrofit decisions. Yamamoto et al. (2024) developed a model-predictive-control-based robust optimization method for smart apartment buildings under uncertainty. Zhang et al. (2025) proposed an interpretable machine-learning method for robust building-energy retrofit assessment.

These studies represent an important transition from nominal optimization towards uncertainty-aware decisions. However, they do not fully address the stability of Pareto-optimal policies when the optimization is repeated under different stochastic conditions. They also generally exclude structural serviceability reliability and do not combine uncertainty propagation with recursive digital-twin updating and explicit IAQ constraints.

Another methodological issue concerns the treatment of rare-event probabilities. Monte Carlo sample sizes that are sufficient for estimating expected annual energy or comfort may be inadequate for estimating small structural exceedance probabilities. Reliability analysis therefore requires methods such as FORM, importance sampling, or other variance-reduction techniques rather than assuming that convergence of annual performance metrics establishes convergence of structural risk.

The present study distinguishes annual-performance uncertainty propagation from rare-event reliability estimation. Carbon and comfort are evaluated using Monte Carlo simulation, whereas structural serviceability reliability is evaluated using FORM and verified using design-point importance sampling.

Structural digital twins have developed rapidly in civil infrastructure research. Torzoni et al. (2024) proposed a digital-twin framework for civil engineering structures in which computational models, monitoring data, uncertainty representation, and model updating were integrated for structural assessment.

Sakr and Sadhu (2024) reviewed recent developments and future directions in digital twins for structural health monitoring. Their review emphasized the importance of data–model integration, condition assessment, damage detection, predictive maintenance, and uncertainty-aware structural evaluation. Li and Brennan (2024) developed a digital-twin-enabled structural-integrity-management framework and highlighted the relevance of continuously updated models for reliability-informed asset management.

Structural reliability research provides methods for quantifying long-term performance under uncertain material properties, loads, environmental exposure, model form, and deterioration. Val et al. (2025) reviewed probabilistic deterioration modeling for reinforced-concrete structures and demonstrated that uncertainty in degradation and resistance must be incorporated explicitly in long-term reliability assessment. Wang et al. (2024) demonstrated probabilistic uncertainty quantification for the strain response of prestressed reinforced-concrete structures using fractile-based sampling.

Structural digital twins and building-energy digital twins have nevertheless developed largely as separate research streams. Structural twins commonly focus on condition monitoring, damage detection, integrity management, and maintenance planning. Building twins commonly focus on energy prediction, HVAC control, occupancy, and thermal performance.

This disciplinary separation prevents a complete evaluation of how operational energy strategies may influence structural serviceability over time. Aggressive energy-saving policies can alter indoor-temperature trajectories, supply-air patterns, and thermal cycling. External temperature and solar exposure simultaneously affect façade and edge components. Where movement is partially restrained, these effects may influence serviceability stress and deformation.

The magnitude and practical significance of such effects depend on structural configuration, material properties, restraint, environmental exposure, and degradation. Therefore, the appropriate structural interpretation is not a generic claim that HVAC control causes structural failure. Instead, operational thermal histories should be propagated through a clearly defined component-level serviceability model with uncertain resistance and degradation.

The present study adopts this limited and physically interpretable scope. It evaluates the probability that the demand of an equivalent reinforced-concrete edge-beam/slab component exceeds its degraded serviceability resistance over a 20-year horizon. The calculated probability is not interpreted as whole-building collapse probability.

Robust optimization is intended to identify decisions that retain acceptable performance when uncertain conditions change. In building applications, robustness is commonly represented through worst-case performance, expected-value penalties, variance measures, chance constraints, or scenario-based evaluation.

These formulations improve upon deterministic optimization because they recognize that a policy selected under nominal conditions may perform poorly under adverse weather, occupancy, equipment, or model conditions. However, objective uncertainty and decision instability are distinct.

Objective instability describes variation in performance when a fixed policy is exposed to uncertainty. A policy may have low expected carbon emissions but high variability in carbon, discomfort, or reliability response.

Decision instability measures cross-run migration of equivalent Pareto policies. Policies are assigned by one-to-one minimum-cost Hungarian matching in the fixed normalized three-objective space (Kuhn 1955) and retained only when the matching distance does not exceed 0.15. M denotes the number of accepted matches, n_x = 4 is the number of decision variables, and the reference vector is the componentwise mean of the accepted policy vectors. One unit of I_d corresponds to an RMS migration equal to 10% of each variable’s admissible range.

Robustness in multi-objective optimization has also been studied beyond the building domain. Classical formulations evaluate the stability of objective values under perturbations of decision variables or uncertain parameters, whereas decision-space robustness examines whether efficient solutions remain implementable or efficient within a neighborhood. The present decision-instability metric addresses a different question: the repeatability and migration of matched Pareto-policy regions across independently seeded stochastic optimization runs. It therefore complements, rather than replaces, established perturbation-based robustness concepts (Deb and Gupta 2006; Stiglmayr et al. 2022).

Previous research has warned that optimization algorithms in building-energy studies may produce ineffective or misleading solutions when uncertainty, algorithm settings, and problem formulation are not handled carefully (Si et al. 2019). Robust building-control and retrofit studies have started to address uncertainty in decision-making (Yamamoto et al. 2024; Zhang et al. 2025), but explicit quantification of Pareto-policy migration remains limited.

A further issue is dimensional consistency. Carbon dispersion, discomfort dispersion, reliability variation, temperature shifts, airflow changes, and controller gains have different units and scales. They cannot be added directly without explicit normalization. A valid instability-aware framework therefore requires dimensionless objective-dispersion measures and normalized decision-space distances.

The present study distinguishes objective instability from decision instability. Objective instability is calculated using stochastic dispersion normalized by a common pooled percentile range. Decision instability is calculated by matching equivalent Pareto policies across independent optimization runs and measuring normalized RMS migration in the four-dimensional operational decision space.

This distinction allows a policy to be evaluated not only according to expected performance but also according to the stability of its response and the repeatability of its operational settings.

Recursive state updating is a defining characteristic of an operational digital twin. Without data assimilation, a building model remains an open-loop simulation that may drift from observed operation as weather, occupancy, equipment performance, or structural response changes.

Kalman-filter-based approaches provide a systematic method for updating predicted states using noisy observations. Linear Kalman filters are appropriate when state-transition and observation models are linear. Extended or unscented Kalman filters are more suitable when thermal, IAQ, HVAC, or structural models are nonlinear.

In the present benchmark, thermal, CO₂, electrical, and structural observations are generated using a numerical truth model and are assimilated using an extended Kalman filter. This provides a controlled demonstration of recursive state updating while avoiding an unsupported claim of field validation.

Multi-objective optimization under uncertainty can require a large number of repeated simulations. Surrogate models are therefore frequently used to reduce computational cost. Gaussian-process models are particularly attractive because they provide predictive means, uncertainty estimates, and automatic relevance determination.

However, surrogate-assisted optimization introduces the risk that an inaccurate approximation may distort the Pareto front or reliability estimate. Surrogate models should therefore be validated using independent data, and final nondominated policies should be re-evaluated using the complete physics-based model.

The present study adopts separate Gaussian-process models for carbon, discomfort, and structural reliability index. The surrogate models accelerate NSGA-II search, but all retained nondominated policies and the final robust policy are re-evaluated using the complete sequentially coupled model.

Within the literature reviewed for this study, six recurring methodological limitations are evident.

First, building digital twins remain predominantly energy-centered. Although recent studies demonstrate substantial progress in energy prediction, occupancy monitoring, white-box modeling, information integration, and predictive control, structural serviceability reliability is rarely included as an operational objective (Cespedes-Cubides and Jradi 2024; Boukaf et al. 2024; Sayed et al. 2025; Hwang et al. 2025; Nuevo-Gallardo et al. 2025).

Second, indoor-air quality is often monitored but not fully coupled with operational optimization. Dynamic CO₂ generation, minimum outdoor-air requirements, ventilation energy, and time-varying carbon-aware ventilation control are not consistently represented within the same framework.

Third, multi-objective building optimization has not fully incorporated structural serviceability reliability as a dynamic co-equal objective. Existing studies provide strong foundations for carbon, comfort, cost, envelope, and retrofit optimization but generally exclude long-term structural resistance and degradation (Chen et al. 2022; Chen et al. 2023; Benaddi et al. 2024; Yang and Jradi 2025; Alexakis et al. 2025).

Fourth, uncertainty is often evaluated as a separate post-processing step. Prior research demonstrates that weather, occupancy, urban climate, equipment, and model parameters influence building performance (Tian et al. 2018; Pang et al. 2020; Kang et al. 2024b; Choi et al. 2025; Xu et al. 2025), but candidate policies are still frequently generated under nominal conditions.

Fifth, robust optimization studies rarely distinguish objective instability from decision instability. Confidence intervals and sensitivity indices describe uncertain performance, but they do not necessarily show whether the selected Pareto policy remains stable across repeated optimization runs.

Sixth, digital-twin updating, uncertainty propagation, IAQ constraints, structural serviceability reliability, surrogate-assisted optimization, and dimensionally consistent instability analysis are rarely combined within a single decision framework.

The literature search focused primarily on peer-reviewed studies published between 2018 and 2026 using combinations of the terms building digital twin, building-energy optimization, robust multi-objective control, indoor-air quality, structural reliability, uncertainty quantification, and decision stability. Earlier sources were retained for established algorithms and standards. The review positions the methodology and is not presented as a formal systematic review (Table 1).

For the mechanically conditioned benchmark considered here, thermal comfort is evaluated using the PMV/PPD procedure in ISO 7730:2025. BS EN 16798-1:2019 is used as contextual guidance for indoor-environmental input parameters; the benchmark is not presented as a formal regulatory compliance assessment (International Organization for Standardization 2025; British Standards Institution 2019).

The present study differs from closely related prior building digital-twin studies in both objective structure and decision methodology (Table 2). Those studies addressed diagnosis-informed climate-response control, stability preservation under non-stationary disturbances, and cumulative health-exposure control. The present framework introduces three-objective carbon–comfort–serviceability optimization, hard dynamic IAQ constraints, separate objective- and decision-instability metrics, and FORM-based serviceability assessment.

2.1. Research Gap and Contribution

The literature demonstrates that a decision-oriented smart-building framework must do more than predict energy demand or generate a nominal Pareto front. It must represent interacting operational domains, enforce indoor-environmental constraints, propagate uncertainty during policy evaluation, quantify long-term serviceability reliability, update predicted states recursively, and determine whether the resulting decisions remain stable.

Existing building digital twins provide strong foundations for sensing, monitoring, and energy management but rarely include structural serviceability reliability or explicit decision-instability analysis (Tahmasebinia et al. 2023; Arsecularatne et al. 2024; Samaei and Riffat 2026b; Samaei and Riffat 2026c; Samaei and Riffat 2026d).

Existing multi-objective building studies provide useful carbon–comfort–cost and design–performance trade-offs but generally neglect degradation-sensitive structural reliability and dynamic IAQ constraints (Chen et al. 2022; Kang et al. 2024a; Yang and Jradi 2025).

Existing uncertainty studies improve the realism of building-performance prediction but do not commonly quantify the migration of equivalent Pareto policies across repeated optimization runs (Tian et al. 2018; Pang et al. 2020; Xu et al. 2025).

Structural digital-twin studies provide strong foundations for monitoring, updating, deterioration assessment, and reliability-informed management but remain largely disconnected from operational carbon, thermal comfort, and ventilation control (Torzoni et al. 2024; Sakr and Sadhu 2024; Li and Brennan 2024; Val et al. 2025).

The present framework addresses these gaps by integrating operational carbon, occupied thermal discomfort, dynamic IAQ constraints, time-dependent structural serviceability reliability, recursive state updating, uncertainty propagation, and separate objective- and decision-instability measures within one operational decision architecture.

3. Methodology

The framework is formulated as a multi-rate, sequentially coupled thermal–energy–moisture–IAQ–structural decision system. Operational states are advanced at 15-min intervals, temperature, CO₂, and electrical-power observations are assimilated hourly, structural strain is assimilated daily, and resistance degradation is updated annually. Figure 1 summarizes an offline policy-design and evaluation cycle with recursive state updating inside the synthetic benchmark; the complete surrogate-assisted optimization and cross-run instability analysis are not claimed to operate online at each control interval. For each candidate operating policy, the computational procedure follows a defined sequence. First, the four control parameters are applied to the reduced thermal–energy–moisture–IAQ model, and the resulting operational states are recursively updated using synthetic observations and the extended Kalman filter. Second, the candidate is evaluated across 1,000 common-random-number uncertainty realizations to estimate annual operational carbon emissions, occupied thermal discomfort, and compliance with the prescribed IAQ constraints. Third, the resulting thermal histories are propagated through the equivalent structural serviceability model, where the 20-year limit-state exceedance probability is estimated using FORM and independently checked by design-point importance sampling. Fourth, Gaussian-process surrogates support the NSGA-II search, while independently trained and frozen instability-response models provide candidate-level screening information. Finally, retained nondominated policies are re-evaluated using the complete numerical model, objective instability is recomputed from the stochastic response distributions, and decision instability is obtained by matching corresponding Pareto-policy regions across the 30 independent optimization runs. Policies satisfying the IAQ, equipment, airflow, and instability requirements are ranked using the normalized robust-selection score.

3.1. State representation and operational decisions

The augmented digital-twin vector contains 21 components: 20 operational components participating in measurement updating and one prediction-only logarithmic resistance component. Six humidity-ratio variables are advanced as auxiliary physics states outside the EKF vector. The numerical truth model contains 32 primary states and six additional humidity-ratio states. The logarithmic resistance component is propagated annually through the degradation model but is not corrected by operational measurements because no direct resistance, inspection, or condition-state observation is assumed. It is retained for synchronized prognostic bookkeeping rather than interpreted as an observable EKF state.

The optimized policy vector contains the occupied heating setpoint, occupied cooling setpoint, high-carbon CO₂ target, and carbon-aware ventilation gain.

Building geometry, envelope construction, structural properties, internal-load densities, and installed HVAC capacities remain fixed. Total supply airflow responds to sensible thermal demand, whereas outdoor airflow responds to occupancy, CO₂ concentration, and grid-carbon intensity. The optimized occupied setpoints are constant policy parameters; the 1 K h⁻¹ slew limit applies only to implemented transitions between unoccupied and occupied setpoints, using a rate-limited ramp at each 15-min interval.

3.2. Thermal, moisture, energy, and IAQ models

The reduced operational twin uses a two-node resistance–capacitance network for each floor. Outdoor and ground boundary terms are evaluated separately. Infiltration is represented consistently as both an outdoor-boundary sensible conductance and an uncontrolled outdoor-air mass flow in the moisture and CO₂ balances.

The outdoor-boundary term includes wall, window, roof, and infiltration conductances, whereas ground-floor conduction is evaluated separately using the time-varying ground temperature defined in Supplementary Equation (S1). Roof and ground conductances are nonzero only for the sixth and first floors, respectively. Mechanical outdoor air passes through sensible heat recovery; infiltration bypasses the AHU and receives no heat-recovery credit.

The uncontrolled infiltration airflow and the internal and solar gains are evaluated as follows. The complete schedule and gain-partition settings are reported in Supplementary Table S5B.

The same daily infiltration multiplier is applied to the uncontrolled airflow and its sensible conductance. Infiltration enters the sensible balance once through H_inf,z(T_out − T_a,z), while its airflow representation is used only in the moisture and CO₂ mass balances. It is excluded from heat recovery because it bypasses the AHU.

Occupant latent gains are represented through a humidity-ratio balance. The cooling-coil bypass-factor model determines condensate and latent cooling load.

The AHU mixed-air humidity ratio is calculated from the mechanical outdoor-air fraction and recirculated zone air. Uncontrolled infiltration bypasses the AHU and is excluded from this mixed-air calculation.

Saturation humidity ratio is evaluated at a benchmark barometric pressure of 101.325 kPa using the Buck saturation-pressure relation (Buck 1981), as detailed in the Supplementary Material.

The reported specific fan power represents the combined supply- and return-fan system and is applied once. Auxiliary electricity comprises circulation pumps equal to 1.2% of heat-pump compressor electricity, AHU controls and actuators equal to 0.5% of fan electricity, and a fixed whole-building control and standby demand of 1.2 kW that remains active throughout the 8,760-h simulation. Baseline heat-pump COP and fan part-load relationships are:

Heating and cooling electricity are calculated from the positive and negative components of the net whole-building HVAC heat input, QHVAC, k = Σz QHVAC, z, k. Whole-building electrical power includes heat-pump, fan, auxiliary, lighting, and equipment electricity.

Indoor CO₂ is advanced analytically over each 15-min interval using a ppm-consistent source term.

The minimum outdoor-air requirement is adopted as a benchmark design rate rather than as a jurisdiction-specific compliance determination. It is calculated at zone level and summed only when whole-building airflow is reported.

χ_vent(t) = 1 for weekdays 07:00–18:00; χ_vent(t) = 0 otherwise

The minimum mechanical outdoor-air requirement is activated only during weekday preconditioning and occupied periods from 07:00 to 18:00. Outside these periods, mechanical outdoor air is set to zero; uncontrolled infiltration remains active continuously. No protective purge event occurs in the benchmark.

The carbon-aware target is:

All feasible policies satisfy the absolute occupied CO₂ constraint of 1,000 ppm and produce zero occupied zone-hours above the limit. Heating and cooling are mutually exclusive. The PI controller uses a 0.20 K mode-switching hysteresis, a minimum mode duration of 15 min, and integral clamping at the active capacity limit. Total supply airflow is increased before outdoor-air clipping whenever necessary to preserve the 0.60 maximum outdoor-air fraction.

Heating and cooling commands are non-negative and mutually exclusive. The net zone heat input is Q_HVAC,z = Q_h,cmd,z − Q_c,cmd,z, so positive values denote heating and negative values denote cooling. The supply-air reset fractions are f_h = Q_h,cmd/Q_h,max and f_c = Q_c,cmd/Q_c,max.

Annual operational carbon intensity is:

3.3. Thermal comfort

Mean radiant temperature is represented by the reduced-order relation T_r = 0.65T_m + 0.35T_a. The relation was calibrated on 20 policy–season combinations and evaluated on a disjoint 10-policy held-out set, for which RMSE was 0.34°C, mean bias was +0.03°C, maximum error was 0.81°C, and R² was 0.984.

At each occupied time step, the lower and upper operative-temperature limits are obtained by solving the ISO 7730 PMV model for PMV = −0.5 and +0.5 using current relative humidity, seasonal clothing, metabolic rate, air speed, and mean radiant temperature. A bracketed bisection search over 15–32°C with a 10⁻⁴°C tolerance is used; if a root lies outside the bracket, the nearest bracket boundary is applied.

Occupied discomfort is calculated as degree-hours outside the PMV-derived operative-temperature band:

Thermal discomfort is summed over occupied zones and 15-min intervals. It is therefore reported as occupied zone-degree-hours per year rather than as a single whole-building clock-hour measure.

3.4. Structural serviceability reliability

The structural model represents an equivalent south-façade reinforced-concrete edge-beam/slab strip. The term $\frac{GHI}{20}$ is an effective calibrated solar-temperature forcing rather than a direct prediction of south-façade surface temperature. The indoor, outdoor, and solar weights and the 6-h time constant were calibrated against the orientation-resolved two-node structural thermal model used in the numerical truth model.

The structural reliability objective represents the probability that the calculated demand of an equivalent reinforced-concrete south-façade edge-beam/slab component exceeds its degradation-dependent serviceability resistance during the 20-year assessment horizon. In engineering terms, this probability quantifies the likelihood of exceeding the prescribed component-level serviceability stress margin under uncertain thermal restraint, material properties, resistance, degradation, and loading. It is not a prediction of whole-building structural failure, nor does it constitute a separate code-compliance assessment of cracking, deflection, fatigue, or ultimate strength. The corresponding reliability index provides a probabilistic representation of this defined limit state. FORM is implemented using the Hasofer–Lind–Rackwitz–Fiessler procedure (Hasofer and Lind 1974; Rackwitz and Fiessler 1978), with probability estimates independently verified using 50,000 design-point importance samples.

Structural strain is retained as a daily constrained latent state rather than as an independently persistent autoregressive state. Its prior mean is re-centered on the restrained thermal-strain relationship in Equation (22), and process noise represents unresolved façade–frame and local-model discrepancy. The corrected strain in Equation (23), rather than an uncorrected temperature-only demand, is propagated into the serviceability limit state. Temporal persistence is carried by the structural-temperature state.

3.5. Uncertainty propagation, FORM, and state updating

The baseline analysis contains 22 independent primary uncertainty factors spanning climate, occupancy, internal loads, grid carbon, envelope conductance, capacitance, HVAC performance, ventilation, material properties, resistance, degradation, loading, and model form. Measurement noise is handled separately. Complete process distributions and measurement-noise specifications are reported in Supplementary Tables S11 and S12.

Annual carbon and discomfort are estimated from 1,000 common-random-number realizations. Global sensitivity indices are estimated using the Saltelli design and first- and total-order estimators (Saltelli et al. 2010). The serviceability probability is evaluated separately by FORM and design-point importance sampling.

The full numerical benchmark includes 22 primary uncertainty factors for operational and structural performance assessment. The 12-dimensional FORM representation is a reduced reliability-analysis space rather than a replacement for the complete uncertainty model. It retains the variables represented explicitly in the structural limit-state calculation and its prescribed time-history construction, while the remaining operational uncertainty factors continue to be represented in the Monte Carlo assessment and in the thermal histories supplied to the reliability model. The influence of excluding variables from the explicit FORM coordinate set is examined through the reliability-oriented screening procedure described in Supplementary Section S7. Accordingly, the reduced FORM dimension should not be interpreted as the removal of those uncertainties from the overall computational framework. The temporal construction, probabilistic transformations, importance-sampling mixture proposal, convergence tolerances, variable-screening procedure, and multi-start design-point checks are documented in Supplementary Section S7. Non-normal variables are transformed using their marginal cumulative distribution functions and the Nataf correlation model where dependence is prescribed (Liu and Der Kiureghian 1986), while the invariant second-moment formulation follows Hasofer and Lind (1974).

Recursive state estimation follows the extended nonlinear form of the Kalman filtering framework (Kalman 1960). Jacobians of the linear RC, CO₂, and structural-temperature transitions are evaluated analytically. Derivatives of nonlinear electrical-power and structural observation functions are evaluated by central finite differences using a perturbation equal to 10⁻⁵ of the corresponding normalized state range.

Weather, occupancy, lighting schedules, equipment schedules, and sampled process multipliers are treated as known exogenous inputs in each synthetic assimilation experiment. This choice isolates state-estimation performance and is more favorable than field deployment, where occupancy, internal loads, and some boundary conditions must themselves be inferred and may be biased or delayed.

Operative temperature is a derived output rather than a directly assimilated observation. Its open-loop and posterior RMSE values are reported without innovation or NIS statistics.

The EKF-corrected structural temperature and strain define the posterior initial condition for prognostic assessment. Future 20-year reliability trajectories are propagated from this posterior state using sampled structural properties, degradation parameters, annual operational histories, and structural process discrepancy. No future measurements are assumed during the reliability projection, and the logarithmic resistance component evolves only through the prescribed degradation model.

3.6. Surrogate-assisted optimization and instability metrics

Uncertainty-aware Gaussian-process regressors with automatic-relevance-determination Matérn-5/2 kernels predict annual carbon, occupied zone-time discomfort, and reliability index β (Rasmussen and Williams 2006). The uncertainty-aware modeling pool contains 1,000 initial scrambled-Sobol policies and 250 adaptively enriched policies, partitioned into 1,000 training and 250 locked-validation cases; a separate 250-policy set is reserved for untouched testing. A second nominal-operational dataset contains 750 complete-model policies: 500 training, 125 locked validation, and 125 untouched test cases. The nominal models predict C_N, T_N, and β_N for B1–B2. Five-fold cross-validation is applied only to training data, and each hyperparameter fit uses 10 maximum-marginal-likelihood restarts.

Objective- and decision-instability response engines are developed from an independent pilot ensemble comprising 20 NSGA-II runs with seeds 20001–20020, a population of 100, 120 generations, and 12,000 evaluations per run. The pilot archives contain 2,724 raw entries, 1,024 unique anchor policies, and 768 matched policy regions with at least 12 accepted matches. Gaussian-process models for IoCT, I_o, and I_d are trained on 614 regions, validated on 77 locked regions, and assessed once on 77 held-out policy regions. The held-out regions are disjoint from the fitting and calibration sets but originate from the same 20-run pilot ensemble; they therefore assess generalization to unseen policy regions within the pilot pool rather than transfer to a separately generated optimizer ensemble. All three instability models are frozen before the 30 paired final framework runs with seeds 21001–21030, and no final-run archive contributes to their fitting.

During the B4–B6 and Proposed searches, the frozen pilot predictors Ĩ_o,pilot(x) and Ĩ_d,pilot(x) are used only as candidate-response engines. After the 30 paired final runs are complete, I_o,final(x) is recomputed from the complete-model trajectory arrays and I_d,final(x) is recomputed from the matched final-run policy archives. These complete-model and archive-derived values replace the pilot predictions in the 146 → 58 → 52 → 44 screening sequence and in the final J_rob ranking. No final-run instability value is used to refit or update a search surrogate.

The final NSGA-II runs use a population of 100, 150 generations, simulated-binary crossover, polynomial mutation, feasibility-first constraint handling, and 30 paired independent seed sets (Deb et al. 2002). Framework-specific candidate ranking uses the frozen response engines defined above. Every retained archive policy is subsequently re-evaluated with the complete operational model and FORM; no surrogate value is retained as a final reported result.

Objective instability is based on robust normalized response dispersion. Carbon and discomfort use the 1%–99% winsorized dispersion of the 1,000 common-random-number trajectories. For reliability, each operational trajectory is paired with a 2,048-point Owen-scrambled Sobol structural bank; 20 independent scrambles verify inner-bank convergence. The non-thermal service-load history is included only in the structural variable set ξs and is not sampled a second time in ξo. The dispersion of conditional serviceability probabilities defines I_o,R. The policy-level P_SLS values used for optimization and comparison remain FORM estimates verified by design-point importance sampling.

Equivalent policies are matched across runs by one-to-one minimum-cost Hungarian assignment in the common normalized three-objective space (Kuhn 1955). The combined pilot cost cmr(λ) = λ||F̃m − F̃r||₂ + (1 − λ)||x̃m − x̃r||₂ with λ = 0.80 is used only to establish stable pilot-region correspondence, with a normalized objective-distance threshold of 0.15 and at least 12 accepted matches per region. After region correspondence is fixed, the GP-I_d training targets are recomputed using objective-only matching with λ = 1.00. The GP-I_d target and the final reported I_d therefore use the same objective-only metric; λ = 0.80 and 0.50 are retained solely as matching-sensitivity cases. The objective-space matching threshold of 0.15 limits correspondence to policy regions that are sufficiently close in the fixed normalized objective space, thereby reducing the risk of interpreting unrelated Pareto solutions as repeated instances of the same policy region. This threshold is an algorithmic correspondence criterion selected for the benchmark and has no universal physical meaning. Its suitability must be interpreted together with the accepted-match requirement and the matching-sensitivity analysis; it should be recalibrated when the objective ranges, policy domain, or stochastic optimization setting changes. Decision instability is the root-mean-square normalized migration of the four policy variables.

The final robust score applies equal weights to min–max-normalized carbon, discomfort, exceedance probability, objective instability, and decision instability. Fixed bounds obtained from the union of all complete-model post-evaluated B1–B6 and Proposed archives are used for every active selection score, policy comparison, and weight-sensitivity case:

3.7. Baseline and ablation design

All frameworks use the same 30 paired final-run seeds, population size of 100, maximum of 150 generations, 15,000 surrogate evaluations per run, archive limit of 200, feasibility-first constraint handling, and complete-model post-evaluation protocol (Table 3). Matched budget refers to optimizer evaluations and random seeds rather than equal wall-clock cost. Base and instability response engines are fitted before the final paired runs and remain frozen. Hypervolume is calculated after complete-model post-evaluation in the common normalized carbon–discomfort–serviceability-risk space using the fixed reference point (1.10, 1.10, 1.10).

For B4, objective instability is restricted to the operational components IoCT = (I_o,C + I_o,T)/2. Neither P_SLS nor I_o,R enters B4 optimization, active-score selection, or feasibility filtering. The full three-component I_o is used only in frameworks that retain structural serviceability reliability. The objective-instability screening threshold of I_o ≤ 0.30 is adopted as an operational decision criterion for the present numerical benchmark rather than as a universal robustness requirement. Its interpretation depends on the fixed response-normalization ranges, prescribed uncertainty distributions, and admissible control domain. The sensitivity analysis considers alternative thresholds between 0.20 and 0.40, providing a direct check of whether the final policy selection depends strongly on this screening choice.

4. Case Study and Model Implementation

The benchmark is a six-story, 10,200 m² office building represented by six equivalent floor zones (Table 4). It uses a synthetic London–Heathrow temperate-maritime climate sequence constructed to reproduce the annual and occupied-hour statistics reported in Supplementary Table S9. The hourly grid-carbon sequence is likewise a synthetic average operational electricity-carbon benchmark constrained to the statistics in Supplementary Table S10. Figure 2 summarizes the implementation, while detailed parameter, verification, sensitivity, and reproducibility registries are provided in Supplementary Sections S1–S12.

4.1. Building, envelope, HVAC, and schedules

One AHU serves each floor. The VAV system is recirculating and uses no fan-coil units or separate DOAS. A reversible air-to-water heat pump supplies the AHU coils. The controller sequence, PI gains, supply-air reset, capacity limits, outdoor-air logic, schedules, and auxiliary-power formulation are reported in Supplementary Section S2.

4.2. Climate, occupancy, IAQ, and comfort

The synthetic climatic sequence has an annual mean outdoor temperature of 11.8°C, an occupied-hour mean of 15.0°C, and annual horizontal irradiation of 1,005 kWh m⁻². The grid-carbon series has a mean of 0.205 kgCO₂e kWh⁻¹ and a 5th–95th percentile range of 0.102–0.318 kgCO₂e kWh⁻¹. The normalized control signal is clipped between these percentiles.

The moisture model uses an occupant moisture generation rate of 2.245×10⁻⁵ kg s⁻¹ person⁻¹, equivalent to 55 W person⁻¹ at a latent heat of vaporization of 2,450 kJ kg⁻¹, a cooling-coil bypass factor of 0.10, and no active humidification. PMV/PPD uses 1.2 met, seasonal clothing of 0.5–1.0 clo, simulated RH, and air speeds of 0.10 and 0.15 m s⁻¹ in heating and cooling modes, respectively.

4.3. Numerical truth model, reduced twin, and EKF

The numerical truth model contains 32 primary states and six auxiliary humidity-ratio states. The reduced physics model contains the 21-component augmented vector in Equation (1) and six auxiliary humidity-ratio states. The truth model uses a three-node thermal network per floor, six first-order HVAC states, and a two-node structural thermal model; the reduced model uses a two-node thermal network, algebraic HVAC response, and a one-node structural thermal model (Table 5).

The baseline $P_0$, $Q_k$, and $R_k$ matrices are diagonal. Variances are obtained by squaring the standard deviations in Supplementary Tables S12 and S14 and are scaled to their hourly, daily, or annual update intervals. Electrical-power measurement variance is updated as ${\left[0.01P_{el,k}\right]}^2$.

4.4. Structural reliability, surrogate modeling, and optimization

At the beginning of each 20-year trajectory, R₀, λD, Ec, αT, χr, θS, mUA, and mC are sampled once and retained throughout the trajectory. Annual weather, occupancy, and equipment histories are regenerated independently between years conditional on those trajectory-level factors while preserving prescribed subannual dependence. The annual non-thermal service-load history is sampled within the structural variable set; in the conditional-risk calculation it is not duplicated in the operational trajectory set.

The FORM analysis uses a deterministic common-random-number mapping from the reduced 12-dimensional reliability vector to the annual time histories, permitting stable central finite-difference gradients while maintaining consistency with the separately defined full-uncertainty operational assessment. The reduction of the explicit reliability coordinates and the associated variable-screening checks are detailed in Supplementary Section S7. Importance sampling uses the mixture proposal reported in Supplementary Section S7 and is accepted when CoV does not exceed 0.10.

The uncertainty-aware surrogate-development and untouched-test dataset contains 1,500 complete-model policies. The nominal-operational response engines use a separate 750-policy dataset. The independent pilot instability ensemble produces 1,024 unique anchor policies and 768 matched regions. Across the seven final ablation frameworks, 26,010 raw archive entries reduce to 19,854 within-framework unique policies and 7,946 unique complete-model evaluations after cross-framework caching. Duplicate vectors are identified using a normalized decision-space tolerance of 10⁻⁶. Full framework-level accounting, convergence statistics, and seed separation are reported in Supplementary Tables S20A–S20G and S23B.

P0-N fixes climatic, occupancy, equipment, HVAC, and grid-carbon inputs at their nominal values while retaining probabilistic material, resistance, restraint, degradation, service-load, and structural model variables. It is therefore a nominal-input reliability-aware evaluation rather than a fully deterministic structural analysis. P0-U applies the full modeled uncertainty set to the identical decision vector (Table 6).

5. Results and Discussion

5.1. Numerical convergence, surrogate accuracy, and state updating

Independent convergence against a separate 10,000-realization reference bank gave mean absolute relative errors of 0.9 ± 0.6% for carbon and 1.0 ± 0.6% for discomfort at N = 1,000. The P4 importance-sampling estimate stabilized at 0.0123 with CoV = 0.040 at 50,000 samples (Figure 3).

NMAE and NRMSE are normalized by the independent-test 5th–95th percentile response range (Table 7). All retained nondominated policies and the final robust policy were re-evaluated using the complete coupled model.

EKF updating reduced RMSE by 62.7–69.1% for directly assimilated operational and structural variables. Air-temperature, CO₂, electrical-power, and structural-strain NIS coverages were 94.6–95.4%. Ljung–Box tests over lags 1–24 did not reject innovation whiteness at the 5% level. Operative temperature was a derived output; its RMSE decreased from 0.81 to 0.29°C, but no innovation or NIS statistic was assigned.

5.2. Pareto structure and the effect of uncertainty

The Pareto set demonstrated a genuine three-way trade-off (Figure 4). Low-carbon policies used wider thermal deadbands and weaker low-carbon overventilation while continuing to satisfy the minimum outdoor-air and occupied CO₂ constraints. Comfort- and reliability-priority policies used narrower thermal deadbands and stronger carbon-aware ventilation, increasing heat-pump and fan electricity.

For the identical P0 vector, carbon increased from 42.6 to 50.4 kgCO₂e m⁻² year⁻¹, discomfort increased from 115 to 128 K h year⁻¹, and P_SLS increased from 0.021 to 0.031 (Figure 5). These changes correspond to increases of 18.4%, 11.3%, and 47.6% relative to the nominal-input results, or nominal-input underestimations of 15.5%, 10.2%, and 32.3% relative to the uncertainty-aware values.

5.3. Representative policy performance and IAQ

EUI is reported in kWh m⁻² year⁻¹, weighted carbon intensity in kgCO₂e kWh⁻¹, carbon in kgCO₂e m⁻² year⁻¹, and discomfort in occupied K h year⁻¹. The P0 policy-region decision-instability value is reported under P0-U because it is estimated from stochastic cross-run re-optimization rather than from the single P0-N evaluation. Main-text serviceability probabilities are FORM estimates used during optimization and policy comparison; Supplementary Table S19 reports independent design-point importance-sampling verification.

P1 achieved the lowest carbon among the optimized stochastic extremes but produced the highest discomfort, serviceability risk, and instability (Tables 8a and 8b). P2 minimized discomfort, and P3 minimized P_SLS. P4 did not minimize any single objective; it provided the lowest balanced robust score.

The ppm-consistent CO₂ mass balance produced occupied concentrations that were physically compatible with the continuously enforced minimum outdoor-air requirement. Table 9 pools occupied 15-min zone-time observations across the six zones and 1,000 stochastic realizations; each occupied zone and interval receives equal weight. The reported mean and 99th percentile are therefore zone-time rather than occupant-exposure-weighted statistics. Supplementary Table S8C confirms that occupant weighting preserves the policy ranking. Carbon reduction arose from temporal shifting of ventilation and HVAC demand rather than from relaxation of IAQ requirements.

5.4. Objective instability, decision migration, and robust selection

Objective instability was largest for the low-carbon policy P1 (I_o = 0.61) and the uncertainty-aware P0 region (I_o = 0.52). P2 and P3 satisfied the 0.30 threshold but remained biased toward comfort and reliability, respectively. P4 provided the lowest balanced representative-policy value, I_o = 0.22 (Figure 6). The screening threshold I_o ≤ 0.30 is a benchmark-specific decision rule tied to the fixed normalization ranges, uncertainty model, and operational decision domain; it is not a universal physical or regulatory limit.

The cross-run decision-instability index decreased from 1.08 for the P0 policy region to 0.37 for P4, a 65.7% reduction. The reduction reflects smaller normalized migration in heating setpoint, cooling setpoint, high-carbon CO₂ target, and ventilation gain (Figure 7).

The reference run contained 146 complete-model objective-nondominated candidates before final screening. Fifty-eight satisfied I_o ≤ 0.30; six additional candidates failed the occupied CO₂ criterion, and eight failed equipment or airflow constraints, leaving 44 fully feasible policies (Figure 8). P4 minimized the robust-ranking score, with J_rob = 0.190. It remained top ranked for all tested weight scenarios and feasible instability thresholds from 0.25 to 0.40; no policy satisfied I_o ≤ 0.20.

Paired two-sided centered bootstrap tests used 10,000 resamples with the +1 correction. Relative to P0-U, P4 changed carbon by −1.601 kgCO₂e m⁻² year⁻¹ (95% CI −1.88 to −1.32; p = 0.0002), discomfort by −30.0 K h year⁻¹ (−31.4 to −28.6; p = 0.0002), P_SLS by −0.019 (−0.021 to −0.017; p = 0.0004), and I_d by −0.71 (−0.78 to −0.64; p = 0.0002). For the probability comparison, common importance-sampling indices are resampled with replacement and the paired weighted contributions (Zn, P4, Zn, P0-U) remain together, preserving the proposal draw and likelihood-ratio pairing rather than resampling unweighted Bernoulli indicators.

5.5. Sensitivity and ablation results

Grid-carbon intensity was the leading carbon factor, followed by outdoor temperature, equipment, and occupancy (Figure 9). Discomfort was governed primarily by outdoor temperature, occupancy, envelope conductance, and thermal capacitance. Reliability was dominated by degradation rate, initial resistance, restraint, thermal expansion, outdoor temperature, and non-thermal service load. The gaps between first- and total-order indices demonstrated material interaction effects. Sobol and Shapley intervals combine estimator uncertainty with surrogate uncertainty from 200 posterior Gaussian-process function draws. The leading rankings were additionally checked with a 64-run complete-model three-level D-optimal screening design; the resulting normalized main-effect statistics are reported only as ranking diagnostics and are not interpreted as independent total-order Sobol estimates. Correlation-aware sensitivity uses the Shapley-effect variance-allocation framework of Song, Nelson and Staum (2016) for ranking comparison rather than as a numerical equivalent of independent-input Sobol indices.

Carbon, discomfort, P_SLS, the active instability metric, and I_d are median selected-policy values across 30 paired final runs; hypervolume is reported as run-level mean ± standard deviation (Table 10). The B4 instability entry is IoCT, whereas the other instability entries are the full three-component I_o. Hypervolume for every framework is calculated after complete-model post-evaluation in the common normalized three-objective space. Paired Wilcoxon signed-rank tests are applied only to the 30 paired run-level hypervolume values, with Holm correction across six comparisons. Framework-specific convergence and archive accounting are reported in the Supplementary Material. These are matched framework comparisons rather than an exact factorial decomposition because each formulation generates its own Pareto archive.

5.6. Scope and limitations

This study is a controlled numerical demonstration rather than a field validation. The six-story office building is represented by six equivalent floor zones and does not resolve room-level airflow, interzonal contaminant transport, occupant adaptation, local draught, or heterogeneous terminal-unit behavior. The climatic and grid-carbon sequences are synthetic benchmark inputs; the grid signal represents average operational carbon intensity rather than marginal emissions. The primary uncertainty analyses are conditional on the adopted marginal distributions and dependence assumptions, while the correlated-input Shapley analysis is a dependence stress test rather than a calibrated building-specific dependence model. Weather, occupancy, schedules, and sampled process multipliers are treated as known exogenous inputs during synthetic assimilation. The primary discomfort and CO₂ summaries are occupied zone-time metrics rather than occupant-exposure-weighted metrics. The structural objective is an equivalent component-level serviceability stress-margin benchmark and is not a code-compliance assessment of cracking, deflection, fatigue, ultimate strength, or whole-building safety. Synthetic observations cannot reproduce all sensor faults, commissioning errors, cyber-physical delays, or model-form discrepancies encountered in operation. These limitations restrict external generalization; field validation and building-specific calibration are required before deployment. The I_o screening threshold and all instability values are benchmark-specific and require recalibration before transfer to another building, policy space, or uncertainty model.

6. Conclusions and Future Research Directions

This study developed and numerically demonstrated an instability-aware digital twin for robust co-optimization of operational carbon emissions, occupied thermal comfort, and structural serviceability reliability under uncertainty. The framework integrates multi-rate thermal, moisture, energy, IAQ, structural, data-assimilation, reliability, surrogate, and evolutionary-optimization components within one decision architecture.

Re-evaluation of the identical P0 decision vector under the full modeled uncertainty set increased carbon from 42.6 to 50.4 kgCO₂e m⁻² year⁻¹, occupied zone-time discomfort from 115 to 128 K h year⁻¹, and P_SLS from 0.021 to 0.031. The selected P4 policy used T_h,set = 20.8°C, T_c,set = 25.4°C, C_high = 925 ppm, and k_γ = 0.75. It produced 48.8 kgCO₂e m⁻² year⁻¹, 98 K h year⁻¹, P_SLS = 0.012, β = 2.257, I_o = 0.22, and I_d = 0.37. Relative to P0-U, P4 reduced carbon by 3.2%, discomfort by 23.4%, P_SLS by 61.3%, and decision instability by 65.7%.

The IAQ model confirmed that P4 maintained a mean occupied CO₂ concentration of 684 ppm, a 99th percentile of 821 ppm, and zero occupied zone-hours above 1,000 ppm. The operational benefit was achieved by carbon-aware timing of ventilation and HVAC demand rather than by weakening ventilation requirements.

The matched ablation results show that reliability coupling, uncertainty propagation, objective-instability screening, and decision-instability evaluation change the retained policy region and selected policy. The independent pilot ensemble and frozen instability surrogates remove circular dependence between I_d estimation and the final paired runs. B5 and B6 provide complementary framework-level controls, but they do not constitute an exact factorial decomposition because each formulation generates its own Pareto archive. The proposed framework achieved the highest hypervolume and the lowest combined instability among the compared formulations.

Future research should extend the framework through field deployment, high-resolution multizone airflow and contaminant models, nonlinear finite-element structural models, measured degradation and inspection data, correlated nonstationary uncertainty, lifecycle cost and embodied carbon, renewable energy and storage, and safety-filtered learning-based control. Field implementation remains necessary before the framework can support building-specific operational decisions.