Ra incorporation into (Ba,Ra)SO4 solid solutions is a key control on Ra mobility in groundwater systems and typically occurs through coprecipitation and recrystallization. The effectiveness and persistence of these mechanisms under long-term reactive transport conditions remain poorly constrained, particularly in fractured crystalline rocks, where Ra migration is controlled by the coupled effects of heterogeneous flow, advective-diffusive transport, fracture-matrix mass exchange, reaction kinetics, and evolving hydrogeochemical conditions. Here, we employ 3D reactive transport modeling using PFLOTRAN to investigate Ra mobility under near-field conditions relevant to geological nuclear waste repositories. The models couple fluid flow, solute transport, and non-ideal solid solution-aqueous solution (SS-AS) interactions, incorporating a regular Guggenheim solid solution model and composition-dependent dissolution-precipitation kinetics for stoichiometric solid solutions. Simulations are conducted in a 10 m & times; 10 m & times; 10 m fracture-matrix domain upscaled from a discrete fracture network, with a constant-flux inflow boundary and fixed Ra concentration representing a sustained Ra source over 10,000 years. The results show that coprecipitation leads to strong but transient Ra immobilization, with substantial Ra uptake within the first similar to 200 years, followed by progressive Ra remobilization as sulfate is depleted and previously formed solid solutions dissolve. Consequently, Ra retention decreases markedly and becomes minimal after similar to 1000 years. In contrast, recrystallization supports persistent Ra immobilization throughout the entire 10,000-year simulation period, provided that sufficient barite remains available within fracture zones. This mechanism is sustained by kinetically controlled coupled dissolution-reprecipitation and produces a characteristic spatial zonation, with Ra-rich solid solutions near inflow regions and progressively Ba-rich compositions downstream. Sensitivity analyses further demonstrate that increasing Ba/Ra ratios in the inflowing fluid can reduce the long-term Ra retention under kinetically controlled reactive transport conditions, in contrast to predictions based solely on thermodynamic equilibrium. Fractures are identified as the dominant domains for long-term Ra immobilization, whereas the low-permeability matrix contributes only minimally due to limited diffusive accessibility. Once fractures lose their retention capacity, aqueous Ra is predominantly flushed from the system rather than retained within the matrix. Overall, these results suggest that equilibrium assumptions commonly adopted in radionuclide safety assessments are insufficient to predict Ra behavior in complex subsurface systems, and thus robust evaluation of long-term Ra mobility requires coupling reaction mechanisms and kinetics with flow and transport in evolving fracture-matrix systems. Although both coprecipitation and recrystallization form (Ba,Ra)SO4 solid solutions, their distinct microscopic mechanisms lead to contrasting long-term behaviors when upscaled to the field scale. These findings have important implications for nuclear waste disposal and managing Ra contamination in geothermal and mining environments.

The call for an “energy transition” undeniably requires finding low-emission energies, among which serpentinization-sourced natural hydrogen (H2) is one potential candidate. However, the knowledge needed to successfully explore for this new energy source, i.e., to understand when and where natural H2 forms during the Wilson cycle, is yet incomplete. Answering these questions is challenging and requires not only a much-improved understanding of serpentinization processes, but also a deep understanding of the geological systems in which serpentinization occurs. The aim of this contribution is not only to focus on how natural H2 forms, migrates and is trapped, but also to investigate the characteristics and evolution of the geological systems that hosted, or host, serpentinization and potential natural H2 production. More specifically, we examine the evolution of mantle rocks during their emplacement in ocean continent transitions (OCTs) and their later reactivation and emplacement in rift-inversion orogens. To do so, we review the evolution of OCT-derived ophiolites located in the Alpine system. We focus on two settings: the Western Pyrenean fossil OCT exposed in SW-France and the Grischun-Malenco fossil OCTs exposed in SE-Switzerland and northern Italy. Their choice is linked to the fact that they are not only the best documented and exposed examples of OCT-derived ophiolites worldwide, but also to the fact that they bear active H2 systems. Compared with ophiolites derived from Mid Ocean Ridges (MOR) or Supra-Subduction Zones (SSZ), OCT-derived ophiolites are substantially different in terms of mantle rock composition, emplacement mechanisms during collision, and potential for natural H2 production. In this study, we review the characteristics of OCTs based on the well calibrated distal Iberian margin and the Grischun-Malenco fossil OCT. We then synthesize the state-of-the art knowledge on these well described OCTs and OCT-derived ophiolites and their potential for natural H2 systems. Finally, we discuss how far this knowledge can be used to explore for native H2 along the Tethyan suture zones in Eastern Europe, and in collisional orogens in general.

