Metadata-Version: 2.4
Name: GalvCalc
Version: 1.0.1
Summary: A computational framework for modeling micro-galvanic corrosion of alloys with coupled anodic dissolution and hydrogen evolution kinetics
Author-email: Gaoning Shi <gaoning_shi@sjtu.edu.cn>
License-Expression: MIT
Project-URL: Homepage, https://github.com/GNscSHI/GalvCalc
Project-URL: Documentation, https://github.com/GNscSHI/GalvCalc#readme
Project-URL: Source, https://github.com/GNscSHI/GalvCalc
Keywords: corrosion,galvanic corrosion,electrochemistry,density functional theory,materials science,hydrogen evolution,pymatgen,polarization
Classifier: Development Status :: 4 - Beta
Classifier: Intended Audience :: Science/Research
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.9
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Topic :: Scientific/Engineering :: Chemistry
Classifier: Topic :: Scientific/Engineering :: Physics
Classifier: Topic :: Scientific/Engineering :: Visualization
Requires-Python: >=3.9
Description-Content-Type: text/markdown
License-File: LICENSE
Requires-Dist: numpy>=1.21
Requires-Dist: scipy>=1.7
Requires-Dist: pandas>=1.3
Requires-Dist: matplotlib>=3.4
Requires-Dist: PyYAML>=5.4
Requires-Dist: sympy>=1.9
Requires-Dist: pymatgen>=2023.10.4
Requires-Dist: tqdm>=4.60
Requires-Dist: joblib>=1.1
Requires-Dist: openpyxl>=3.0
Provides-Extra: defects
Requires-Dist: pymatgen<2027,>=2024.10; extra == "defects"
Requires-Dist: pymatgen-analysis-defects<2027,>=2024.10; extra == "defects"
Requires-Dist: scikit-image>=0.19; extra == "defects"
Requires-Dist: mp-pyrho>=0.4.4; extra == "defects"
Provides-Extra: ml
Requires-Dist: torch>=1.10; extra == "ml"
Requires-Dist: scikit-learn>=1.0; extra == "ml"
Requires-Dist: tabpfn>=0.1.9; extra == "ml"
Provides-Extra: dev
Requires-Dist: pytest>=7.0; extra == "dev"
Requires-Dist: build>=0.10; extra == "dev"
Requires-Dist: twine>=4.0; extra == "dev"
Dynamic: license-file

<p align="center">
  <img src="GalvCalc/assets/logo.png" alt="GalvCalc" width="400"/>
</p>

<p align="center">
  <em>From atomic-scale descriptors to macroscopic corrosion polarization curves</em>
</p>

<p align="center">
  <a href="https://pypi.org/project/GalvCalc/"><img alt="PyPI" src="https://img.shields.io/pypi/v/GalvCalc?color=2b6cb0&label=PyPI"></a>
  <a href="https://pypi.org/project/GalvCalc/"><img alt="Python" src="https://img.shields.io/pypi/pyversions/GalvCalc?color=2b6cb0"></a>
  <a href="LICENSE"><img alt="License" src="https://img.shields.io/github/license/GNscSHI/GalvCalc?color=2b6cb0"></a>
</p>

**GalvCalc** models the micro-galvanic corrosion of alloys with coupled
anodic dissolution and hydrogen-evolution kinetics. Starting from a bulk
crystal structure, it builds the surface models, estimates the
electrochemical descriptors (equilibrium potentials, work functions, surface
energies and hydrogen adsorption free energies), and assembles multi-phase
polarization curves — with built-in anode/cathode area-ratio optimization and
alloy-content scans of the corrosion current and potential.

The package is developed alongside the manuscript *"GalvCalc: A Framework for
Modeling Micro-galvanic Corrosion of Alloys with Coupled Anodic Dissolution
and Hydrogen Evolution Kinetics"*.

---

## Features

| Module | Purpose |
| --- | --- |
| `GalvCalc.core` | `Bulk` / `Surface` structure classes, Nernst-equation equilibrium potentials (bundled aqueous-ion thermodynamic data), Wulff-aware plotting helpers |
| `GalvCalc.cathode` | Surface-property manager (surface energies, work functions, Wulff shapes), hydrogen-adsorption analysis, facet-weighted exchange-current estimation |
| `GalvCalc.anode` | Substitutional surface doping of anode materials and electrochemical descriptor tables |
| `GalvCalc.polarization` | Butler–Volmer polarization curves, multi-anode/multi-cathode plots, area-ratio optimization, alloy-content scans |
| `GalvCalc.predictor` | ML predictors (CGCNN / TabPFN) for surface properties and hydrogen adsorption free energies |

**Highlights**

- Full pipeline: bulk crystal → slab generation → Wulff shape →
  facet-weighted descriptors → corrosion polarization curves.
- Nernst-equation equilibrium potentials for arbitrary ion combinations (Mg,
  Fe, Zn, Al, ...), backed by a bundled database of aqueous-ion formation
  free energies.
