"""Power-law melt viscosity for a two-component blend.

eta = K * shear_rate ** (n - 1)
Fold a function, open the symbol list, or press Ctrl+F.
"""

from dataclasses import dataclass


@dataclass
class PowerLaw:
    """Consistency K is Pa*s^n. Index n is dimensionless."""

    name: str
    consistency: float
    index: float

    def viscosity(self, shear_rate: float) -> float:
        if shear_rate <= 0:
            raise ValueError("shear rate must be positive")
        return self.consistency * shear_rate ** (self.index - 1)


def blend_viscosity(parts: list[tuple[PowerLaw, float]], shear_rate: float) -> float:
    total = sum(weight for _, weight in parts)
    if total <= 0:
        raise ValueError("weights must sum above zero")
    mixed = 0.0
    for model, weight in parts:
        mixed += model.viscosity(shear_rate) * (weight / total)
    return mixed


def describe(model: PowerLaw, shear_rate: float) -> str:
    eta = model.viscosity(shear_rate)
    return f"{model.name}: eta({shear_rate:g} 1/s) = {eta:,.1f} Pa*s"


if __name__ == "__main__":
    hdpe = PowerLaw("HDPE", consistency=8500, index=0.42)
    ldpe = PowerLaw("LDPE", consistency=3200, index=0.55)
    shear = 100.0
    print(describe(hdpe, shear))
    print(describe(ldpe, shear))
    mix = blend_viscosity([(hdpe, 0.7), (ldpe, 0.3)], shear)
    print(f"70/30 blend: {mix:,.1f} Pa*s")