feat: initial commit on computing density

Author: yucongalicechen

news/compute-density.rst

@@ -0,0 +1,23 @@
+**Added:**
+
+* Function that computes theoretical density from a given CIF metadata file.
+
+**Changed:**
+
+* <news item>
+
+**Deprecated:**
+
+* <news item>
+
+**Removed:**
+
+* <news item>
+
+**Fixed:**
+
+* <news item>
+
+**Security:**
+
+* <news item>

src/diffpy/utils/tools.py

@@ -4,6 +4,8 @@
 from pathlib import Path
 
 import numpy as np
+from periodictable import formula
+from scipy.constants import Avogadro as AVOGADRO_NUMBER
 from scipy.optimize import dual_annealing
 from scipy.signal import convolve
 from xraydb import material_mu
@@ -190,8 +192,9 @@ def get_package_info(package_names, metadata=None):
     Package info stored in metadata as
     {'package_info': {'package_name': 'version_number'}}.
 
+    Parameters
     ----------
-    package_name : str or list
+    package_names : str or list
         The name of the package(s) to retrieve the version number for.
     metadata : dict
         The dictionary to store the package info. If not provided, a new
@@ -214,6 +217,49 @@ def get_package_info(package_names, metadata=None):
     return metadata
 
 
+def compute_density_from_cif(sample_composition, cif_data):
+    """Compute the theoretical density from given CIF metadata.
+
+    Parameters
+    ----------
+    sample_composition : str
+        The chemical formula of the material, e.g. "NaCl".
+    cif_data : dict
+        The dictionary containing CIF metadata,
+        typically parsed from a JSON file retrieved
+        from the Crystallography Open Database (COD).
+
+    Returns
+    -------
+    density : float
+        The material density in g/cm^3.
+    """
+    molar_mass = formula(sample_composition).mass
+    a, b, c = (float(cif_data[k]) for k in ("a", "b", "c"))
+    alpha_deg, beta_deg, gamma_deg = (
+        float(cif_data[k]) for k in ("alpha", "beta", "gamma")
+    )
+    Z = int(float(cif_data["Z"]))
+    alpha_rad, beta_rad, gamma_rad = map(
+        np.radians, [alpha_deg, beta_deg, gamma_deg]
+    )
+    volume = (
+        a
+        * b
+        * c
+        * np.sqrt(
+            1
+            - np.cos(alpha_rad) ** 2
+            - np.cos(beta_rad) ** 2
+            - np.cos(gamma_rad) ** 2
+            + 2 * np.cos(alpha_rad) * np.cos(beta_rad) * np.cos(gamma_rad)
+        )
+    )
+    volume_cm3 = volume * 1e-24
+    density = (Z * molar_mass) / (volume_cm3 * AVOGADRO_NUMBER)
+    return density
+
+
 def get_density_from_cloud(sample_composition, mp_token=""):
     """Function to get material density from the MP or COD database.
 

tests/test_tools.py

@@ -9,6 +9,7 @@
 from diffpy.utils.tools import (
     _extend_z_and_convolve,
     check_and_build_global_config,
+    compute_density_from_cif,
     compute_mu_using_xraydb,
     compute_mud,
     get_package_info,
@@ -270,6 +271,28 @@ def test_get_package_info(monkeypatch, inputs, expected):
     assert actual_metadata == expected
 
 
+@pytest.mark.parametrize(
+    "inputs, expected_density",
+    [
+        (
+            {
+                "sample_composition": "NaCl",
+                "cif_data_filename": "cif_data.json",
+            },
+            2.187,
+        ),
+    ],
+)
+def test_compute_density_from_cif(inputs, expected_density):
+    path = Path("testdata") / inputs["cif_data_filename"]
+    with open(path) as f:
+        cif_data = json.load(f)
+    actual_density = compute_density_from_cif(
+        inputs["sample_composition"], cif_data
+    )
+    assert actual_density == pytest.approx(expected_density, rel=0.01, abs=0.1)
+
+
 @pytest.mark.parametrize(
     "inputs",
     [

tests/testdata/cif_data.json

@@ -0,0 +1,75 @@
+{
+  "file": "1000041",
+  "a": "5.62",
+  "siga": null,
+  "b": "5.62",
+  "sigb": null,
+  "c": "5.62",
+  "sigc": null,
+  "alpha": "90",
+  "sigalpha": null,
+  "beta": "90",
+  "sigbeta": null,
+  "gamma": "90",
+  "siggamma": null,
+  "vol": "177.5",
+  "sigvol": null,
+  "celltemp": null,
+  "sigcelltemp": null,
+  "diffrtemp": null,
+  "sigdiffrtemp": null,
+  "cellpressure": null,
+  "sigcellpressure": null,
+  "diffrpressure": null,
+  "sigdiffrpressure": null,
+  "thermalhist": null,
+  "pressurehist": null,
+  "compoundsource": null,
+  "nel": "2",
+  "sg": "F m -3 m",
+  "sgHall": "-F 4 2 3",
+  "sgNumber": "225",
+  "commonname": null,
+  "chemname": "Sodium chloride",
+  "mineral": null,
+  "formula": "- Cl Na -",
+  "calcformula": "- Cl Na -",
+  "cellformula": "- Cl4 Na4 -",
+  "Z": "4",
+  "Zprime": "0.0208333",
+  "acce_code": null,
+  "authors": "Abrahams, S C; Bernstein, J L",
+  "title": "Accuracy of an automatic diffractometer. measurement of the sodium chloride structure factors",
+  "journal": "Acta Crystallographica (1,1948-23,1967)",
+  "year": "1965",
+  "volume": "18",
+  "issue": null,
+  "firstpage": "926",
+  "lastpage": "932",
+  "doi": "10.1107/S0365110X65002244",
+  "method": null,
+  "radiation": null,
+  "wavelength": null,
+  "radType": null,
+  "radSymbol": null,
+  "Rall": "0.022",
+  "Robs": null,
+  "Rref": null,
+  "wRall": null,
+  "wRobs": null,
+  "wRref": null,
+  "RFsqd": null,
+  "RI": null,
+  "gofall": null,
+  "gofobs": null,
+  "gofgt": null,
+  "gofref": null,
+  "duplicateof": null,
+  "optimal": null,
+  "status": null,
+  "flags": "has coordinates",
+  "svnrevision": "130149",
+  "date": "2020-10-21",
+  "time": "18:00:00",
+  "onhold": null
+}