- Calibrated Butler–Volmer kinetics for both Mg- and Fe-based systems.
- DFT-grounded second-phase exchange currents for nine common Mg-alloy
  intermetallics (Mg<sub>17</sub>Al<sub>12</sub>, Mg<sub>2</sub>Al<sub>3</sub>, MgZn<sub>2</sub>, LaMg<sub>12</sub>, CaMg<sub>2</sub>, Y<sub>5</sub>Mg<sub>24</sub>, NdMg<sub>3</sub>, Mg<sub>2</sub>Si, CeMg<sub>12</sub>).
- CGCNN prediction of surface energies and work functions, and
  TabPFN-based prediction of hydrogen adsorption free energies.

---

## Installation

```bash
pip install GalvCalc
```

Core dependencies are installed by default. The `ml` extra adds the
pre-trained predictors (CGCNN / TabPFN), and the `defects` extra (Python ≥
3.10) enables the substitution-based doping workflow in `GalvCalc.anode`.

**Requirements:** Python >= 3.9; `numpy`, `scipy`, `pandas`, `matplotlib`,
`pymatgen`, `PyYAML`, `sympy`, `tqdm`, `joblib`, `openpyxl`. The `ml` extra
adds `torch`, `scikit-learn`, `tabpfn`.

---

## Quick start

All structure files used below (`Mg2Si.vasp`, `Mg.vasp`, `surface/*.vasp`)
and the demo DFT adsorption dataset (`adsorption_analysis.csv`) ship with the
package under `GalvCalc/examples/`, so the snippets run fully offline from
that folder (as in the bundled notebook).

### 1. Equilibrium potential

```python
from GalvCalc.core.structures import Bulk

# Load a bulk structure from a POSCAR / CIF / .vasp file
bulk = Bulk.from_file("Mg2Si.vasp")

# Single-ion reaction: Mg -> Mg[2+] + 2 e-
Ee = Bulk.get_equilibrium_potential(
    ions="Mg[2+]",
    ion_numbers=[1],
    energy_formation=0.0,
)
print(f"Mg E_eq = {Ee:.3f} V vs. SHE")

# Multi-ion compound dissolution: MgZn2 -> Mg[2+] + 2 Zn[2+] + 6 e-
Ee2 = Bulk.get_equilibrium_potential(
    ions=["Mg[2+]", "Zn[2+]"],
    ion_numbers=[1, 2],
    energy_formation=-0.24,
)
print(f"MgZn2 E_eq = {Ee2:.3f} V vs. SHE")
```

### 2. Surfaces, ML-predicted facet properties and Wulff construction

```python
from GalvCalc.core.structures import Bulk
from GalvCalc.cathode.surfaces import SurfaceProperties

bulk = Bulk.from_file("Mg2Si.vasp")

# Batch-generate all surface terminations up to max_index = 1
surface_props = SurfaceProperties.from_bulk_structure(
    bulk,
    max_index=1,
    min_slab_size=15.0,
    min_vacuum_size=15.0,
    center_slab=True,
    symmetrize=True,
    max_normal_search=1,
)

# CGCNN-predicted facet surface energies / work functions (needs the `ml` extra)
df_pred = surface_props.get_predicted_surface_properties()

# Wulff shape from the predicted facet properties
wulff = surface_props.wulff_construct(
    surface_energies=surface_props.surface_energies_dict,
    work_functions=surface_props.work_functions_dict,
    output_dir="demo_output",
    save_plot=True,
    save_csv=True,
)
```

### 3. Hydrogen adsorption and exchange-current estimation

The demo DFT adsorption dataset is bundled as
`adsorption_analysis.csv`; the facet-weighted hydrogen-evolution exchange
current follows from the adsorption free energy through the BEP-type relation
of the manuscript:

```python
import pandas as pd
from GalvCalc.cathode import ic0_mg

df_ads = pd.read_csv("adsorption_analysis.csv")
mg2si = df_ads[df_ads["formula"] == "Mg2Si"]
print(mg2si[["miller_index", "termination", "ads_position", "Eads", "workfunction"]].head())

# Mg2Si (110) H6 site: G_ads = 0.028 eV, work function = 3.7378 eV
i0 = ic0_mg(delta_G=0.028, wf=3.7378)
print(f"ic0_mg = {i0:.3e} A/cm^2")
```

```python
from GalvCalc.cathode.hydrogen import AdsorptionManager

adsorption_manager = AdsorptionManager(surface_props.surfaces)
results = adsorption_manager.H_adsorption_analysis(
    adsorbate="H",
    site_indices="all",
    output_dir="H_adsorption_analysis",
    include_summary=True,
    include_dataframe=True,
    include_visualization=True,
)

formatted_data = adsorption_manager.get_formatted_data()
adsorption_manager.save_formatted_data("adsorption_data.json")
```

`AdsorptionManager.H_adsorption_analysis` locates the adsorption sites on
every surface, builds the H-adsorbed slabs, predicts the hydrogen adsorption free energies
with the bundled TabPFN model and exports POSCAR files plus CSV/JSON
summaries. This path needs the `ml` extra.

### 4. Substitutional doping of anode surfaces

