Ensemble Simulations Research Articles

North Atlantic Subtropical Mode Water properties: intrinsic and atmospherically forced interannual variability

Abstract. This study investigates the contributions of the ocean's chaotic intrinsic variability (CIV) and atmospherically forced variability to the interannual fluctuations in the North Atlantic Subtropical Mode Water (STMW) properties. Utilizing a 0.25° regional 50-member ocean–sea-ice ensemble simulation driven by an original surface forcing method and perturbed initially, the forced variability in STMW properties is estimated from ensemble mean fluctuations, while CIV is determined from deviations around the ensemble mean within each member. The model successfully captures the main features of STMW, showing correct agreement with observation-based ARMOR3D data in terms of location, seasonality, mean temperature and volume, as well as interannual variance of STMW's main properties. CIV significantly impacts STMW, explaining 10 %–13 % and 28 %–44 % of the interannual variance of its geometric and thermohaline mean properties, respectively, with a maximum imprint on STMW temperature. Observation-based and simulated intrinsic-to-total variance ratios are mostly consistent, dispelling concerns about a signal-to-noise paradox. This study also illustrates the advantages of ensemble simulations over single simulations in understanding oceanic fluctuations and attributing them to external drivers, while also cautioning against overreliance on individual simulation assessments.

Contrasting Domestic and Global Impacts of Emission Reductions in China on Tropospheric Ozone

AbstractSustained near‐surface ozone pollution despite stringent emission controls in China has drawn wide attention. However, how such rapid regional emission reductions affect global tropospheric ozone remains unexplored. Here, using an ensemble of measurements and global model simulations, we find Chinese emission reductions during 2013–2020 have increased domestic surface ozone; particularly in North China, they have decreased the global tropospheric ozone burden by 2.7 Tg (54% of China's anthropogenic contribution and 7% of global anthropogenic contribution in 2019) and ozone radiative forcing by 8.8 mW m−2 (64% of China's anthropogenic contribution and 9% of global anthropogenic contribution in 2019). Nitrogen oxides (NOx) associated nonlinear chemistry and height‐dependent atmospheric transport patterns primarily cause the contrasting impacts. Future Chinese NOx emission controls would continue to decrease global tropospheric ozone, while regional surface ozone may still increase unless NOx reductions exceed ∼53% of the 2019 level. Our results reveal the substantial benefit of regional emission reductions in China on global tropospheric ozone mitigation.

Quantum Simulations of Radiation Damage in a Molecular Polyethylene Analog.

An atomic-level understanding of radiation-induced damage in simple polymers like polyethylene is essential for determining how these chemical changes can alter the physical and mechanical properties of important technological materials such as plastics. Ensembles of quantum simulations of radiation damage in a polyethylene analog are performed using the Density Functional Tight Binding method to help bind its radiolysis and subsequent degradation as a function of radiation dose. Chemical degradation products are categorized with a graph theory approach, and occurrence rates of unsaturated carbon bond formation, crosslinking, cycle formation, chain scission reactions, and out-gassing products are computed. Statistical correlations between product pairs show significant correlations between chain scission reactions, unsaturated carbon bond formation, and out-gassing products, though these correlations decrease with increasing atom recoil energy. The results present relatively simple chemical descriptors as possible indications of network rearrangements in the middle range of excitation energies. Ultimately, the work provides a computational framework for determining the coupling between nonequilibrium chemistry in polymers and potential changes to macro-scale properties that can aid in the interpretation of future radiation damage experiments on plasticmaterials.

Constrained cosmological simulations of the Local Group using Bayesian hierarchical field-level inference

