diff --git a/TrimSPlib.js b/TrimSPlib.js
index 9dfb54e..81529d1 100644
--- a/TrimSPlib.js
+++ b/TrimSPlib.js
@@ -3080,28 +3080,40 @@ function rho_fun()
 	if (isFinite(rhos[chem])) {
 	    // Found in densities list of compositions
 	    rho = rhos[chem];
-	} else {
-	    // suggest a material density based on composition (when it's undefined)
-            let layer_formula = parse_formula(chem);
+		} else {
+		    // suggest a material density based on composition (when it's undefined)
+	            let layer_formula = parse_formula(chem);
 
-            // determine the stoichiometry sum (for normalization)
-            let stoichiometry_sum = 0;
-	    for (key of Object.keys(layer_formula)) {
-	        stoichiometry_sum = stoichiometry_sum + layer_formula[key];
+	            // determine the stoichiometry sum (for normalization)
+	            let stoichiometry_sum = 0;
+		    for (const key of Object.keys(layer_formula)) {
+		        const amount = Number(layer_formula[key]);
+		        if (isFinite(amount) && amount > 0) {
+			    stoichiometry_sum += amount;
+			}
+		    }
+
+		    // determine the density using a weighted average of the elemental densities
+		    let density_estimate = 0;
+		    for (const key of Object.keys(layer_formula)) {
+		        const amount = Number(layer_formula[key]);
+		        if (!isFinite(amount) || amount <= 0) {
+			    continue;
+			}
+			if (!elemPars[key] || !isFinite(elemPars[key].rho)) {
+			    continue;
+			}
+		        if (stoichiometry_sum > 0) {
+			    density_estimate += (amount / stoichiometry_sum) * elemPars[key].rho;
+			}
+		    }
+
+	            rho = density_estimate;
+		    if (rho > 0) {
+			alert("Warning: Estimated density from elemental constituents; please verify manually.");
+		    }
+		}
 	    }
-
-	    // determine the density using on a weighted average of the elemental densities
-	    let density_estimate = 0;
-	    for (key of Object.keys(layer_formula)) {
-	        if (layer_formula[key] != null) {
-	            density_estimate = density_estimate + (layer_formula[key] / stoichiometry_sum) * elemPars[key].rho;
-	        }
-	    }
-
-            rho = density_estimate;
-	    alert("Warning: The density for this layer is only an estimate!")
-	}
-    }
     // Set value in appropriate cell
     document.getElementById(rhoLi).value = rho;    
 }
@@ -3769,3 +3781,123 @@ function amIWeb() {
     }
 }
 
+async function getDensity(formula) {
+    function normalizeCodFormula(codFormula) {
+	return codFormula.replace(/-/g, ' ').replace(/\s+/g, '').trim();
+    }
+
+    function sameFormula(formulaA, formulaB) {
+	const a = parse_formula(formulaA);
+	const b = parse_formula(formulaB);
+	const aKeys = Object.keys(a).sort();
+	const bKeys = Object.keys(b).sort();
+	if (aKeys.length !== bKeys.length) {
+	    return false;
+	}
+	for (let i = 0; i < aKeys.length; i++) {
+	    const key = aKeys[i];
+	    if (key !== bKeys[i] || Math.abs(a[key] - b[key]) > 1e-8) {
+		return false;
+	    }
+	}
+	return true;
+    }
+
+    function scaleFormula(formulaObj, factor) {
+	const scaled = {};
+	for (const key of Object.keys(formulaObj)) {
+	    scaled[key] = formulaObj[key] * factor;
+	}
+	return scaled;
+    }
+
+    function sameFormulaObject(formulaA, formulaB) {
+	const aKeys = Object.keys(formulaA).sort();
+	const bKeys = Object.keys(formulaB).sort();
+	if (aKeys.length !== bKeys.length) {
+	    return false;
+	}
+	for (let i = 0; i < aKeys.length; i++) {
+	    const key = aKeys[i];
+	    if (key !== bKeys[i] || Math.abs(formulaA[key] - formulaB[key]) > 1e-8) {
+		return false;
+	    }
+	}
+	return true;
+    }
+
+    function molarMass(formulaObj) {
+	let mass = 0;
+	for (const key of Object.keys(formulaObj)) {
+	    if (!elemPars[key] || !elemPars[key].A) {
+		return null;
+	    }
+	    mass += formulaObj[key] * elemPars[key].A;
+	}
+	return mass;
+    }
+
+    try {
+	const parsedFormula = parse_formula(formula);
+	const formulaMass = molarMass(parsedFormula);
+	if (formulaMass == null) {
+	    return null;
+	}
+
+	const encodedFormula = encodeURIComponent(formula);
+	const codRes = await fetch(`https://www.crystallography.net/cod/result?formula=${encodedFormula}&format=json`);
+	if (!codRes.ok) {
+	    throw new Error(`COD lookup failed with status ${codRes.status}`);
+	}
+
+	let records = await codRes.json();
+	if ((!Array.isArray(records) || records.length === 0) && Object.keys(parsedFormula).length > 0) {
+	    const params = new URLSearchParams();
+	    const elements = Object.keys(parsedFormula).sort();
+	    elements.forEach((element, index) => {
+		params.append(`el${index + 1}`, element);
+	    });
+	    params.append('strictmin', elements.length);
+	    params.append('strictmax', elements.length);
+	    params.append('format', 'json');
+
+	    const fallbackRes = await fetch(`https://www.crystallography.net/cod/result?${params.toString()}`);
+	    if (!fallbackRes.ok) {
+		throw new Error(`COD fallback lookup failed with status ${fallbackRes.status}`);
+	    }
+	    records = await fallbackRes.json();
+	}
+	if (!Array.isArray(records) || records.length === 0) {
+	    return null;
+	}
+
+	for (const record of records) {
+	    const z = parseFloat(record.Z);
+	    const vol = parseFloat(record.vol);
+	    if (!isFinite(z) || !isFinite(vol) || z <= 0 || vol <= 0) {
+		continue;
+	    }
+
+	    const recordFormula = normalizeCodFormula(record.formula || '');
+	    if (!sameFormula(formula, recordFormula)) {
+		continue;
+	    }
+
+	    const cellFormula = normalizeCodFormula(record.cellformula || '');
+	    if (cellFormula !== '') {
+		const parsedCellFormula = parse_formula(cellFormula);
+		const expectedCellFormula = scaleFormula(parsedFormula, z);
+		if (!sameFormulaObject(parsedCellFormula, expectedCellFormula)) {
+		    continue;
+		}
+	    }
+
+	    // Convert unit-cell data (A^3, Z) into density in g/cm^3.
+	    return 1.66053906660 * z * formulaMass / vol;
+	}
+    } catch (error) {
+	console.error("Error fetching data from COD:", error);
+    }
+
+    return null;
+}