Tc-99 has drawn widespread concern because of its long half-life, high fission yield and high mobility in research of radioactive waste disposal and environmental remediation. TcO4− through-diffusion experiments in Opalinus Clay (OPA) were performed under air and under argon No noticeable Tc-breakthrough was observed for over one year while the total Tc concentration in the source reservoir was steadily decreasing under both air and argon atmosphere. The total Tc activity distribution in the clay sample along the diffusion direction was obtained by slicing the OPA clay samples retrieved from the diffusion cells, using the abrasive peeling technique. In the case of diffusion under air atmosphere, almost no Tc was measured in that part of the sample close to source reservoir, while much more Tc was measured under argon atmosphere. A reasonable explanation for this observation is that the reductive retention of Tc plays a significant role during transport. A reactive transport model was constructed to simulate the diffusion process whereby the diffusion of Tc was coupled with redox reactions. Even though reduction of TcO4− by aqueous Fe2+ is thermodynamically feasible, it was not observed in the experiment. Furthermore, Fe2+ associated with solid phase was demonstrated to be more active than aqueous Fe2+. Therefore, surface complexation redox reaction was proposed. Dissolution rate of pyrite, equilibrium constant and diffusion coefficient of TcO4− were considered as possible factors controlling the redox reaction. Modeling results showed that TcO4− diffused into the clay and was partially reduced into surface complexed Tc(IV) by pyrite. When TcO4− transported under air atmosphere, O2 competitively consumed Fe2+ and pyrite, resulting in no Tc immobilization in the related zone.

The evolution of steel/bentonite interfaces in engineered barrier systems of radioactive waste disposal is governed by the corrosion of the steel and the interaction of corrosion products with the bentonite. In the early post-closure phase of the repository, the transition from aerobic to anaerobic corrosion in combination with evolving temperature conditions and chemical gradients due to the hydration of the bentonite take place. In situ tests in underground laboratories provide a high degree of representativeness for the early transient evolution of the interfaces. Complex Fe-bentonite interaction patterns have been observed, which include a concave Fe accumulation front and distinct coloured halos around the corroding steel. A reactive transport model for the FEBEX in situ test is presented, which describes the transition from aerobic to anaerobic corrosion, the transport of O2 in the gas and liquid phase and the chemical evolution of the steel/bentonite interface during the 18 years of the experiment. The model successfully reproduces the major steps of the interface evolution: an aerobic corrosion phase characterized by accumulation of goethite in the corrosion layer, and an anaerobic corrosion phase, where Fe(II) and Fe(II)/Fe(III) corrosion products form, including a gradual re-dissolution of the aerobic corrosion product. In the anaerobic corrosion phase, some of the Fe(II) diffuses into the bentonite, where it partly reacts with remaining O2 to form Fe(III) precipitates or interacts with the clay minerals via surface complexation or cation exchange. Two different approaches for the implementation of a potential electron transfer from sorbed Fe(II) to structural Fe(III) are presented, but the lack of experimental data does not allow evaluation of their representativeness for the FEBEX in situ experiment. The model calculations qualitatively reproduce the concave shape of the Fe accumulation front in the bentonite and a distinct zonation with an iron-oxide precipitation dominated zone adjacent to the steel and a Fe(II) sorption dominated zone further into the bentonite. Sensitivity cases with respect to the hydration of the bentonite and the transport parameters point to the importance of the fast O2 diffusion in the gas phase for the formation of these characteristic interface features. The calculated thickness of the Fe-bentonite interaction zone varies from 1 to 4 cm, which is smaller than observed in some areas of the experiment, where locally this zone extended > 10 cm into the bentonite. It is possible that this difference is due to a more complex flow pattern in the in situ experiment that is not captured by the model or due to an additional decoupled electron transport across oxide rich layers.

