#!/usr/bin/env python3
"""Reproduce the fictional migration/service example with the Python standard library.
All positions, travel times, capacities, and arrivals are synthetic.
Coordinates are in a fictional local kilometre grid, not a geographic CRS.
Run from anywhere; outputs are written relative to this script's package.
"""
from __future__ import annotations
import heapq
import json
import math
from html import escape
from pathlib import Path
CONTENT = Path(__file__).resolve().parent
DATA = CONTENT / 'data'
ASSETS = CONTENT / 'assets'
def shortest_times(nodes, edges, source, season):
"""Undirected shortest paths in minutes; disconnected nodes remain None."""
graph = {n['id']: [] for n in nodes}
for e in edges:
weight = e[season + '_minutes']
if weight is None:
continue
if weight < 0:
raise ValueError('Travel times must be non-negative')
graph[e['from']].append((e['to'], weight))
graph[e['to']].append((e['from'], weight))
dist = {k: math.inf for k in graph}
dist[source] = 0
queue = [(0, source)]
while queue:
value, node = heapq.heappop(queue)
if value != dist[node]:
continue
for neighbour, weight in graph[node]:
new = value + weight
if new < dist[neighbour]:
dist[neighbour] = new
heapq.heappush(queue, (new, neighbour))
return {k: None if not math.isfinite(v) else v for k, v in dist.items()}
def calculate(example):
params = example['parameters']
travel = {s: shortest_times(example['nodes'], example['edges'], 'clinic', s)
for s in ['dry', 'wet']}
rows = []
for node in example['nodes']:
if node['kind'] != 'settlement':
continue
base = node['population_initial'] * (1 + params['background_growth_fraction'])
q = params['water_litres_per_person_day']
capacity = node['water_capacity_litres_day']
background_gap = max(0, base * q - capacity)
for scenario in ['lower', 'central', 'higher']:
arrivals = node['additional_residents'][scenario]
pop = base + arrivals
demand = pop * q
gap = max(0, demand - capacity)
rows.append({
'id': node['id'], 'name': node['name'], 'scenario': scenario,
'population_without_project': round(base, 6),
'additional_residents': arrivals,
'population_with_project': round(pop, 6),
'demand_litres_day': round(demand, 6),
'capacity_litres_day': capacity,
'background_gap_litres_day': round(background_gap, 6),
'gap_litres_day': round(gap, 6),
'gap_increase_litres_day': round(gap - background_gap, 6),
'gap_percent_capacity': round(100 * gap / capacity, 6),
'dry_clinic_minutes': travel['dry'][node['id']],
'wet_clinic_minutes': travel['wet'][node['id']],
})
return rows
def svg_open(width, height, title, description):
return [f'']
return '\n'.join(out)
def deficit_chart(rows):
central = [r for r in rows if r['scenario'] == 'central']
out = svg_open(940, 412, 'Separate existing shortages from additional demand',
'Synthetic daily water deficits. Gate has a 14400 litre additional deficit. Junction has a 9600 litre additional deficit. Riverside has an existing 1920 litre deficit plus a 9600 litre increase, totalling 11520 litres. Hill has no deficit. Sum of local deficits is 35520 litres; pooled regional shortfall is 28800 litres.')
out += [text_svg(24,34,'Separate existing shortages from additional demand','heading'),
text_svg(24,59,'Fictional central scenario • daily water shortfall • units: litres per day','sub')]
for i,r in enumerate(central):
y=98+i*51
x0=138
scale=.042
existing=r['background_gap_litres_day']*scale
added=r['gap_increase_litres_day']*scale
out.append(text_svg(24,y+18,r['name'],'label'))
if existing:
out.append(f'')
if added:
out.append(f'')
out.append(text_svg(x0+existing+added+12,y+18,f"{r['gap_litres_day']:,.0f}",'value'))
out += [f'',text_svg(45,334,'Shortfall without project','small'),
f'',text_svg(289,334,'Increase in shortfall in central scenario','small'),
text_svg(24,366,'Local deficits total 35,520 L/day. Regional pooling gives 28,800 L/day, but assumes capacity can be transferred.','sub'),
text_svg(24,391,'All inputs are synthetic. This calculation does not establish water adequacy, quality, affordability, or engineering capacity.','small'),'']
return '\n'.join(out)
def main():
example=json.loads((DATA/'synthetic-example.json').read_text())
rows=calculate(example)
# Verify the published figures and a genuinely important spatial distinction.
central={r['id']:r for r in rows if r['scenario']=='central'}
assert central['gate']['gap_litres_day']==14400
assert central['junction']['gap_litres_day']==9600
assert central['riverside']['gap_litres_day']==11520
assert central['riverside']['background_gap_litres_day']==1920
assert central['hill']['wet_clinic_minutes']==66
assert central['riverside']['wet_clinic_minutes']==56
assert sum(r['gap_litres_day'] for r in central.values())==35520
pooled=max(0,sum(r['demand_litres_day'] for r in central.values())-sum(r['capacity_litres_day'] for r in central.values()))
assert pooled==28800
results={'synthetic':True,'rows':rows,'central_summary':{
'additional_residents':sum(r['additional_residents'] for r in central.values()),
'sum_local_gap_litres_day':35520,'pooled_gap_litres_day':pooled,
'note':'Pooled capacity is not assumed transferable in the local calculations.'}}
DATA.mkdir(exist_ok=True,parents=True); ASSETS.mkdir(exist_ok=True,parents=True)
(DATA/'synthetic-results.json').write_text(json.dumps(results,indent=2)+'\n')
(ASSETS/'spatial-priorities.svg').write_text(map_portfolio(example,rows))
(ASSETS/'water-shortfalls.svg').write_text(deficit_chart(rows))
print(json.dumps({'verified':True,'settlements':len(central),'scenario_records':len(rows),'central_summary':results['central_summary']},indent=2))
if __name__=='__main__':
main()