# How Lighthouses Work 3, Exercise 2: what six times the light buys
# Ranges of an assumed 2,000 cd light and one six times as bright, and the factor needed to double the range.
import math

NM = 1852.0

def transmissivity(visibility_nm):
    return 0.05 ** (1 / visibility_nm)

def luminous_range_nm(intensity_cd, visibility_nm, threshold_lx=2e-7):
    t = transmissivity(visibility_nm)
    lo, hi = 0.01, 200.0
    for _ in range(100):
        mid = (lo + hi) / 2
        e = intensity_cd * t ** mid / (mid * NM) ** 2
        if e > threshold_lx:
            lo = mid
        else:
            hi = mid
    return (lo + hi) / 2

def geographic_range_nm(light_m, eye_m):
    return 3.86 * (math.sqrt(light_m) + math.sqrt(eye_m)) / 1.852

BASE_CD = 2000.0          # assumed starting intensity
for vis in (5, 10, 20):
    base = luminous_range_nm(BASE_CD, vis)
    new = luminous_range_nm(BASE_CD * 6, vis)
    print(vis, round(base, 1), round(new, 1), round(new / base, 2))

target = 2 * luminous_range_nm(BASE_CD, 10)   # double the nominal range
lo, hi = 1.0, 1e6
for _ in range(100):
    mid = math.sqrt(lo * hi)
    if luminous_range_nm(BASE_CD * mid, 10) < target:
        lo = mid
    else:
        hi = mid
print(round(hi, 1), "x the intensity doubles the range")

# Output when run:
# 5 6.9 9.0 1.31
# 10 10.8 14.7 1.36
# 20 16.1 23.2 1.44
# 100.6 x the intensity doubles the range
#
# Each of the first three lines gives: visibility, range of the 2,000 cd light, range of the
# 12,000 cd light, and the ratio. Six times the intensity gave 1.31 times the range in poor
# visibility (5 nm), 1.36 times at 10 nm and 1.44 times in clear air (20 nm). The last line says
# the intensity has to be multiplied by about 100.6 to double the 10 nm range. The inverse square
# law alone would need only 4 times (range goes as the square root of intensity), so the rest of
# the factor is the light lost in the air, which grows with the distance travelled.