Reactive‐transport models (RTMs) have been an integral part of the safety case and design optimization for the nuclear waste repository in Finland. Such highly resolved, internally complex RTMs tend to be computationally demanding, thus their design is often a compromise between representativeness and computational efficiency. A modeling strategy is presented, demonstrating how a 2D simplified model can be constructed based upon and calibrated against a geometrically realistic 3D RTM. Key geometrical parameters of the 2D model are consistent with the 3D model, and the 2D model is verified by a series of benchmarks for various designs and conditions of the repository. Overall, the 2D model is proven suitable for fast scoping calculations, as (1) it reduces the mesh elements by 16 times and shortens the run time by about 40 times; (2) it has much more flexible remeshing possibilities; and (3) it provides highly consistent results with the 3D model.

The experience gained in modeling the evolution, from past to present, of natural tracer profiles in geologic media can help support safety assessments of disposal concepts for radioactive wastes in deep geologic repositories because the assessments will be based on predictions using similar models of conditions that could evolve, from the present to the future, in the repository's near-field and far-field environments. Solute-transport models were developed in the present study using a forward modeling approach constrained by boundary conditions inferred from the paleo-hydrogeological evolution of the Horonobe area in Hokkaido, Japan and ranges in model parameter values determined in laboratory and field-based investigations at two deep borehole locations in the area. The models were calibrated by adjusting individual parameter values within limits imposed by the site-characterization data to optimize agreement between model predictions and observations. Results from the calibrated models were consistent with results from other studies suggesting that observed variations with sampling depth in Cl concentrations and 818O and 8D values resulted primarily from an advective transport regime that existed for a relatively brief period during the most recent glacial-interglacial cycle. Advective transport was possible during this time because thawing of a discontinuous permafrost zone starting about 18 ka allowed meteoric water to recharge the flow system and hydraulic gradients were larger than at present because sea levels were significantly lower than the current levels reached at the beginning of the Holocene (i.e., since about 12 ka). Similar climatic controls on the generation and duration of advective/diffusive transport regimes likely existed during nine previous glacial-interglacial cycles that occurred in the Horonobe area since regional uplift began about 1 Ma. If so, alternative versions of the transport models suggest the tracer profiles observed in this area may have evolved over a much longer cumulative period of advective transport and at lower Darcy velocities than assumed in the calibrated models but in a manner that is still compatible with ranges in model parameter values determined in laboratory and field investigations. Apparent differences in transport behavior at the two borehole locations considered in this study, which were situated only about 1 km apart, appear to have resulted from relatively small differences in accessible porosity and hydraulic conductivity, which in turn may have been controlled by local differences in fracture density and fracture connectivity.

Copper canisters are a central component in the safety of the Finnish spent fuel repository concept (KBS-3), where the main corrodent potentially affecting the canister integrity is sulfide. In this study, a 3D numerical model is developed to assess the evolution of sulfide fluxes and the spatially resolved canister corrosion depths for the Finnish spent nuclear fuel repository concept. The backfilled tunnel and the disposal hole are implemented using repository geometries, with sulfide being produced at their interface with the rock (excavation damaged zone) by sulfate reducing bacteria (SRB). Recent experimental findings regarding the microbial sulfate reduction process as well as the scavenging of sulfide via iron (oxy)hydroxides are incorporated in the reactive transport model. Long-term simulations are performed, predicting a heterogeneous corrosion of the canister with a max. corrosion depth of 1.3 mm at the bottom corner after one million years. The evolution of sulfide fluxes shows two main phases, depending on the source of sulfate: first sulfate is supplied by the dissolution of gypsum from the bentonite barriers, followed by a steady, low-level supply from the groundwater. Sensitivity cases demonstrate that both the organic carbon and Fe(III) oxide contents in the bentonite are critical to the corrosion evolution, by being the main electron donor for SRB activities and the major sulfide scavenger in the bentonite, respectively. The backfilled tunnel contributes little to the flux of corrosive sulfide to the canister due to the attenuation by Fe(III)-oxides/hydroxides but induces a notable flux of sulfate into the disposal hole.

Abstract Numerical modeling is used to understand the regional scale flow dynamics of the fault‐hosted orogenic geothermal system at the Grimsel Mountain Pass in the Swiss Alps. The model is calibrated against observations from thermal springs discharging in a tunnel some 250 m underneath Grimsel Pass to derive estimates for the bulk permeability of the fault. Simulations confirm that without the fault as a hydraulic conductor the thermal springs would not exist. Regional topography alone drives meteoric water in a single pass through the fault plane where it penetrates to depths exceeding 10 km and acquires temperatures in excess of 250°C. Thermal constraints from the thermal springs at Grimsel Pass suggest bulk fault permeabilities in the range of 2e−15 m2–4.8e−15 m2. Reported residence times of >30,000 and 7 years for the deep geothermal and shallow groundwater components in the thermal spring water, respectively, suggest fault permeabilities of around 2.5e−15 m2. We show that the long residence time of the deep geothermal water is likely a consequence of low recharge rates during the last glaciation event in the Swiss Alps, which started some 30,000 years ago. Deep groundwater discharging at Grimsel Pass today thus infiltrated the Grimsel fault prior to the last glaciation event. The range of permeabilities estimated from observational constraints is fully consistent with a subcritical single‐pass flow system in the fault plane.

