Add arithemtic std as a cli param

Author: arunoruto

.gitignore

@@ -1,6 +1,7 @@
 .envrc
 .direnv
 RESULTS/
+*.speedscope.json
 
 *.csv
 

docs/source/usage.md

@@ -30,7 +30,8 @@ pyfracval -n 512 --df 1.9 --kf 1.4 --rp-g 50 --rp-gstd 1.25
 - `--df`: Target fractal dimension.
 - `--kf`: Target fractal prefactor.
 - `--rp-g`: Geometric mean radius of primary particles.
-- `--rp-gstd`: Geometric standard deviation of radii (`1.0` for monodisperse).
+- `--rp-gstd`: **Geometric** standard deviation of radii (`>= 1.0`). If provided, **takes precedence** over `--rp-std`. If neither is provided, defaults to [Default Value, e.g., 1.5].
+- `--rp-std`: Approximate **arithmetic** standard deviation of radii. Used to _estimate_ `--rp-gstd` via `exp(std/mean)` **only if `--rp-gstd` is not given**. A warning will show the estimated geometric value being used.
 - `--tol-ov`: Overlap tolerance (e.g., `1e-5`).
 - `--n-subcl-perc`: Target fraction for PCA subcluster size (e.g., `0.1`).
 - `--num-aggregates`: Generate multiple aggregates sequentially.
@@ -41,6 +42,25 @@ pyfracval -n 512 --df 1.9 --kf 1.4 --rp-g 50 --rp-gstd 1.25
 - `--log-file`: Redirect log output to a file.
 - `-h`, `--help`: Show all available options and their defaults.
 
+**Example: Using Arithmetic Standard Deviation**
+
+Generate an aggregate with N=200, Df=1.9, kf=1.2, geometric mean radius 20, and an approximate _arithmetic_ standard deviation of 5. `pyfracval` will estimate the geometric standard deviation and use that.
+
+```bash
+pyfracval -n 200 --df 1.9 --kf 1.2 --rp-g 20 --rp-std 5 -vv
+```
+
+_(Check the log output for a WARNING indicating the calculated `rp_gstd` value being used)._
+
+**Example: Geometric Standard Deviation Takes Precedence**
+
+If you provide both, `--rp-gstd` will be used:
+
+```bash
+# rp_gstd=1.3 will be used, rp-std=5 will be ignored (with a warning)
+pyfracval -n 100 --df 1.8 --rp-gstd 1.3 --rp-std 5
+```
+
 **Example: Generating Multiple Aggregates with Plotting**
 
 Generate 3 aggregates, each with N=100, Df=1.7, kf=1.1, and show the plots:

pyfracval/app.py

@@ -5,7 +5,9 @@
 from pathlib import Path
 
 import numpy as np
-import polars as pl
+
+# import polars as pl
+import pandas as pd
 import pyvista as pv
 import streamlit as st
 from stpyvista import stpyvista
@@ -49,7 +51,7 @@
 
 match Path(file).suffix:
     case ".csv":