```python
from GalvCalc.core.structures import Bulk, Surface
from GalvCalc.anode import SurfaceDopingManager

mg_bulk = Bulk.from_file("Mg.vasp")
mg_surface = Surface.from_bulk(bulk_structure=mg_bulk, miller_index=(0, 0, 1))
manager = SurfaceDopingManager(mg_surface, "Mg_001")

# Substitute surface Mg by Zn, Al and Y (second layer by default)
manager.batch_dope(["Zn", "Al", "Y"], save_to_file=False)
for name, info in manager.doping_info.items():
    print(f"{name}: host {info.host_element} -> {info.dopant_element}, "
          f"site {info.site_index}, layer {info.layer}, "
          f"multiplicity {info.multiplicity}, depth {info.depth:.2f} A")

manager.set_property("Mg_001", "work_function", 3.72)
manager.set_property("Mg_001", "vacancy_energy", 0.84)
manager.set_property("Mg_001_Zn", "work_function", 3.75)
manager.set_property("Mg_001_Zn", "vacancy_energy", 0.79)
manager.set_property("Mg_001_Al", "work_function", 3.74)
manager.set_property("Mg_001_Al", "vacancy_energy", 0.84)
manager.set_property("Mg_001_Y", "work_function", 3.29)
manager.set_property("Mg_001_Y", "vacancy_energy", 0.78)

df = manager.calculate_electrochemical_properties(E00=-2.37, ia00=1e-5)
print(df[["surface_name", "E0_calculated", "ia0_calculated", "dopant_element"]])
```

Set `save_to_file=True` to export the doped structures as POSCAR files into
`doped_surfaces/`; the full doping record (host/dopant, site, layer, symmetry
multiplicity) is available from `manager.doping_info`. This workflow builds on
pymatgen's defect framework and needs the `defects` extra (Python ≥ 3.10).

### 5. Polarization curves and area-ratio optimization

```python
from GalvCalc.polarization import (
    ElectrodeParameters, Composition, plot_single_polarization,
    mg_second_phase_example, plot_mg_second_phases,
)
from GalvCalc.polarization.area_ratio import (
    AreaRatioAnalyzer, create_example_parameters,
)

# Fe-based system: n = 1, exponent (alpha_a + 1) * n, prefactor 2
anode = ElectrodeParameters(4.1e-8, -0.44, 0.5, name="Fe/Fe2+", kinetic_form="fe")
cathode = ElectrodeParameters(7.9e-8, -0.059, 0.5, name="H+/H2", kinetic_form="fe")
comp = Composition("Fe", anode=anode, cathodes=[cathode], area_ratios=[1, 1])
fig = plot_single_polarization(comp, reference_electrode="SHE")

# Mg-based system: nine intermetallic second-phase cathodes on one plot
fig_mg = plot_mg_second_phases(mg_second_phase_example())

# Anode/cathode area-ratio optimization
params = create_example_parameters()
fig, optimal_ratio, i_max = AreaRatioAnalyzer(reference_electrode="SCE").plot_area_ratio_analysis(params)
print(f"optimal anode area ratio = {optimal_ratio:.2%}")
```

The second-phase cathodic exchange currents in `mg_second_phase_example()`
are DFT-calibrated for nine Mg intermetallics; supply your own values to
study other phases.

### 6. Alloy-content scan

```python
import numpy as np
from GalvCalc.polarization import (
    wt_to_vol_fraction, two_phase_kinetics, scan_corrosion_vs_content,
)

# Any intermetallic second phase works: swap the phase formula, densities
# and kinetic parameters. Mg-Nd / Mg3Nd is the example used here.
rho_matrix, rho_phase = 1.8, 3.58  # Mg matrix, Mg3Nd second phase
nd_wt = np.linspace(0, 6, 13)      # Nd content (wt%)

def kinetics_at_content(w):
    vol_fraction = wt_to_vol_fraction(
        w, "Mg3Nd", rho_matrix=rho_matrix, rho_phase=rho_phase
    )
    return two_phase_kinetics(
        vol_fraction=vol_fraction,
        anode_i0=6e-23 if w <= 0 else 3e-23,
        anode_equilibrium_potential=-2.37,
        cathode_i0s=[10 ** -8.1, 10 ** -8.732],
        cathode_equilibrium_potential=-0.61,
    )

df_scan = scan_corrosion_vs_content(
    nd_wt,
    params_fn=kinetics_at_content,
    domain=(0.0, 1.0),
    n_ratios=51,
    n_potentials=1500,
)
df_scan[["content", "max_log10_i_corr", "E_corr_at_max"]].round(4)
```

---

## Example notebooks

Runnable notebooks are shipped under `GalvCalc/examples/`:

- `galvcalc_demo.ipynb` — a step-by-step walkthrough covering equilibrium
  potentials, slab generation, adsorption-site schematics, Wulff construction,
  facet-weighted exchange currents, anode doping, polarization curves,
  area-ratio optimization and alloy-content scans.

Example structures (`*.vasp`, `surface/*.vasp`) and the demo DFT adsorption
dataset (`adsorption_analysis.csv`) are bundled in the same folder.

---

## License

Distributed under the [MIT License](LICENSE).

## Contact

Gaoning Shi — gaoning_shi@sjtu.edu.cn