Supporting information for "Reaction mechanism and water/rock ratios involved in epidosite alteration of the oceanic crust" by Samuel Weber, Larryn W. Diamond, Peter Alt-Epping and Alannah C. Brett-Adams. Figure S1: Molar Mg/(Mg + Fe) of chlorite across spilite–epidosite reaction fronts. Plots transition from least epidosite-altered (left) to most altered (right) part of thin-sections. Figure S2: Mineral composition of sample “Median” at W/R ratio 40, 75, 3000 and 30000. Calculations at 350 °C, 50 MPa. Table S1: Bulk-rock chemical compositions of spilites (spl) and epidosites (epi), their sampling locations, rock types and densities. Table S2: Chlorite EMPA spot analyses, reported in wt% oxides. Reported analyses of each sample start in least epidotized area of thin-section and increase in alteration grade Table S3: Electron microprobe analyses of epidote in spilitized Geotimes lavas (Sample AB17-32-spl) Table S4: Dimensionless element enrichment factors (Equation 1) from reactive-transport simulations at 350 °C Table S5: Minerals and aqueous species used in the reactive-transport simulations Table S6: Electron microprobe analyses of relict igneous plagioclase in spilitized Geotimes lavas (Sample AB17-38-spl)

Many meteoric‐recharged, fault‐hosted geothermal systems in amagmatic orogenic belts have been active through the Pleistocene glacial/interglacial climate fluctuations. The effects of climate‐induced recharge variations on fluid flow patterns and residence times of the thermal waters are complex and may influence how the geothermal and mineralization potential of the systems are evaluated. We report systematic thermal‐hydraulic simulations designed to reveal the effects of recharge variations, using a model patterned on the orogenic geothermal system at Grimsel Pass in the Swiss Alps. Previous studies have shown that fault‐bounded circulation of meteoric water is driven to depths of ∼10 km by the high alpine topography. Simulations suggest that the current single‐pass flow is typical of interglacial periods, during which (a) meteoric recharge into the fault is high (above tens of centimeters per year), (b) conditions are at or somewhat below the critical Rayleigh number, and (c) hydraulic connectivity along the fault plane is extensive (an extent of at least 10 km into increasingly higher terrain is required to explain the 10 km penetration depth). The subcritical condition constrains the bulk fault permeability to <1e‐14 m 2 . In contrast, the limited recharge during the numerous Pleistocene glaciation events likely induced a layered flow system, with single‐pass flow confined to shallow depths while non‐Rayleigh convection occurred deeper in the fault. The same layering can be observed at low aspect ratios (length/depth) of the fault plane, when the available recharge area limits flux through the fault.

Epidosites are a prominent type of subseafloor hydrothermal alteration of basalts in ophiolites and greenstone belts, showing an end‐member mineral assemblage of epidote + quartz + titanite + Fe‐oxide. Epidosites are known to form within crustal‐scale upflow zones and their fluids have been proposed as deep equivalents of black‐smoker seafloor vent fluids. Proposals of the mass of fluid per mass of rock ( W / R ratio) needed to form epidosites are contradictory, varying from 20 (Sr isotopes) to > 1,000 (Mg mobility). To test these proposals we have conducted a petrographic, geochemical and reactive‐transport numerical simulation study of the chemical reaction that generates km 3 ‐size epidosite zones within the lavas and sheeted dike complex of the Samail ophiolite, Oman. At 250–400°C the modeled epidosite‐forming fluid has near‐neutral pH (∼ 5.2), high f O 2 , low sulfur and very low Fe (10 −6 mol/kg) contents. These features argue against a genetic link with black‐smoker fluids. Chemical buffering by the epidosite fluid enriches the precursor spilites in Ca and depletes them in Na and Mg. Completion of the spilite‐to‐epidosite reaction requires enormous W / R ratios of 700–∼40,000, depending on initial Mg content and temperature. Collectively, the variably altered rocks in the Samail epidosite zones record flow of ∼10 15 kg of fluid through each km 3 of precursor spilite rock. This fluid imposed on the epidosite an Sr‐isotope signature inherited from the previous rock‐buffered chemical evolution of the fluid through the oceanic crust, thereby explaining the apparently contradictory low W/R ratios based on Sr isotopes.