We present a novel approach based on Bayesian field-level inference that provides representative initial conditions for simulation of the Local Group (LG) of galaxies and its neighbourhood, constrained by present-day observations. We extended the Bayesian Origin Reconstruction from Galaxies ( algorithm with a multi-resolution approach, allowing us to reach the smaller scales needed to apply the constraints. Our data model simultaneously accounts for observations of mass tracers within the dark haloes of the Milky Way (MW) and M31, for their observed separation and relative velocity, and for the quiet surrounding Hubble flow, represented by the positions and velocities of 31 galaxies at distances between one and four megaparsec. Our approach delivers representative posterior samples of realisations that are statistically and simultaneously consistent with all of these observations, leading to significantly tighter mass constraints than found if the individual datasets are considered separately. In particular, we estimate the virial masses of the MW and M31 to be $ (M_ 200c /M_ and $12.33 respectively, their sum to be $ 200c /M_ and the enclosed mass within spheres of radius $R$ to be $ (M(R)/M_ and $12.96 for $R=1$ Mpc and $3$ Mpc, respectively. The M31-MW orbit is nearly radial for most of our realisations, and most of them feature a dark matter sheet aligning approximately with the supergalactic plane, despite the surrounding density field not being used explicitly as a constraint. High-resolution, high-fidelity resimulations from initial conditions identified using the approximate simulations of our inference scheme continue to satisfy the observational constraints, demonstrating a route to future high-resolution, full-physics simulations of ensembles of LG look-alikes, all of which closely mirror the observed properties of the real system and its immediate environment.

Nonlinear changes in urban heat island intensity, urban breeze intensity, and urban air pollutant concentration with roof albedo

This study systematically examines how the urban heat island (UHI) and urban breeze circulation (UBC) respond to an increase in roof albedo (αr) and its influence on urban air pollutant dispersion. For this, idealized ensemble simulations are performed using the Weather Research and Forecasting (WRF) model. The increase in αr from 0.20 to 0.65 decreases the UHI intensity, UBC intensity, and urban planetary boundary layer (PBL) height in the daytime (from 1200 to 1700 LST) by 47%, 36%, and 6%, respectively. As both UBC intensity and urban PBL height decrease, the daytime urban near-surface passive tracer concentration increases by 115%. The daytime UHI intensity, UBC intensity, and urban tracer concentration nonlinearly change with αr: For 0.10 ≤ αr < 0.80, the rates of changes in the UHI intensity, UBC intensity, and urban tracer concentration with αr overall increase as αr increases. For αr ≥ 0.80, the daytime roof surface temperature is notably lower than the daytime urban near-surface air temperature, the UHI intensity, UBC intensity, and urban tracer concentration very slightly changing with αr. This study provides insights into the associations between changes in roof surface temperature and roof surface energy fluxes with αr and those in UHI intensity.

A computational framework to support probabilistic criticality modelling for the geological disposal of radioactive waste

Nuclear Waste Services is tasked with disposal of the UK’s higher-activity radioactive waste in a Geological Disposal Facility. The disposal of fissile nuclides requires a demonstration that there is no significant concern from criticality, i.e. a fission chain reaction. While waste packages will initially be emplaced in a subcritical configuration, over the long timescales following closure there is potential for waste packages to degrade and for nuclides to be dispersed in the subsurface by groundwater, leading to the potential for a critical system forming. To facilitate modelling, a codebase has been developed which interfaces a probabilistic simulation tool (GoldSim) with a neutron transport code (MONK/MCNP). This allows large ensemble simulations to be run iteratively to determine limiting fissile masses which satisfy a criticality safety criterion. This paper documents the main algorithms and methodologies implemented within this framework, and provides background and example results illustrating the application to post-closure criticality modelling.

Computationally efficient and error aware surrogate construction for numerical solutions of subsurface flow through porous media

Limiting the injection rate to restrict the pressure below a threshold at a critical location can be an important goal of simulations that model the subsurface pressure between injection and extraction wells. The pressure is approximated by the solution of Darcy’s partial differential equation for a given permeability field. The subsurface permeability is modeled as a random field since it is known only up to statistical properties. This induces uncertainty in the computed pressure. Solving the partial differential equation for an ensemble of random permeability simulations enables estimating a probability distribution for the pressure at the critical location. These simulations are computationally expensive, and practitioners often need rapid online guidance for real-time pressure management. An ensemble of numerical partial differential equation solutions is used to construct a Gaussian process regression model that can quickly predict the pressure at the critical location as a function of the extraction rate and permeability realization. The Gaussian process surrogate analyzes the ensemble of numerical pressure solutions at the critical location as noisy observations of the true pressure solution, enabling robust inference using the conditional Gaussian process distribution.Our first novel contribution is to identify a sampling methodology for the random environment and matching kernel technology for which fitting the Gaussian process regression model scales as O(nlogn) instead of the typical O(n3) rate in the number of samples n used to fit the surrogate. The surrogate model allows almost instantaneous predictions for the pressure at the critical location as a function of the extraction rate and permeability realization. Our second contribution is a novel algorithm to calibrate the uncertainty in the surrogate model to the discrepancy between the true pressure solution of Darcy’s equation and the numerical solution. Although our method is derived for building a surrogate for the solution of Darcy’s equation with a random permeability field, the framework broadly applies to solutions of other partial differential equations with random coefficients.

A framework for visual comparison of scalar fields with uncertainty

AbstractScientists working with uncertain data, such as climate simulations, medical images, or ensembles of physical simulations, regularly confront the problem of comparing observations, e.g., to identify similarities, differences, or patterns. Current approaches in comparative visualization of uncertain scalar fields mainly rely on juxtaposition of both data and uncertainties, where each is represented using, e.g., color mapping or volume rendering. While interpretation of uncertain scalar data from visual encodings is already cognitively challenging, comparison of uncertain fields without explicit visualization support adds a further layer of complexity. In this paper, we present a theoretical framework to devise and describe a class of techniques that directly visualize differences between two or more uncertain scalar fields in a single image. We model each such technique as a combination of one or more interpolation stages, with the application of distance measures on random variables to the resulting distributions, and an appropriate visual encoding. Our framework captures existing methods and lends itself well to formulating new comparative visualization techniques for uncertain data for different visualization scenarios. Furthermore, by modeling uncertain scalar field differences as random variables themselves, we enable additional opportunities for comparison. We demonstrate the usefulness of our framework and its properties by applying it to effective comparative visualization techniques for several synthetic and real-world data sets.

Symmetry-adapted Markov state models of closing, opening, and desensitizing in α 7 nicotinic acetylcholine receptors.

α7 nicotinic acetylcholine receptors (nAChRs) are homopentameric ligand-gated ion channels with critical roles in the nervous system. Recent studies have resolved and functionally annotated closed, open, and desensitized states of these receptors, providing insight into ion permeation and lipid binding. However, the process by which α7 nAChRs transition between states remains unclear. To understand gating and lipid modulation, we generated two ensembles of molecular dynamics simulations of apo α7 nAChRs, with or without cholesterol. Using symmetry-adapted Markov state modeling, we developed a five-state gating model. Free energies recapitulated functional behavior, with the closed state dominating in absence of agonist. Open-to-nonconducting transition rates corresponded to experimental open durations. Cholesterol relatively stabilized the desensitized state, and reduced open-desensitized barriers. These results establish plausible asymmetric transition pathways between states, define lipid modulation effects on the α7 nAChR conformational cycle, and provide an ensemble of structural models applicable to rational design of lipidic pharmaceuticals.

Exploration of Cryptic Pockets Using Enhanced Sampling Along Normal Modes: A Case Study of KRAS G12D.

Identification of cryptic pockets has the potential to open new therapeutic opportunities by discovering ligand binding sites that remain hidden in static apo structures of a target protein. Moreover, allosteric cryptic pockets can become valuable for designing target-selective ligands when the natural ligand binding sites are conserved in variants of a protein. For example, before an allosteric cryptic pocket was discovered, KRAS was considered undruggable due to its smooth surface and conservation of the GDP/GTP binding pocket across the wild type and oncogenic isoforms. Recent identification of the Switch-II cryptic pocket in the KRASG12C mutant and FDA approval of anticancer drugs targeting this site underscores the importance of cryptic pockets in solving pharmaceutical challenges. Here, we present a newly developed approach for the exploration of cryptic pockets using weighted ensemble molecular dynamics simulations with inherent normal modes as progress coordinates applied to the wild type KRAS and the G12D mutant. We performed extensive all-atomic simulations (>400 μs) with and without several cosolvents (xenon, ethanol, benzene), and analyzed trajectories using three distinct methods to search for potential binding pockets. These methods have been applied as a proof-of-concept to KRAS and have shown they can predict known cryptic binding sites. Furthermore, we performed ligand-binding simulations of a known inhibitor (MRTX1133) to shed light on the nature of cryptic pockets in KRASG12D and the role of conformational selection vs induced-fit mechanism in the formation of these cryptic pockets.

Exploring the potential of forest snow modeling at the tree and snowpack layer scale

Abstract. Boreal and sub-alpine forests host seasonal snow for multiple months per year; however, snow regimes in these environments are rapidly changing due to rising temperatures and forest disturbances. Accurate prediction of forest snow dynamics, relevant for ecohydrology, biogeochemistry, cryosphere, and climate sciences, requires process-based models. While snow schemes that track the microstructure of individual snow layers have been proposed for avalanche research, so far, tree-scale processes resolving canopy representations only exist in a few snow-hydrological models. A framework that enables layer- and microstructure-resolving forest snow simulations at the meter scale is lacking to date. To fill this research gap, this study introduces the forest snow modeling framework FSMCRO, which combines two detailed, state-of-the art model components: the canopy representation from the Flexible Snow Model (FSM2) and the snowpack representation of the Crocus ensemble model system (ESCROC). We apply FSMCRO to discontinuous forests at boreal and sub-alpine sites to showcase how tree-scale forest snow processes affect layer-scale snowpack properties. Simulations at contrasting locations reveal marked differences in stratigraphy throughout the winter. These arise due to different prevailing processes at under-canopy versus gap locations and due to variability in snow metamorphism dictated by a spatially variable snowpack energy balance. Ensemble simulations allow us to assess the robustness and uncertainties of simulated stratigraphy. Spatially explicit simulations unravel the dependencies of snowpack properties on canopy structure at a previously unfeasible level of detail. Our findings thus demonstrate how hyper-resolution forest snow simulations can complement observational approaches to improve our understanding of forest snow dynamics, highlighting the potential of such models as research tools in interdisciplinary studies.

Bayesian estimates for changes of the Russian river runoff in the 21st century as based on the CMIP6 model ensemble simulations

Based on ensemble calculations with the CMIP6 (Coupled Model Intercomparison Project, phase 6) climate models and using Bayesian averaging, an analysis was conducted on the changes in the 21st century runoff of several Russian rivers – the Volga, Ob, Yenisei, Lena, Amur, and Selenga. Bayesian weights considered the quality of models’ reproduction of runoff (long-term average runoff, linear runoff trend over the time interval with available runoff observations, interannual and interdecadal variability). The quality of runoff characteristics reproduction by individual models in the CMIP6 ensemble varies most significantly for the long-term average runoff, runoff trend, and, to a lesser extent, for interannual variability. In the 21st century, the ensemble average runoff increases for most of the analyzed rivers, except for the Volga. This increase is more pronounced under scenarios with larger anthropogenic impacts. It is especially significant for the SSP5-8.5 scenario (Shared Socioeconomic Pathways, 5-8.5), under which the runoff increase trend from 2015 to 2100 relative to its current long-term average is (10 ± 4)% for the Ob, (16 ± 3)% for the Yenisei, (39 ± 7)% for the Lena, (36 ± 7)% for the Amur, and (18 ± 6)% for the Selenga. The primary reason for the change in ensemble average runoff in the 21st century in models under all SSP scenarios is the change in precipitation. Accounting for differences in model quality in reproducing river runoff on average for 2015–2100 reduces inter-model deviations relative to the corresponding values with uniform weighting of model results by 6–26%, depending on the SSP scenario and river basin.

Anthropogenic influence on precipitation in Aotearoa New Zealand with differing circulation types

This study quantifies the influences of anthropogenic forcing to date on precipitation over Aotearoa New Zealand (ANZ). Large ensembles of simulations from the weather@home regional climate model experiments are analysed under two scenarios, a natural (NAT) or counter-factual scenario which excludes human-induced changes to the climate system and an anthropogenic (ANT) or factual scenario. The impacts of anthropogenic forcing on precipitation are analysed in the context of large-scale circulation types characterized using an existing Self Organizing Map classification. The combined effect of both thermodynamics and dynamics are compared with values expected from the Clausius–Clapeyron (C–C) relation. Changes in the precipitation intensity due to greenhouse gas-forced temperature rise are lower than the expected C–C value. However extreme precipitation changes approach the C–C value for some circulation types. Specifically westerly flows enhance precipitation change across ANZ relative to the C–C rate, particularly over the West Coast. Conversely, northwesterly flows reduce the change over the North Island relative to the C–C value. Moreover, the wet day frequency generally reduces in the ANT scenario relative to NAT, reductions are largest on the West Coast of the South Island for westerly flows. Additionally, the frequency of days with extreme precipitation rises over ANZ for most circulation patterns, except in Northland and for northwesterly flows. This underscores the combined influence of dynamics and thermodynamics in shaping both precipitation intensity and frequency patterns across ANZ.

Atomistic simulations of polysaccharide materials for insights into their crystal structure, nanostructure, and dissolution mechanism

AbstractCrystalline polysaccharides are abundant in nature and can be transformed into highly functional materials. However, the molecular basis for the formation of higher-order structures remains unclear. Computer simulation is an advanced tool for modeling macromolecular structures, and the atomistic simulations provide valuable information on the crystalline polysaccharides. Fiber deformation, crystalline transition, and novel nanostructures of cellulose were characterized through molecular dynamics simulations and density functional theory calculations of models of molecular chain sheets extracted from the crystal structure of the cellulose polymorphs. Extended ensemble molecular dynamics simulations were applied to analyze the artificial crystal structure of non-natural amylose analog polysaccharides, revealing the hexagonal packing of double helices through the self-assembly of molecular chains dispersed in aqueous solution. Dissolution simulations of the cellulose and chitin crystalline fibers revealed that the anions of ionic liquids, with their solvation power, played a key role in the cleavage of intermolecular hydrogen bonds in the crystal structure, whereas the cations contributed to irreversible molecular chain dispersion. The good correlation between the actual solubility of polysaccharides and the predicted number of intermolecular hydrogen bonds prompted the development of a platform that combined simulations and machine learning for high-throughput screening of solvents for cellulose and chitin.

Fast generation of high-dimensional spatial extremes

Widespread extreme climate events cause many fatalities, economic losses and have a huge impact on critical infrastructure. It is therefore of utmost importance to estimate the frequency and associated consequences of spatially concurrent extremes. Impact studies of climate extremes are severely hampered by the lack of extreme observations, and even large ensembles of climate simulations often do not include enough extreme or record-breaking climate events for robust analysis. On the other hand, weather generators specifically fitted to extreme observations can quickly generate many physically or statistically plausible extreme events, even with intensities that have never been observed before. We propose a Fourier-based algorithm for generating high-resolution synthetic datasets of rare events, using essential concepts of classical modelling of (spatial) extremes. Here, the key feature is that the stochastically generated datasets have the same spatial dependence as the observed extreme events. Using high-resolution gridded precipitation and temperature datasets, we show that the new algorithm produces realistic spatial patterns, and is particularly attractive compared to other existing methods for spatial extremes. It is exceptionally fast, easy to implement, scalable to high dimensions and, in principle, applicable for any spatial resolution. We generated datasets with 10000 gridpoints, a number that can be increased without difficulty. Since current impact models often require high-resolution climate inputs, the new algorithm is particularly useful for improved impact and vulnerability assessment.

Bridging Hydrological Ensemble Simulation and Learning Using Deep Neural Operators

AbstractEnsemble‐based simulation and learning (ESnL) has long been used in hydrology for parameter inference, but computational demands of process‐based ESnL can be quite high. To address this issue, we propose a deep neural operator learning approach. Neural operators are generic machine learning algorithms that can learn functional mappings between infinite‐dimensional spaces, providing a highly flexible tool for scientific machine learning. Our approach is built upon DeepONet, a specific deep neural operator, and is designed to address several common problems in hydrology, namely, model parameter estimation, prediction at ungaged locations, and uncertainty quantification. Here we demonstrate the effectiveness of our DeepONet‐based workflow using an existing large model ensemble created for an eastern U.S. watershed that is instrumented with 10 streamflow gages. Results suggest DeepONet achieves high efficiency in learning an ML surrogate model from the model ensemble, with the modified Kling‐Gupta Efficiency exceeding 0.9 on holdout test sets. Parameter inference, carried out using the trained DeepONet surrogate model and genetic algorithm, also yields robust results. Additionally, we formulate and train a separate DeepONet model for physics‐informed, seq‐to‐seq streamflow forecasting, which further reduces biases in the pre‐trained DeepONet surrogate model. While this study focuses primarily on a single watershed, our approach is general and may be extended to enable learning from model ensembles across multiple basins or models. Thus, this research represents a significant contribution to the application of hybrid machine learning in hydrology.

Semi-empirical supported, Ab initio derived thermodynamic properties for ClO2 and its sub and extended species, applied in water treatment cycles

This paper describes a group of sixty (60) sub and extended chlorine oxide species with the general formulae of ClxOy (with x ≤ 2, y ≤ 8). Their role in water treatment cycles, behaving as key reactive species, is represented by a complex sequence of chemical inter-dependencies, exposed as a cohesive set of chemical reactions to demonstrate their cyclic role in aqueous media. An empirical/semi-empirical computational approach, supported by Ab Initio simulations, in accordance with open-shell character, has been followed to determine their optimum molecular geometries, to obtain their thermochemical properties. Besides a single molecular analysis, Grand Canonical Ensemble simulations, supported by a revised library of force field parameters, constituted a core component of the computational approach and proved to be invaluable in confirming thermochemical properties. This approach also offered finite estimates of optimum model sizes, a benefit with wider modelling application. Extended molecular species of ClO2 display a complex sequence of bonding character, with a variable charge dissipation (reported as partial charges), which complicates selection of basis sets in optimizing molecular geometries, during Ab Initio analyses.Optimum molecular geometries were obtained using Gaussian and MOPAC, which in turn resulted in reliable Heats of Formation. These correlated well with energies extracted from the open literature. Thermodynamic Analysis of the reaction of selected chlorine oxides with water using FactSage, predicted the production of known and two previously undetected chlorine species, [ClO4]- and [ClOH2]+.

Initiation and evolution of urban-induced precipitation under different background wind speeds: Roles of urban breeze circulation and cold pool

Urban-enhanced precipitation is often observed far downwind of cities, suggesting important roles of background wind. This study examines the impacts of background wind on the initiation and evolution of urban-induced precipitation through idealized ensemble large-eddy simulations with different initial background wind speeds of 0, 1, 2, 3, 4, and 5 m s−1. Without background wind, a strong urban breeze circulation (UBC) initiates a relatively large number of cumulus clouds in the urban area, which exhibit greater cloud top heights than those in the rural area. The early-initiated clouds do not produce precipitation but moisten the middle troposphere, assisting later-initiated clouds to develop into deep convective clouds and produce concentrated precipitation in the urban area. As the background wind speed increases, the UBC weakens. Also, the high-humidity area generated by early-initiated clouds is continuously advected away from the urban area. Consequently, the precipitation initiation is delayed as the background wind speed increases. Precipitating clouds are advected downwind and generate a cold pool, and updrafts are repeatedly produced at the downwind boundary of the cold pool. The downwind-advected clouds are fed by the newly produced updrafts, being a long-lived organized precipitation system that produces a considerable amount of precipitation in the downwind area. This feeding effect is most prominent for a moderate background wind (3 m s−1), causing the precipitation amount and duration to be maximized. For stronger background winds, precipitation is initiated by the mechanical lifting of the warm air of the urban boundary layer over the cool rural boundary layer around sunset.

The Dual Nature of Entrainment-Mixing Signatures Revealed Through Large-Eddy Simulations of a Convection-Cloud Chamber

Abstract Entrainment of subsaturated air into a cloud can influence its optical and microphysical properties in various ways, depending on the droplet evaporation and turbulent mixing time scales. Previous experiments in the Pi Convection-Cloud Chamber have revealed that, given a fixed entrained air property, the mixing of entrained subsaturated air results in complete evaporation of some cloud droplets, with the rest remaining unchanged. This is a signature of inhomogeneous mixing. While comparing the results of entrainment with varying air properties, the mixing signature appears as if the subsaturated air is well mixed with the cloud to evenly reduce the droplets’ size. In other words, taken together, the experiments appear to have the signature of homogeneous mixing. To explore these results in a greater depth, we conduct large-eddy simulations combined with a bin microphysics scheme. Our results reproduce the similar signatures of inhomogeneous and homogeneous mixing, implying that LES can resolve the inhomogeneous mixing when the grid spacing is smaller than the entrained air parcel. Additionally, we observe that increasing aerosol injection rate enhances the signature of inhomogeneous mixing, while coarser grid spacing diminishes it. Finally, the change in wall fluxes in response to various entrained air properties confirms that the homogeneous signature seen in the analysis of an ensemble of simulations is the result of various equilibrium states. This further strengthens the suggestion that the homogeneous mixing signature found in aircraft observations near the cloud top may result from combining entrainment events of different intensities, possibly caused by various-sized eddies.

Statistical Thermodynamic Description of Self-Assembly of Large Inclusions in Biological Membranes.

Recent studies revealed anomalous underscreening in concentrated electrolytes, and we suggest that the underscreened electrostatic forces between membrane proteins play a significant role in the process of self-assembly. In this work, we assumed that the underscreened electrostatic forces compete with the thermodynamic Casimir forces induced by concentration fluctuations in the lipid bilayer, and developed a simplified model for a binary mixture of oppositely charged membrane proteins with different preference to liquid-ordered and liquid-disordered domains in the membrane. In the model, like macromolecules interact with short-range Casimir attraction and long-range electrostatic repulsion, and the cross-interaction is of the opposite sign. We determine energetically favored patterns in a system in equilibrium with a bulk reservoir of the macromolecules. Different patterns consisting of clusters and stripes of the two components and of vacancies are energetically favorable for different values of the chemical potentials. Effects of thermal flutuations at low temperature are studied using Monte Carlo simulations in grand canonical and canonical ensembles. For fixed numbers of the macromolecules, a single two-component cluster with a regular pattern coexists with dispersed small one-component clusters, and the number of small clusters depends on the ratio of the numbers of the molecules of the two components. Our results show that the pattern formation is controlled by the shape of the interactions, the density of the proteins, and the proportion of the components.