-        data = pl.read_csv(file).to_numpy()
+        data = pd.read_csv(file).to_numpy()
         information = re.search(
             r"N(\d+)-D(\d+_\d+)-K(\d+_\d+)-(\d+)_(\d+)_(\d+)",
             str(file),
@@ -76,8 +78,28 @@
     case _:
         st.error("File type not supported")
 
-pl = plot_particles(data[:, :3], data[:, 3])
-stpyvista(pl)
+plotter = plot_particles(data[:, :3], data[:, 3])
+stpyvista(plotter)
+
+gmean = np.exp(np.mean(np.log(data[:, 3])))
+gstd = np.exp(np.std(np.log(data[:, 3])))
+
+col1, col2 = st.columns([1, 4])
+with col1:
+    st.table(
+        pd.DataFrame(
+            dict(
+                Arithemtic=[np.mean(data[:, 3]), np.std(data[:, 3])],
+                Geometric=[gmean, gstd],
+            ),
+            index=["Mean", "STD"],  # pyright: ignore
+        )
+    )
+with col2:
+    st.write(
+        "Approximate Geometric STD: ",
+        np.exp(np.std(data[:, 3]) / np.mean(data[:, 3])),
+    )
 
 st.write(metadata)
 with st.expander("Full file path"):

pyfracval/cli.py

@@ -4,6 +4,7 @@
 from pathlib import Path
 
 import click
+import numpy as np
 
 from pyfracval import config as default_config
 from pyfracval.logs import TRACE_LEVEL_NUM, create_logger
@@ -65,9 +66,19 @@
 @click.option(
     "--rp-gstd",
     type=float,
-    default=DEFAULT_SIGMA,
-    show_default=True,
-    help="Geometric standard deviation of primary particle radii (>= 1.0).",
+    default=None,
+    show_default=f"Calculated from --rp-std or defaults to {DEFAULT_SIGMA}",  # Show calculated default
+    help="Geometric standard deviation of primary particle radii (>= 1.0). "
+    "If provided, this value takes precedence over --rp-std.",  # Added precedence info
+)
+@click.option(
+    "--rp-std",
+    type=float,
+    default=None,
+    help="Approximate arithmetic standard deviation of primary particle radii. "
+    "If --rp-gstd is NOT provided, this value will be used to estimate "
+    "a geometric standard deviation using the heuristic exp(std/mean). "
+    "A warning will be shown with the calculated geometric value.",
 )
 @click.option(
     "--ext-case",
@@ -145,6 +156,9 @@ def cli(ctx, **kwargs) -> None:
     by Cluster-Cluster Aggregation (CCA) approach described by
     Moran et al. (2019) to generate aggregates with tunable fractal
     dimension (Df) and prefactor (kf) from polydisperse primary particles.
+
+    Allows specifying target Df, kf, and particle size distribution
+    (via geometric mean/std dev or estimated from arithmetic std dev).
     """
     if ctx.invoked_subcommand:
         return
@@ -162,11 +176,6 @@ def cli(ctx, **kwargs) -> None:
     logger = create_logger(log_level, kwargs["log_file"])
 
     # --- Validate Inputs ---
-    if kwargs["rp_gstd"] < 1.0:
-        raise click.BadParameter(
-            "Geometric standard deviation (rp_gstd) must be >= 1.0.",
-            param_hint="--rp-gstd",
-        )
     if kwargs["rp_g"] <= 0:
         raise click.BadParameter(
             "Geometric mean radius (rp_g) must be > 0.", param_hint="--rp-g"
@@ -175,14 +184,45 @@ def cli(ctx, **kwargs) -> None:
         raise click.BadParameter(
             "Number of particles (n) must be at least 2.", param_hint="-n"
         )
+    if kwargs["rp_gstd"] is not None:
+        # Geometric STD provided, use it directly (takes precedence)
+        if kwargs["rp_gstd"] < 1.0:
+            raise click.BadParameter(
+                "Geometric standard deviation (--rp-gstd) must be >= 1.0.",
+                param_hint="--rp-gstd",
+            )
+        final_rp_gstd = kwargs["rp_gstd"]
+        if kwargs["rp_std"] is not None:
+            logger.warning("Both --rp-gstd and --rp-std provided. Using --rp-gstd.")
+    elif kwargs["rp_std"] is not None:
+        # Arithmetic STD provided, Geometric STD not provided
+        if kwargs["rp_std"] < 0:
+            raise click.BadParameter(
+                "Arithmetic standard deviation (--rp-std) cannot be negative.",
+                param_hint="--rp-std",
+            )
+        # Apply heuristic: sigma_g = exp(sigma_a / mu_g)
+
+        final_rp_gstd = np.exp(kwargs["rp_std"] / kwargs["rp_g"])
+        logger.warning(
+            f"Using heuristic to calculate geometric standard deviation from arithmetic std: "
+            f"exp(rp_std / rp_g) = exp({kwargs['rp_std']:.2f} / {kwargs['rp_g']:.2f}) = {final_rp_gstd:.3f}. "
+            f"Targeting rp_gstd = {final_rp_gstd:.3f} for generation."
+        )
+    else:
+        # Neither provided, use the default geometric STD
+        final_rp_gstd = DEFAULT_SIGMA
+        logger.info(
+            f"No geometric or arithmetic standard deviation provided. Using default rp_gstd = {final_rp_gstd:.3f}"
+        )
 
     # --- Prepare Configuration for Runner ---
     sim_config = {
         "N": kwargs["num_particles"],
         "Df": kwargs["df"],
         "kf": kwargs["kf"],
         "rp_g": kwargs["rp_g"],
-        "rp_gstd": kwargs["rp_gstd"],
+        "rp_gstd": final_rp_gstd,
         "tol_ov": kwargs["tol_ov"],
         "n_subcl_percentage": kwargs["n_subcl_perc"],
         "ext_case": kwargs["ext_case"],

pyproject.toml

@@ -35,7 +35,10 @@ docs = [
   "sphinxcontrib-bibtex>=2.6.3",
   "sphinxcontrib-napoleon>=0.7",
 ]
-test = ["pytest>=8.3.3"]
+test = [
+    "py-spy>=0.4.0",
+    "pytest>=8.3.3",
+]
 
 [project.scripts]
 pyfracval = "pyfracval.cli:cli"