In order to assess the thermo-hydraulic modelling capabilities of various geothermal simulators, a comparative test suite was created, consisting of a set of cases designed with conditions relevant to the low-enthalpy range of geothermal operations within the European HEATSTORE research project. In an effort to increase confidence in the usage of each simulator, the suite was used as a benchmark by a set of 10 simulators of diverse origin, formulation, and licensing characteristics: COMSOL, MARTHE, ComPASS, Nexus-CSMP++, MOOSE, SEAWATv4, CODE_BRIGHT, Tough3, PFLOTRAN, and Eclipse 100. The synthetic test cases (TCs) consist of a transient pressure test verification (TC1), a well-test comparison (TC2), a thermal transport experiment validation (TC3), and a convection onset comparison (TC4), chosen to represent well-defined subsets of the coupled physical processes acting in subsurface geothermal operations. The results from the four test cases were compared among the participants, to known analytical solutions, and to experimental measurements where applicable, to establish them as reference expectations for future studies. A basic description, problem specification, and corresponding results are presented and discussed. Most participating simulators were able to perform most tests reliably at a level of accuracy that is considered sufficient for application to modelling tasks in real geothermal projects. Significant relative deviations from the reference solutions occurred where strong, sudden (e.g. initial) gradients affected the accuracy of the numerical discretization, but also due to sub-optimal model setup caused by simulator limitations (e.g. providing an equation of state for water properties).

Epidosites (basalts altered to epidote + quartz with minor titanite and Fe-oxide) have been proposed to form along the deep segments of hydrothermal upwelling zones in oceanic crust, possibly leading up to VMS deposits at the seafloor. During epidosite alteration of the precursor rock, Na and Mg are leached completely and Ca is added. This metasomatic exchange requires substantial fluid flux but the exact amount is unclear. We use reactive-transport modeling to calculate the mass of fluid per mass of rock (water- rock ratio) needed to form complete epidosites. Our simulated mineralogical changes during epidosite alteration closely match changes measured in samples from the Oman ophiolite, supporting the validity of our approach. Calculated water-rock ratios for epidosite alteration range from 1200 up to as high as 80,000, strongly depending on the reactant rock composition. Such high water-rock ratios support the proposal that epidosites form along focused flow paths of discharging hydrothermal fluid.

Tables and figures containing the data used in the manuscript called "Quantification of 3D thermal anomalies from surface observations of an orogenic geothermal system (Grimsel Pass, Swiss Alps)", which was submitted to JGR:Solid Earth

Geothermal systems in amagmatic orogens involve topography‐driven infiltration of meteoric water up to 10 km deep into regional‐scale faults and exfiltration of the heated water in surface springs. The thermal anomalies along the upflow zones have not been quantified, yet they are key to estimating the geothermal exploitation potential of such systems. Here we quantify the three‐dimensional heat anomaly below the orogenic geothermal system at Grimsel Pass, Swiss Alps, where warm springs emanate from an exhumed, fossil hydrothermal zone. We use discharge rates and temperatures of the springs, temperature measurements along a shallow tunnel, and the formation temperature and depth of the fossil system to constrain coupled thermal–hydraulic numerical simulations of the upflow zone. The simulations reveal that upflow rates act as a first‐order control on the temperature distribution and that the site is underlain by an ellipsoidal thermal plume enclosing 10 2 –10 3 PJ of anomalous heat per km depth. When the fossil system was active (3.3 Ma), the thermal plume was double its present size, corresponding to a theoretical petrothermal power output of 30–220 MW, with the 120 °C threshold for geothermal electricity production situated at less than 2‐km depth. We conclude that mountainous orogenic belts without igneous activity and even with only low background geothermal gradients typical of waning orogens are surprisingly promising plays for petrothermal power production. Our study implies exploration should focus on major valley floors because there the hydraulic head gradients and thus upflow rates and heat anomalies reach maximum values.

