import math, random

# 1. Fitting the power law τ = K·γ̇ⁿ to a flow curve. The "measurements" here are simulated
#    (K = 0.80 Pa·sⁿ, n = 0.55, 3 % noise), to show the method; real values must be measured.
random.seed(7)
rates = [1, 2, 5, 10, 20, 50, 100, 200]
stress = [0.80 * g ** 0.55 * random.gauss(1, 0.03) for g in rates]
x = [math.log(g) for g in rates]
y = [math.log(s) for s in stress]
mx, my = sum(x) / len(x), sum(y) / len(y)
n = sum((a - mx) * (b - my) for a, b in zip(x, y)) / sum((a - mx) ** 2 for a in x)
K = math.exp(my - n * mx)
print(f"fitted K = {K:.2f} Pa·sⁿ, n = {n:.2f} (n < 1: shear-thinning)")

# 2. Apparent viscosity η = K·γ̇ⁿ⁻¹, with K scaling as c^3.3 above the overlap concentration
#    (typical of random-coil polysaccharides; Morris et al. 1981) and an Arrhenius temperature
#    dependence (activation energy 20 kJ/mol assumed).
def eta(c, rate, t_c, K1=0.80, n1=0.55, Ea=20e3):
    Kc = K1 * (c / 1.0) ** 3.3 * math.exp(Ea / 8.314 * (1 / (t_c + 273.15) - 1 / 293.15))
    nc = max(0.3, n1 - 0.15 * math.log(c))          # more concentrated, more shear-thinning
    return Kc * rate ** (nc - 1)

print("\napparent viscosity, mPa·s (water at 20 °C ≈ 1 mPa·s)")
print("polysaccharide   20 °C at 10 / 50 s⁻¹     45 °C at 10 / 50 s⁻¹")
for c in (0.25, 0.5, 1.0):
    print(f"   {c:4.2f} %        {1000 * eta(c, 10, 20):6.0f} / {1000 * eta(c, 50, 20):5.0f}"
          f"          {1000 * eta(c, 10, 45):6.0f} / {1000 * eta(c, 50, 45):5.0f}")
# fitted K = 0.79 Pa·sⁿ, n = 0.55 (n < 1: shear-thinning)
#
# apparent viscosity, mPa·s (water at 20 °C ≈ 1 mPa·s)
# polysaccharide   20 °C at 10 / 50 s⁻¹     45 °C at 10 / 50 s⁻¹
#    0.25 %             5 /     3               2 /     2
#    0.50 %            37 /    21              19 /    11
#    1.00 %           284 /   138             149 /    72

import math

# 1. How much polysaccharide is in a cup? Assumed: fresh stems 80 % water; polysaccharide
#    30 % of dry matter (the Chinese Pharmacopoeia requires at least 25 % in the dried herb).
fresh_g, water_frac, poly_dry = 12.0, 0.80, 0.30
poly_g = fresh_g * (1 - water_frac) * poly_dry
for cup_ml in (150, 250):
    print(f"{fresh_g:.0f} g fresh stem in {cup_ml} mL: up to {poly_g:.2f} g polysaccharide, "
          f"{100 * poly_g / cup_ml:.2f} % if all of it dissolved")

# 2. How fast does it come out? First-order release, X = 1 − exp(−k·t), with the rate rising
#    with temperature (Arrhenius). Illustrative: k = 0.10 per minute at 100 °C, Ea = 40 kJ/mol.
def released(t_min, temp_c, k100=0.10, Ea=40e3):
    k = k100 * math.exp(-Ea / 8.314 * (1 / (temp_c + 273.15) - 1 / 373.15))
    return 1 - math.exp(-k * t_min)

print("\nshare released from crushed stem")
print("minutes    60 °C    80 °C   100 °C")
for t in (5, 10, 20, 40):
    print(f"  {t:3d}     {100 * released(t, 60):5.0f} %  {100 * released(t, 80):5.0f} %  "
          f"{100 * released(t, 100):5.0f} %")
# 12 g fresh stem in 150 mL: up to 0.72 g polysaccharide, 0.48 % if all of it dissolved
# 12 g fresh stem in 250 mL: up to 0.72 g polysaccharide, 0.29 % if all of it dissolved
#
# share released from crushed stem
# minutes    60 °C    80 °C   100 °C
#     5        10 %     21 %     39 %
#    10        19 %     38 %     63 %
#    20        35 %     62 %     86 %
#    40        57 %     85 %     98 %