Porosity changes due to mineral dissolution–precipitation reactions in porous media and the resulting impact on transport parameters influence the evolution of natural geological environments or engineered underground barrier systems. In the absence of long-term experimental studies, reactive transport codes are used to evaluate the long-term evolution of engineered barrier systems and waste disposal in the deep underground. Examples for such problems are the long-term fate of CO2 in saline aquifers and mineral transformations that cause porosity changes at clay–concrete interfaces. For porosity clogging under a diffusive transport regime and for simple reaction networks, the accuracy of numerical codes can be verified against analytical solutions. For clogging problems with more complex chemical interactions and transport processes, numerical benchmarks are more suitable to assess model performance, the influence of thermodynamic data, and sensitivity to the reacting mineral phases. Such studies increase confidence in numerical model descriptions of more complex, engineered barrier systems. We propose a reactive transport benchmark, considering the advective–diffusive transport of solutes; the effect of liquid-phase density on liquid flow and advective transport; kinetically controlled dissolution–precipitation reactions causing porosity, permeability, and diffusivity changes; and the formation of a solid solution. We present and analyze the results of five participating reactive transport codes (i.e., CORE2D, MIN3P-THCm, OpenGeoSys-GEM, PFLOTRAN, and TOUGHREACT). In all cases, good agreement of the results was obtained.

Diffusion is the main transport mechanism in many argillaceous formations. In this study, tracer and ion profiles in a 280 m thick clay-rich sequence were simulated by single component and multicomponent diffusion modelling. Drillcores from this sequence originating from a deep borehole in Schlattingen (NE Switzerland) had been previously extensively analysed in terms of porewater chemistry, mineralogy and diffusion parameters. In particular, data from high-pressure core squeezing had enabled to obtain depth profiles of major solutes and water tracers over the entire sequence. The hydrogeological conditions at the site were constrained in the model by the analogy of the nearby site at Benken. In a first step, a simple single component diffusion (SCD) model was set up to simulate the profiles of conservative tracers (delta H-2, delta O-18, Cl-), to check reasonable boundary conditions for the adjacent aquifers and to estimate characteristic diffusion times. Based on these findings, a multicomponent diffusion (MCD) model considering explicitly diffusion in the electrical double layer (EDL) and the "free" water and a chemical equilibrium model was used to simulate the diffusion of major cations (Na+, Ca2+, Mg2+, K+, Sr2+) and anions (Cl-, SO42-, HCO3-). The SCD modelling resulted in a good match of the measured water tracer and chloride profiles in spite of the uncertainty in the conditions regarding the surrounding aquifers. Diffusion times of 0.5-1 Ma were deduced which are in the same range as those postulated previously for the Benken site. Using the same type of boundary conditions, a reasonably good fit of the measured major cation and anion data could be obtained with the MCD model. The results were not sensitive to uncertainties inherent in the MCD model, such as the extent of surface charge screening by fixed cations or the thickness of the EDL. This supports the robustness of the model approach as long as key features such as anion exclusion are captured. Overall, the suitability of the MCD model for simulating cation and anion fluxes in argillaceous rocks over large distances and long timescales could be established. The results also support the validity of squeezing data from drillcores as proxy for in-situ porewater data.

Understanding ion transport through clays and clay membranes is important for many geochemical and environmental applications. Ion transport is affected by electrostatic forces exerted by charged clay surfaces. Anions are partly excluded from pore water near these surfaces, whereas cations are enriched. Such effects can be modeled by the Donnan approach. Here we introduce a new, comparatively simple way to represent Donnan equilibria in transport simulations. We include charged surfaces as immobile ions in the balance equation and calculate coupled transport of all components, including the immobile charges, with the Nernst-Planck equation. This results in an additional diffusion potential that influences ion transport, leading to Donnan ion distributions while maintaining local charge balance. The validity of our new approach was demonstrated by comparing Nernst-Planck simulations using the reactive transport code Flotran with analytical solutions available for simple Donnan systems. Attention has to be paid to the numerical evaluation of the electrochemical migration term in the Nernst-Planck equation to obtain correct results for asymmetric electrolytes. Sensitivity simulations demonstrate the influence of various Donnan model parameters on simulated anion accessible porosities. It is furthermore shown that the salt diffusion coefficient in a Donnan pore depends on local concentrations, in contrast to the aqueous salt diffusion coefficient. Our approach can be easily implemented into other transport codes. It is versatile and facilitates, for instance, assessing the implications of different activity models for the Donnan porosity.