simulation.py

"""
Predator-Prey Simulation

Designed for client-server architecture:
- Agents have stable unique IDs
- Spatial indexing via grid cells for O(1) neighbor queries
- Dirty tracking for delta updates
- Clean separation of state vs logic
"""

import random
import math
import numpy as np
from dataclasses import dataclass, field
from enum import Enum
from typing import Dict, Set, List, Tuple, Optional
from collections import defaultdict
import logging

TWOPI = 2 * math.pi

class AgentType(Enum):
    PREY = 0
    PREDATOR = 1

@dataclass
class Agent:
    """Minimal agent state - everything needed to render + simulate"""
    id: int
    agent_type: AgentType
    x: float  # Float position (cell = int(x), int(y))
    y: float
    energy: float
    sex: int = 0            # 0 = male, 1 = female
    born: int = 0           # Tick when agent was spawned
    parent: List[int] = field(default_factory=lambda: [-1, -1])  # IDs of parents ([-1,-1] for founders)
    offspring: List[int] = field(default_factory=list)  # IDs of offspring
    
    # Genetic parameters (inherited with mutation)
    vision: float = 0.0      # Effective vision range (capped by type max)
    brutpflege: float = 0.0  # Maternal investment (extra energy to offspring)
    speed: float = 0.0       # Sprint speed modifier (added to sprint ratio)
    heat: float = 0.0        # Mating priority exponent (predators only), 'heat' misnomer, as lower = more mating focus
    aggression: float = 0.0  # Base attack probability (prey only for now)
    genes: List[float] = field(default_factory=lambda: [5.0, 5.0, 5.0])  # Speciation genes [0-10]
    
    # Movement state
    bearing: float = 0.0     # Current heading in radians [0, 2π)
    
    @property
    def cell(self) -> Tuple[int, int]:
        """Current cell for food/interactions (floor of position)"""
        return (int(self.x), int(self.y))
    
    # Class-level definition of phenotype fields (sent in streaming)
    # Everything else is "internal" (only sent on inspect)
    PHENOTYPE_FIELDS = ('id', 'type', 'x', 'y')
    
    @property
    def position(self) -> tuple:
        """
        Agent's position in the simulation's coordinate system.
        
        Returns an opaque position tuple - the server should not
        interpret the contents, only pass them through.
        For this 2D simulation, returns (x, y).
        """
        return (self.x, self.y)
    
    def to_display_dict(self) -> dict:
        """
        Minimal data for visualization - streamed continuously.
        Only includes phenotype fields for efficient bandwidth.
        Positions rounded to 2 decimal places.
        """
        return {
            'id': self.id,
            'type': self.agent_type.value,
            'sex': self.sex,
            'x': round(self.x, 2),
            'y': round(self.y, 2),
        }

    def to_viewport_dict(self) -> dict:
        """
        Intermediate data for agents inviewport.
        Positions rounded to 2 decimal places.
        """
        return {
            'id': self.id,
            'type': self.agent_type.value,
            'sex': self.sex,
            'x': round(self.x, 2),
            'y': round(self.y, 2),
            'e': round(self.energy, 2),
            'gn': [round(g, 2) for g in self.genes],
        }

    
    def to_full_dict(self) -> dict:
        """
        Complete state - only sent on inspect request.
        Includes phenotype + internal state + genetics.
        Positions rounded to 2 decimal places, bearing (in radians) to 3 places.
	name conventions to save bandwidth (but imposes some load on client to translate back)
        const renameMap = {
	  a: 'bearing',
          e: 'energy',
          b: 'born',
          p: 'parent',
          o: 'offspring',
          v: 'vision',
          pf: 'brutpflege',
          s: 'speed',
          h: 'heat',
          ag: 'aggression'
          gn: 'genes',          
        };

        """
        return {
            'id': self.id,
            'type': self.agent_type.value,
            'sex': self.sex,
            'x': round(self.x, 2),
            'y': round(self.y, 2),
            #above: streamed values, below: "internal state" only sent on request
            'a': round(self.bearing, 3),
            'e': round(self.energy, 2),
            'b': self.born,
            'p': self.parent.copy(),  # Copy to avoid mutation
            'o': self.offspring.copy(),  # Copy to avoid mutation
            'v': self.vision,
            'pf': self.brutpflege,
            's': self.speed,
            'h': self.heat,
            'ag': self.aggression,
            'gn': [round(g, 2) for g in self.genes],
        }
    
    def to_dict(self) -> dict:
        """Serializable representation (full state) - backward compatible"""
        return self.to_full_dict()


class SpatialGrid:
    """
    Spatial index for O(1) queries of "which agents are in cell (x,y)?"
    
    Critical for scalability: without this, finding nearby agents is O(n).
    With this, it's O(agents_in_cell).
    
    Accepts float positions - converts to int cell coordinates internally.
    
    Maintains separate indices by agent type for efficient type-specific queries.
    """
    
    def __init__(self, width: int, height: int):
        self.width = width
        self.height = height
        # cell -> set of agent IDs in that cell (all types)
        self._grid: Dict[Tuple[int, int], Set[int]] = defaultdict(set)
        # Type-specific grids: type -> cell -> set of agent IDs
        self._type_grids: Dict[int, Dict[Tuple[int, int], Set[int]]] = {
            0: defaultdict(set),  # PREY
            1: defaultdict(set),  # PREDATOR
        }
        # Cells that contain each type (for fast iteration)
        self._type_cells: Dict[int, Set[Tuple[int, int]]] = {
            0: set(),  # cells with PREY
            1: set(),  # cells with PREDATOR
        }
    
    def add(self, agent_id: int, x: float, y: float, agent_type: int = 0):
        cell = (int(x), int(y))
        self._grid[cell].add(agent_id)
        self._type_grids[agent_type][cell].add(agent_id)
        self._type_cells[agent_type].add(cell)
    
    def remove(self, agent_id: int, x: float, y: float, agent_type: int = 0):
        cell = (int(x), int(y))
        self._grid[cell].discard(agent_id)
        type_cell_set = self._type_grids[agent_type][cell]
        type_cell_set.discard(agent_id)
        if not type_cell_set:
            self._type_cells[agent_type].discard(cell)
    
    def move(self, agent_id: int, old_x: float, old_y: float, new_x: float, new_y: float, agent_type: int = 0):
        """Update agent's cell assignment. Only modifies grid if cell changed."""
        old_cell = (int(old_x), int(old_y))
        new_cell = (int(new_x), int(new_y))
        if old_cell != new_cell:
            self._grid[old_cell].discard(agent_id)
            self._grid[new_cell].add(agent_id)
            # Update type-specific grid
            old_type_set = self._type_grids[agent_type][old_cell]
            old_type_set.discard(agent_id)
            if not old_type_set:
                self._type_cells[agent_type].discard(old_cell)
            self._type_grids[agent_type][new_cell].add(agent_id)
            self._type_cells[agent_type].add(new_cell)
    
    def agents_at(self, x: int, y: int) -> Set[int]:
        return self._grid[(x, y)]
    
    def agents_in_region(self, x1: int, y1: int, x2: int, y2: int) -> Set[int]:
        """Get all agent IDs in a rectangular region (for AOI queries)"""
        result = set()
        for x in range(x1, x2 + 1):
            for y in range(y1, y2 + 1):
                result.update(self._grid[(x, y)])
        return result
    
    def agents_of_type_in_region(self, x1: int, y1: int, x2: int, y2: int, agent_type: int) -> Set[int]:
        """
        Get agent IDs of specific type in a rectangular region.
        
        Optimized: only checks cells that actually contain agents of this type.
        For sparse types (like predators), this is MUCH faster than agents_in_region.
        """
        result = set()
        type_cells = self._type_cells[agent_type]
        
        # Only check cells that have this type
        for cell in type_cells:
            cx, cy = cell
            if x1 <= cx <= x2 and y1 <= cy <= y2:
                result.update(self._type_grids[agent_type][cell])
        
        return result


class World:
    """
    The simulation world state and logic.
    
    Scalability considerations:
    - Food is stored per-cell (not per-agent), so memory is O(cells) not O(agents)
    - Spatial grid gives O(1) cell lookups
    - Dirty set tracks which agents changed for delta updates
    
    Biomes:
    - 0: baseline (no modifiers)
    - 1: desert (no seasons, reduced food, 1.2x metabolism/movement)
    - 2: arctic (extreme seasons, reduced food, 1.5x metabolism/movement)
    - 3: tropical (no seasons, bonus food)
    - 7: wall (impassable, no food) - reserved for future
    """
    
    # Biome constants
    BIOME_BASELINE = 0
    BIOME_DESERT = 1
    BIOME_ARCTIC = 2
    BIOME_TROPICAL = 3
    BIOME_WALL = 7
    
    # Biome modifiers: {biome: (food_mult, cost_mult, season_mult)}
    # food_mult: multiplier for food regen and ceiling
    # cost_mult: multiplier for metabolism and movement costs
    # season_mult: multiplier for season effect (0 = no seasons, >1 = extreme seasons)
    BIOME_MODIFIERS = {
        0: (1.0, 1.0, 1.0),   # baseline
        1: (0.75, 1.1, 0.25),   # desert: less food, higher cost, weak seasons
        2: (0.75, 1.1, 1.25),   # arctic: less food, higher cost, extreme seasons
        3: (1.25, 1.0, 0.0),   # tropical: more food, normal cost, no seasons
        7: (0.0, 999, 0.0),   # wall: no food, impassable
    }

    # Foraging logic
    BEARING_ADJUST = 0.3  # 0.1 rad ~6 degrees, for prey (adjust towards more attractive pasture)
    PREY_BEARING_DRIFT_STDDEV = 0.5
    PREDATOR_BEARING_DRIFT_STDDEV = 0.1

    _QUADRANT_CELLS = (
        ((1, 0), (1, 1), (1, -1)),   # Q0: E, NE, SE
        ((0, 1), (1, 1), (-1, 1)),   # Q1: N, NE, NW
        ((-1, 0), (-1, 1), (-1, -1)), # Q2: W, NW, SW
        ((0, -1), (1, -1), (-1, -1)), # Q3: S, SE, SW
    )
    
    # Simulation parameters
    # Tuned for stable oscillating predator-prey dynamics
    FOOD_REGEN_PER_TURN = 0.2
    FOOD_CEILING = 4.0
    SEASON_LENGTH = 40
    SEASON_STRENGTH = 0.5 # in [0,1]
    
    PREY_INITIAL_ENERGY = 2.0
    PREY_MOVE_COST = 0.02        # Base cost for normal movement (foraging)
    PREY_MOVE_SPEED = 0.1      # Base speed for normal movement
    PREY_SPRINT_RATIO = 4.0     # Speed multiplier when fleeing
    PREY_SPRINT_COST = 6.0      # Cost multiplier when fleeing (actual = MOVE_COST * SPRINT_COST)
    PREY_METABOLISM_COST = 0.05
    PREY_REPRODUCTION_THRESHOLD = 2*PREY_INITIAL_ENERGY
    PREY_REPRODUCTION_COST = PREY_INITIAL_ENERGY
    PREY_EAT_RATE = 0.8  # max food eaten per turn
    PREY_VISION_RANGE = 5  # max vision range for prey
    PREY_VISION_COST = 0.01
    PREY_DEFAULT_VISION = 0.3  # founder population vision range (0: start blind)
    PREY_SENESCENCE = 160
    PREY_AGE_PENALTY = 0.1
    PREY_MAX_ENERGY = 10.0      # Stop foraging above this (satiated)
    PREY_MALE_INVESTMENT = 0.3  # Male pays this fraction of reproduction cost
    
    PREDATOR_INITIAL_ENERGY = 6.0
    PREDATOR_MOVE_COST = 0.01        # Base cost for normal movement (patrolling)
    PREDATOR_MOVE_SPEED = 0.1      # Base speed for normal movement
    PREDATOR_SPRINT_RATIO = 4.5     # Speed multiplier when pursuing
    PREDATOR_SPRINT_COST = 5.0      # Cost multiplier when pursuing
    PREDATOR_METABOLISM_COST = 0.01
    PREDATOR_REPRODUCTION_THRESHOLD = 2*PREDATOR_INITIAL_ENERGY
    PREDATOR_REPRODUCTION_COST = PREDATOR_INITIAL_ENERGY*2
    PREDATOR_ENERGY_FROM_PREY = 3.5  # fixed value, or have it depend on prey energy?
    PREDATOR_VISION_RANGE = 5  # max vision range for predators
    PREDATOR_VISION_COST = 0.01
    PREDATOR_DEFAULT_VISION = 3  # founder population vision range
    PREDATOR_SENESCENCE = 500
    PREDATOR_AGE_PENALTY = 0.1
    PREDATOR_MAX_ENERGY = 30.0      # Stop hunting above this (satiated)
    PREDATOR_MALE_INVESTMENT = 0.48  # Male pays this fraction of reproduction cost
    PREDATOR_HEAT = 10.0            # Baseline mating priority exponent (lower = more mating focus)
    
    MATING_DISTANCE = 1.1  # Max distance for mating (needs vision to see mates outside own cell)
    CATCH_DIST_SQ = 0.15 * 0.15 # Predator catches prey within this distance (square)
    
    # Speciation and combat (prey only for now)
    PREY_SPECIATION_THRESHOLD_SQ = 1.0  # Genetic distance^2 beyond which agents are different species
    PREY_GENE_MUTATION_SD = 0.4        # Mutation rate for speciation genes
    PREY_AGGRESSION_MUTATION_SD = 0.2  # Mutation rate for aggression gene
    PREY_ATTACK_COST = 1.0              # Energy spent to attack
    PREY_ATTACK_DAMAGE = 3.0            # Energy dealt to target
    PREY_INTRASPECIES_MODIFIER = 0.3	# modifies attack probability/cost/damage for intraspecies aggression
    PREY_FEMALE_AGGRESSION_MULT = 0.5   # Females are less aggressive
    PREY_DEFAULT_AGGRESSION = 0.3       # Founders start non-aggressive
    PREY_ATTACK_REPULSION = 0.15        # Distance attacker is pushed back after interspecies attack
    CONSENSUS_STRENGTH = 0.8		# foraging prey herd behavior
    
    # Default board size and founder populations
    DEFAULT_WIDTH = 100
    DEFAULT_HEIGHT = 100
    DEFAULT_INITIAL_PREY = 1000
    DEFAULT_INITIAL_PREDATORS = 15
    
    # Statistics logging
    STATS_LOG_INTERVAL = 10  # Print stats every N ticks (0 to disable)
    HIST_BINS = 10 #number of bins in histograms

    def __init__(self, width: int, height: int, seed: Optional[int] = None):
        self.width = width
        self.height = height
        self.tick = 0
        
        if seed is not None:
            random.seed(seed)
        
        # Biome per cell (generated after construction via generate_biomes)
        self.biome: Dict[Tuple[int, int], int] = {}
        for x in range(width):
            for y in range(height):
                self.biome[(x, y)] = self.BIOME_BASELINE
        
        # Food per cell (initialized after biomes are set)
        self.food: Dict[Tuple[int, int], float] = {}
        for x in range(width):
            for y in range(height):
                food_mult, _, _ = self.BIOME_MODIFIERS[self.biome[(x, y)]]
                self.food[(x, y)] = random.uniform(0, self.FOOD_CEILING * food_mult)
        
        # Agent management
        self._next_id = 0
        self.agents: Dict[int, Agent] = {}
        self.spatial_grid = SpatialGrid(width, height)
        
        # Delta tracking: which agents have changed since last "mark_clean()"
        self._dirty_agents: Set[int] = set()
        self._spawned_agents: Set[int] = set()  # new since last clean
        self._despawned_agents: Set[int] = set()  # removed since last clean
        self._mated_this_tick: Set[int] = set()  # agents who have mated this tick
        self._fled_this_tick: Set[int] = set()   # agents who fled this tick (skip combat)
        
        # Per-tick statistics (reset each tick, accumulated for logging)
        self._stats = {
            'attacks_interspecies': 0,
            'attacks_intraspecies': 0,
            'attack_kills': 0,
            'matings': 0,
            'deaths_predation': 0,
            'deaths_starvation': 0,
        }

        #statistics histograms, set default maxima
        self._hist_max = {
            'prey_age': 100,      # initial defaults
            'prey_energy': self.PREY_MAX_ENERGY,
            'prey_speed': 1.0,
            'prey_vision': self.PREY_VISION_RANGE,
            'prey_brutpflege': 1.0,
            'prey_aggression': -1,  # sentinel value -1 means "it's a percentage", max fixed at 1.0
            'predator_age': 100,
            'predator_energy': self.PREDATOR_MAX_ENERGY,
            'predator_speed': 1.0,
            'predator_vision': self.PREDATOR_VISION_RANGE,
            'predator_brutpflege': 1.0,
            'predator_heat': 10.0,
        }
        self._histograms = {}  # filled by _update_histograms
        #how to report each property as histogram value
        self.hist_value_getters = {
            'age': lambda a: self.tick - a.born,
            'energy': lambda a: a.energy,
            'vision': lambda a: a.vision,
            'brutpflege': lambda a: a.brutpflege,
            'speed': lambda a: a.speed + (self.PREY_SPRINT_RATIO - 1 if a.agent_type == AgentType.PREY else self.PREDATOR_SPRINT_RATIO - 1),
            'aggression': lambda a: a.aggression,
            'heat': lambda a: a.heat,
        }
        
        # Simulation state
        self.halted = False  # Set to True when termination condition is met
        self.halt_reason: Optional[str] = None  # Why simulation halted
        
        # Deferred command execution (for commands received during step())
        self._stepping = False  # True while inside step()
        self._deferred_commands: List[Tuple[str, any]] = []  # Queue of (name, value) to execute after step
    
    def generate_biomes(self):
        """
        Generate biome map based on latitude (y coordinate).
        
        Default layout for 100-height world:
        - y in [0, 19]: arctic (2)
        - y in [20, 49]: baseline (0) 
        - y in [50, 69]: desert (1)
        - y in [70, 99]: tropical (3)
        
        Scales proportionally for other world heights.
        
        Additionally, longitudinal walls are placed at:
        - x % 30 == 10 and y % 25 < 10
        """
        for x in range(self.width):
            for y in range(self.height):
                # Check for wall first (longitudinal barriers)
                if x % 30 == 10 and y % 25 < 10:
                    biome = self.BIOME_WALL
                else:
                    # Latitude-based biome assignment
                    y_norm = y / self.height
                    
                    if y_norm < 0.2:
                        biome = self.BIOME_ARCTIC
                    elif y_norm < 0.5:
                        biome = self.BIOME_BASELINE
                    elif y_norm < 0.7:
                        biome = self.BIOME_DESERT
                    else:
                        biome = self.BIOME_TROPICAL
                
                self.biome[(x, y)] = biome
        
        # Reinitialize food based on new biomes
        self._init_food_from_biomes()
    
    def load_biome_map(self, filename: str):
        """
        Load biome map from a text file.
        
        File format:
        - Each line is a row (line 0 = y=0)
        - Each character is a column (char 0 = x=0)
        - Characters '0'-'7' map to biome codes
        
        The world dimensions are set from the map file.
        
        Args:
            filename: Path to the map file
            
        Raises:
            ValueError: If map dimensions don't match or invalid characters found
        """
        with open(filename, 'r') as f:
            lines = [line.rstrip('\n\r') for line in f.readlines()]
        
        # Filter out empty lines
        lines = [line for line in lines if line]
        
        if not lines:
            raise ValueError("Empty map file")
        
        map_height = len(lines)
        map_width = len(lines[0])
        
        # Verify all lines have same width
        for i, line in enumerate(lines):
            if len(line) != map_width:
                raise ValueError(f"Line {i} has width {len(line)}, expected {map_width}")
        
        # Check if dimensions match world (resize if needed)
        if map_width != self.width or map_height != self.height:
            # Resize world to match map
            self.width = map_width
            self.height = map_height
            self.spatial_grid = SpatialGrid(map_width, map_height)
            # Reinitialize biome and food dicts
            self.biome = {}
            self.food = {}
            for x in range(map_width):
                for y in range(map_height):
                    self.biome[(x, y)] = self.BIOME_BASELINE
                    self.food[(x, y)] = 0.0
        
        # Parse map
        for y, line in enumerate(lines):
            for x, char in enumerate(line):
                if char not in '01234567':
                    raise ValueError(f"Invalid character '{char}' at ({x}, {y})")
                self.biome[(x, y)] = int(char)
        
        # Initialize food based on biomes
        self._init_food_from_biomes()
    
    def _init_food_from_biomes(self):
        """Initialize food values based on biome modifiers."""
        for x in range(self.width):
            for y in range(self.height):
                food_mult, _, _ = self.BIOME_MODIFIERS[self.biome[(x, y)]]
                self.food[(x, y)] = random.uniform(0, self.FOOD_CEILING * food_mult)
    
    def is_wall(self, x: int, y: int) -> bool:
        """Check if a cell is a wall (impassable)."""
        return self.biome.get((x, y), self.BIOME_BASELINE) == self.BIOME_WALL
    
    def is_passable(self, x: float, y: float) -> bool:
        """Check if a position is passable (not a wall, within bounds)."""
        ix, iy = int(x), int(y)
        if not (0 <= ix < self.width and 0 <= iy < self.height):
            return False
        return self.biome[(ix, iy)] != self.BIOME_WALL
    
    def get_biome_modifier(self, x: float, y: float) -> Tuple[float, float, float]:
        """
        Get biome modifiers for a position.
        
        Returns (food_mult, cost_mult, season_mult) for the cell containing (x, y).
        """
        cell = (int(x), int(y))
        biome = self.biome.get(cell, self.BIOME_BASELINE)
        return self.BIOME_MODIFIERS.get(biome, self.BIOME_MODIFIERS[self.BIOME_BASELINE])
    
    @classmethod
    def create(cls, **kwargs) -> 'World':
        """
        Create and initialize a world with agents.
        
        This is the main entry point for creating a simulation world.
        The server passes an opaque config dict - only this method knows
        what parameters are valid and what they mean.
        
        Args (all optional, with defaults):
            width: World width in cells (default: 100, ignored if map_file provided)
            height: World height in cells (default: 100, ignored if map_file provided)
            initial_prey: Number of prey to spawn (default: 1000)
            initial_predators: Number of predators to spawn (default: 15)
            seed: Random seed for reproducibility (default: None)
            biomes: Whether to generate biomes (default: True, ignored if map_file provided)
            map_file: Path to biome map file (default: None)
        
        Returns:
            Initialized World with agents spawned and dirty state cleared
        """
        # Extract and apply defaults - only this method knows the parameter names
        width = kwargs.get('width', cls.DEFAULT_WIDTH)
        height = kwargs.get('height', cls.DEFAULT_HEIGHT)
        initial_prey = kwargs.get('initial_prey', cls.DEFAULT_INITIAL_PREY)
        initial_predators = kwargs.get('initial_predators', cls.DEFAULT_INITIAL_PREDATORS)
        seed = kwargs.get('seed', None)
        biomes = kwargs.get('biomes', True)
        map_file = kwargs.get('map_file', 'maps/map2.txt') # 'maps/map1.txt'
        
        # Create the world
        world = cls(width, height, seed=seed)
        
        # Load biome map from file, or generate procedurally
        if map_file:
            world.load_biome_map(map_file)
            # Update width/height in case map resized the world
            width = world.width
            height = world.height
        elif biomes:
            world.generate_biomes()
        
        # Spawn initial prey at random float positions (retry if on wall)
        for _ in range(initial_prey):
            for attempt in range(100):  # Max retries to avoid infinite loop
                x = random.random() * width
                y = random.random() * height
                if world.is_passable(x, y):
                    world.spawn_agent(AgentType.PREY, x, y, cls.PREY_INITIAL_ENERGY)
                    break
        
        # Spawn initial predators as breeding pairs (male + female per cell)
        # initial_predators specifies number of pairs, so total = initial_predators * 2
        for _ in range(initial_predators):
            for attempt in range(100):
                # Pick a random cell
                cell_x = random.randint(0, width - 1)
                cell_y = random.randint(0, height - 1)
                if world.is_passable(cell_x + 0.5, cell_y + 0.5):
                    # Spawn male at random position in cell
                    male_x = cell_x + random.random()
                    male_y = cell_y + random.random()
                    male = world.spawn_agent(AgentType.PREDATOR, male_x, male_y, cls.PREDATOR_INITIAL_ENERGY)
                    male.sex = 0  # Force male
                    
                    # Spawn female at different random position in same cell
                    female_x = cell_x + random.random()
                    female_y = cell_y + random.random()
                    female = world.spawn_agent(AgentType.PREDATOR, female_x, female_y, cls.PREDATOR_INITIAL_ENERGY)
                    female.sex = 1  # Force female
                    break
        
        # Clear dirty state so first delta is clean
        world.mark_clean()
        
        return world
    
    def _allocate_id(self) -> int:
        """Get a unique agent ID. IDs are never reused (important for client sync)."""
        agent_id = self._next_id
        self._next_id += 1
        return agent_id
    
    def spawn_agent(self, agent_type: AgentType, x: float, y: float, energy: float, 
                    parent_ids: Tuple[int, int] = (-1, -1)) -> Agent:
        """
        Create a new agent in the world.
        
        Args:
            agent_type: Type of agent (PREY or PREDATOR)
            x, y: Position (float coordinates)
            energy: Initial energy
            parent_ids: Tuple of parent IDs ((-1, -1) for founders)
        
        Returns:
            The created Agent
        
        Genetic inheritance:
            - Founders (parent_ids == (-1, -1)) get type defaults
            - Offspring inherit average of both parents' genes + Gaussian noise (sd=0.2)
            - Vision and brutpflege clamped >= 0
            - Speed clamped to minimum that keeps sprint ratio >= 1
        
        Sex is assigned randomly (50/50).
        """
        agent_id = self._allocate_id()
        
        # Determine minimum speed based on agent type
        if agent_type == AgentType.PREY:
            min_speed = -(self.PREY_SPRINT_RATIO - 1)
        else:
            min_speed = -(self.PREDATOR_SPRINT_RATIO - 1)
        
        # Determine genetic parameters
        parent1_id, parent2_id = parent_ids
        parent1 = self.agents.get(parent1_id) if parent1_id >= 0 else None
        parent2 = self.agents.get(parent2_id) if parent2_id >= 0 else None
        
        if parent1 and parent2:
            # Sexual reproduction: average parents' genes, then add mutation
            avg_vision = (parent1.vision + parent2.vision) / 2
            avg_brutpflege = (parent1.brutpflege + parent2.brutpflege) / 2
            avg_speed = (parent1.speed + parent2.speed) / 2
            
            vision = max(0.0, avg_vision + random.gauss(0, 0.2))
            brutpflege = max(0.0, avg_brutpflege + random.gauss(0, 0.2))
            speed = max(min_speed, avg_speed + random.gauss(0, 0.2))
            
            # Heat: only mutate for predators, prey always 0
            if agent_type == AgentType.PREDATOR:
                avg_heat = (parent1.heat + parent2.heat) / 2
                # Heat mutates multiplicatively, clamped to [0.005, 2]
                heat = max(0.01, min(200.0, avg_heat * (1 + random.gauss(0, 0.1))))
                aggression = 0  # Predators don't use aggression (yet)
                genes = []  # Predators don't use speciation genes (yet)
            else:
                heat = 0
                # Aggression: only for prey, average + mutation, clamped to [0, 1]
                avg_aggression = (parent1.aggression + parent2.aggression) / 2
                aggression = max(0.0, min(1.0, avg_aggression + random.gauss(0, self.PREY_AGGRESSION_MUTATION_SD)))
                # Speciation genes: average + mutation, clamped to [0, 10]
                genes = []
                for i in range(3):
                    avg_gene = (parent1.genes[i] + parent2.genes[i]) / 2
                    genes.append(max(0.0, min(10.0, avg_gene + random.gauss(0, self.PREY_GENE_MUTATION_SD))))
        else:
            # Founders get type-specific defaults
            if agent_type == AgentType.PREY:
                vision = self.PREY_DEFAULT_VISION
                heat = 0  # Prey don't use heat
                aggression = self.PREY_DEFAULT_AGGRESSION
		# species by map quarter (pending more versatile initialization options)
                genes = [5.0+3*(1 if x<50 else -1), 5.0+3*(1 if y<50 else -1), 5.0+2*(1 if x<50 else -1)]  
            else:
                vision = self.PREDATOR_DEFAULT_VISION
                heat = self.PREDATOR_HEAT  # Predators start at baseline heat
                aggression = 0  # Not used for predators yet
                genes = []  # Not used for predators yet
            brutpflege = 0  # founders start at brutpflege=0
            speed = 0       # founders start at baseline speed
        
        # Random initial bearing and sex
        bearing = (random.random() - 0.5) * TWOPI
        sex = random.randint(0, 1)  # 0 = male, 1 = female
        
        agent = Agent(
            id=agent_id,
            agent_type=agent_type,
            x=float(x),
            y=float(y),
            energy=energy,
            sex=sex,
            born=self.tick,
            parent=list(parent_ids),
            offspring=[],
            vision=vision,
            brutpflege=brutpflege,
            speed=speed,
            heat=heat,
            aggression=aggression,
            genes=genes,
            bearing=bearing
        )
        self.agents[agent_id] = agent
        self.spatial_grid.add(agent_id, x, y, agent_type.value)
        self._spawned_agents.add(agent_id)
        self._dirty_agents.add(agent_id)
        
        # Record this agent as offspring of both parents
        if parent1:
            parent1.offspring.append(agent_id)
        if parent2:
            parent2.offspring.append(agent_id)
        
        return agent
    
    def remove_agent(self, agent_id: int):
        """Remove an agent from the world."""
        if agent_id in self.agents:
            agent = self.agents[agent_id]
            self.spatial_grid.remove(agent_id, agent.x, agent.y, agent.agent_type.value)
            del self.agents[agent_id]
            self._despawned_agents.add(agent_id)
            self._dirty_agents.discard(agent_id)
            self._spawned_agents.discard(agent_id)
    
    def _mark_dirty(self, agent_id: int):
        """Mark an agent as modified (for delta updates)."""
        self._dirty_agents.add(agent_id)
    
    def get_dirty_state(self) -> dict:
        """
        Get all changes since last mark_clean().
        This is what we'll send to clients as delta updates.
        
        Uses to_display_dict() for phenotype-only data (efficient streaming).
        Full agent data is available via inspect command.
        """
        return {
            'tick': self.tick,
            'spawned': [self.agents[aid].to_display_dict() for aid in self._spawned_agents if aid in self.agents],
            'updated': [self.agents[aid].to_display_dict() for aid in self._dirty_agents - self._spawned_agents if aid in self.agents],
            'despawned': list(self._despawned_agents)
        }
    
    def get_dirty_ids(self) -> dict:
        """
        Get just the IDs of changed agents (no data copying).
        Much more efficient for streaming - actual data fetched lazily at send time.
        """
        return {
            'tick': self.tick,
            'spawned_ids': list(self._spawned_agents),
            'updated_ids': list(self._dirty_agents - self._spawned_agents),
            'despawned_ids': list(self._despawned_agents)
        }
    
    def get_full_state(self) -> dict:
        """
        Get complete world state (for initial sync or reconnection).
        
        Map cells are returned as entities alongside agents, with type=-1 (CELL).
        This unified approach allows the same streaming infrastructure to handle
        both agents and map data.
        """
        # Agents as entities (display dict only for streaming efficiency)
        agent_entities = [a.to_display_dict() for a in self.agents.values()]
        
        # Map cells as entities (type=-1 indicates CELL)
        cell_entities = [
            {
                'id': -(x * self.height + y + 1),  # Negative IDs for cells
                'type': -1,  # CELL type
                'x': x,
                'y': y,
                'food': round(self.food[(x, y)], 2),
                'biome': self.biome[(x, y)]
            }
            for x in range(self.width)
            for y in range(self.height)
        ]
        
        return {
            'tick': self.tick,
            'width': self.width,
            'height': self.height,
            'agents': agent_entities,
            'cells': cell_entities,
            'total_entities': len(agent_entities) + len(cell_entities)
        }
    
    def describe_map(self) -> List[List[int]]:
        """
        Get static map data for initial transmission to client.
        
        Returns only cells with non-zero biome codes (client assumes biome 0 by default).
        Each entry is [x, y, biome].
        
        Returns:
            List of [x, y, biome] for all cells where biome != 0
        """
        return [
            [x, y, biome]
            for (x, y), biome in self.biome.items()
            if biome != 0
        ]
    
    def mark_clean(self):
        """Clear dirty tracking (call after sending delta to clients)."""
        self._dirty_agents.clear()
        self._spawned_agents.clear()
        self._despawned_agents.clear()
    
    def inspect_agent(self, agent_id: int) -> Optional[dict]:
        """
        Get full state of a specific agent (for inspect command).
        Returns None if agent doesn't exist.
        """
        if agent_id in self.agents:
            return self.agents[agent_id].to_full_dict()
        return None
    
    @staticmethod
    def distance(p1, p2) -> float:
        """
        Distance between two positions in the simulation's topology.
        
        For this 2D grid simulation, uses Chebyshev distance (L∞ norm).
        Other simulations might use Euclidean, Manhattan, graph distance, etc.
        
        Args:
            p1: First position (tuple/list of coordinates)
            p2: Second position (tuple/list of coordinates)
        
        Returns:
            Distance as a float
        """
        return max(abs(a - b) for a, b in zip(p1, p2))
    
    @staticmethod
    def is_in_viewport(position, center, radius) -> bool:
        """
        Check if a position is within a viewport.
        
        This is simulation-specific - defines what "within radius" means
        for this simulation's topology.
        
        Args:
            position: Agent's position (tuple/list)
            center: Viewport center (tuple/list)
            radius: Viewport radius (float)
        
        Returns:
            True if position is within viewport
        """
        return World.distance(position, center) <= radius
    
    def get_viewport(self, center, radius) -> dict:
        """
        Get all entities (agents and map cells) within a viewport region.
        
        This is the simulation's implementation of viewport queries.
        The server calls this without knowing the topology details.
        
        Args:
            center: Viewport center position (tuple/list of coordinates)
            radius: Viewport radius
        
        Returns:
            {
                'tick': current tick,
                'center': center position as list,
                'radius': radius,
                'agents': list of entity dicts (agents + cells),
                'count': number of entities
            }
        
        Entity types:
            type >= 0: Agents (prey=0, predator=1, etc.)
            type = -1: Map cells (food)
        """
        entities = []
        
        # Add agents in viewport
        for agent in self.agents.values():
            if self.is_in_viewport(agent.position, center, radius):
                entities.append(agent.to_viewport_dict())
        
        # Add map cells in viewport
        for (x, y), food in self.food.items():
            if self.is_in_viewport((x, y), center, radius):
                entities.append({
                    'id': -(x * self.height + y + 1),  # Negative IDs for cells
                    'type': -1,  # CELL type
                    'x': x,
                    'y': y,
                    'food': round(food, 2)
                })
        
        return {
            'tick': self.tick,
            'center': list(center),
            'radius': radius,
            'agents': entities,  # Named 'agents' for backward compatibility
            'count': len(entities)
        }
    
    def _move_agent(self, agent: Agent, new_x: float, new_y: float) -> bool:
        """
        Move an agent to a new position (float coordinates).
        
        Returns:
            True if move succeeded, False if blocked (wall or out of bounds)
        """
        if new_x == agent.x and new_y == agent.y:
            return True
        
        # Clamp to world bounds (keep slightly inside to avoid edge issues)
        new_x = max(0.0, min(float(self.width) - 0.001, new_x))
        new_y = max(0.0, min(float(self.height) - 0.001, new_y))
        
        # Check if destination cell is a wall
        if not self.is_passable(new_x, new_y):
            return False
        
        self.spatial_grid.move(agent.id, agent.x, agent.y, new_x, new_y, agent.agent_type.value)
        agent.x = new_x
        agent.y = new_y
        self._mark_dirty(agent.id)
        return True
    
    def _move_toward(self, agent: Agent, target_x: float, target_y: float, speed: float) -> float:
        """
        Move agent toward target position at given speed.
        
        Args:
            agent: Agent to move
            target_x, target_y: Target position
            speed: Maximum distance to move
        
        Returns:
            Actual distance moved (0 if already at target)
        """
        dx = target_x - agent.x
        dy = target_y - agent.y
        dist = math.sqrt(dx * dx + dy * dy)
        
        if dist < 0.001:  # Already there
            return 0.0
        
        if dist <= speed:
            # Can reach target this tick
            new_x, new_y = target_x, target_y
            actual_dist = dist
        else:
            # Move speed units toward target
            new_x = agent.x + speed * dx / dist
            new_y = agent.y + speed * dy / dist
            actual_dist = speed
        
        self._move_agent(agent, new_x, new_y)
        return actual_dist
    
    def _move_away_from(self, agent: Agent, threat_x: float, threat_y: float, speed: float) -> float:
        """
        Move agent away from threat position at given speed.
        
        Args:
            agent: Agent to move
            threat_x, threat_y: Position to flee from
            speed: Distance to move
        
        Returns:
            Actual distance moved
        """
        dx = agent.x - threat_x
        dy = agent.y - threat_y
        dist = math.sqrt(dx * dx + dy * dy)
        
        if dist < 0.001:
            # Exactly on threat - pick random direction
            angle = (random.random() - 0.5) * TWOPI
            dx, dy = math.cos(angle), math.sin(angle)
            dist = 1.0
        
        # Move speed units away from threat
        new_x = agent.x + speed * dx / dist
        new_y = agent.y + speed * dy / dist
        
        self._move_agent(agent, new_x, new_y)
        return speed
    
    def _get_neighbors(self, x: int, y: int) -> List[Tuple[int, int]]:
        """Get valid neighboring cells (8-directional), excluding walls."""
        neighbors = []
        for dx in [-1, 0, 1]:
            for dy in [-1, 0, 1]:
                if dx == 0 and dy == 0:
                    continue
                nx, ny = x + dx, y + dy
                if 0 <= nx < self.width and 0 <= ny < self.height:
                    # Exclude wall cells
                    if self.biome[(nx, ny)] != self.BIOME_WALL:
                        neighbors.append((nx, ny))
        random.shuffle(neighbors)
        return neighbors
    
    def _expand_ring(self, neighbor_ids: Set[int], cx: int, cy: int, r: int):
        """
        Add ALL agents in the outer ring at radius r (in-place).
        Ring consists of cells at Chebyshev distance exactly r from (cx, cy).
        """
        # Top row (full width): y = cy - r
        if cy - r >= 0:
            for x in range(max(0, cx - r), min(self.width, cx + r + 1)):
                neighbor_ids.update(self.spatial_grid.agents_at(x, cy - r))
        
        # Bottom row (full width): y = cy + r
        if cy + r < self.height:
            for x in range(max(0, cx - r), min(self.width, cx + r + 1)):
                neighbor_ids.update(self.spatial_grid.agents_at(x, cy + r))
        
        # Left column (excluding corners): x = cx - r
        if cx - r >= 0:
            for y in range(max(0, cy - r + 1), min(self.height, cy + r)):
                neighbor_ids.update(self.spatial_grid.agents_at(cx - r, y))
        
        # Right column (excluding corners): x = cx + r
        if cx + r < self.width:
            for y in range(max(0, cy - r + 1), min(self.height, cy + r)):
                neighbor_ids.update(self.spatial_grid.agents_at(cx + r, y))

    def _expand_neighborhood_ring(self, neighbor_ids: Set[int], cx: int, cy: int, r: int):
        """
        Add agents in the outer ring at radius r to neighbor_ids (in-place).
        Ring consists of cells at Chebyshev distance exactly r from (cx, cy).
        """
        # Top row (full width): y = cy - r
        if cy - r >= 0:
            for x in range(max(0, cx - r), min(self.width, cx + r + 1)):
                neighbor_ids.update(self.spatial_grid.agents_at(x, cy - r))
        
        # Bottom row (full width): y = cy + r
        if cy + r < self.height:
            for x in range(max(0, cx - r), min(self.width, cx + r + 1)):
                neighbor_ids.update(self.spatial_grid.agents_at(x, cy + r))
        
        # Left column (excluding corners): x = cx - r
        if cx - r >= 0:
            for y in range(max(0, cy - r + 1), min(self.height, cy + r)):
                neighbor_ids.update(self.spatial_grid.agents_at(cx - r, y))
        
        # Right column (excluding corners): x = cx + r
        if cx + r < self.width:
            for y in range(max(0, cy - r + 1), min(self.height, cy + r)):
                neighbor_ids.update(self.spatial_grid.agents_at(cx + r, y))
    
    def _expand_type_ring(self, neighbor_ids: Set[int], cx: int, cy: int, r: int, agent_type: int):
        """
        Add agents of specific type in the outer ring at radius r (in-place).
        Uses type-specific grid for efficiency with sparse types.
        """
        type_grid = self.spatial_grid._type_grids[agent_type]
        
        # Top row (full width): y = cy - r
        if cy - r >= 0:
            for x in range(max(0, cx - r), min(self.width, cx + r + 1)):
                neighbor_ids.update(type_grid[(x, cy - r)])
        
        # Bottom row (full width): y = cy + r
        if cy + r < self.height:
            for x in range(max(0, cx - r), min(self.width, cx + r + 1)):
                neighbor_ids.update(type_grid[(x, cy + r)])
        
        # Left column (excluding corners): x = cx - r
        if cx - r >= 0:
            for y in range(max(0, cy - r + 1), min(self.height, cy + r)):
                neighbor_ids.update(type_grid[(cx - r, y)])
        
        # Right column (excluding corners): x = cx + r
        if cx + r < self.width:
            for y in range(max(0, cy - r + 1), min(self.height, cy + r)):
                neighbor_ids.update(type_grid[(cx + r, y)])
    
    def _get_required_radius(self, agent: Agent) -> int:
        """Get the cell search radius needed for this agent's vision."""
        if agent.agent_type == AgentType.PREY:
            max_vision = self.PREY_VISION_RANGE
        else:
            max_vision = self.PREDATOR_VISION_RANGE
        
        effective_vision = min(agent.vision, max_vision)
        if effective_vision <= 0:
            return 0
        return int(effective_vision) + 1

    def _find_nearest_agent(self, searcher: Agent, target_type: AgentType) -> Optional[Agent]:
        """
        Find nearest agent of target_type within searcher's effective vision.
        
        Uses Euclidean distance with float positions.
        Vision cap determined by searcher's type.
        """
        # Get vision cap based on searcher type
        if searcher.agent_type == AgentType.PREY:
            max_vision = self.PREY_VISION_RANGE
        else:
            max_vision = self.PREDATOR_VISION_RANGE
        
        effective_vision = min(searcher.vision, max_vision)
        if effective_vision <= 0:  # Blind agents can't see
            return None
        
        r = int(effective_vision) + 1  # Cell search radius (ceiling)
        vision_sq = effective_vision * effective_vision
        
        best_target = None
        best_dist_sq = float('inf')
        
        # Search cells in range - use type-specific query for efficiency
        cx, cy = int(searcher.x), int(searcher.y)
        for agent_id in self.spatial_grid.agents_of_type_in_region(
            max(0, cx - r), max(0, cy - r),
            min(self.width - 1, cx + r), min(self.height - 1, cy + r),
            target_type.value
        ):
            if agent_id not in self.agents:
                continue
            agent = self.agents[agent_id]
            
            # Euclidean distance squared (float positions)
            dx = agent.x - searcher.x
            dy = agent.y - searcher.y
            dist_sq = dx * dx + dy * dy
            
            if dist_sq <= vision_sq and dist_sq < best_dist_sq:
                best_dist_sq = dist_sq
                best_target = agent
        
        return best_target
    
    # =========================================================================
    # MAIN SIMULATION STEP
    # =========================================================================
    
    def step(self):
        """
        Execute one simulation tick using cell-based processing.
        
        PHASE 1: Cell-based local interactions (no movement)
          - Predators eat prey (by ID order, eat closest catchable prey)
          - Surviving prey share food fairly
          
        PHASE 2: Movement with type-specific neighborhood queries
          - Prey search for predators (sparse - fast query)
          - Predators search for prey (dense - grouped by radius)
          
        PHASE 3: Death checks and reproduction (integrated into Phase 2)
        
        Any agent-modifying commands (spawn/kill) received during step() are
        deferred and executed after the step completes.
        """
        self._stepping = True
        try:
            self._step_inner()
        finally:
            self._stepping = False
            self._process_deferred_commands()
    
    def _step_inner(self):
        """Internal step logic - called by step() with _stepping flag managed."""

        #bookkeeping, resetting counters etc.
        self.tick += 1
        # Reset per-tick statistics
        for key in self._stats:
            self._stats[key] = 0

        # initialize histograms
        hist_counts = {k: [0] * self.HIST_BINS for k in self._hist_max}
        hist_new_max = {k: 0.0 for k in self._hist_max}
        
        # Food regeneration with biome modifiers
        for pos in self.food:
            biome = self.biome[pos]
            food_mult, _, season_mult = self.BIOME_MODIFIERS[biome]
            
            if self.SEASON_STRENGTH == 0 or season_mult == 0:
                current_ceiling = self.FOOD_CEILING * food_mult
                current_regen = self.FOOD_REGEN_PER_TURN * food_mult
            else:
                effective_season = self.SEASON_STRENGTH * season_mult
                season_effect = 1 - 0.5 * effective_season * (1 + math.sin((self.tick % self.SEASON_LENGTH) / self.SEASON_LENGTH * TWOPI))
                current_ceiling = self.FOOD_CEILING * food_mult * season_effect
                current_regen = self.FOOD_REGEN_PER_TURN * food_mult * season_effect
            
            self.food[pos] = min(current_ceiling, self.food[pos] + current_regen)
        
        # PHASE 1: Cell-based local interactions
        cells_with_agents = set()
        for agent in self.agents.values():
            cells_with_agents.add(agent.cell)
        
        for cell in cells_with_agents:
            self._process_cell_interactions(cell)

        self._fled_this_tick.clear()   # Reset fleeing tracker for new tick
        self._mated_this_tick.clear()  # Reset mating tracker for new tick
                
        # PHASE 2: Movement with shared neighborhood per cell
        # Optimizations:
        # - Pre-compute predators in neighborhood ONCE (skip threat check if none)
        # - Pre-compute fertile mates in cell (check cell first before wider scan)
        
        cell_agents = defaultdict(list)
        for agent in self.agents.values():
            cell_agents[agent.cell].append(agent)
        
        for (cx, cy), agents_here in cell_agents.items():
            # Separate valid agents by type
            valid_prey = []
            valid_predators = []
            for agent in agents_here:
                if agent.id not in self.agents:
                    continue
                if agent.agent_type == AgentType.PREY:
                    valid_prey.append(agent)
                else:
                    valid_predators.append(agent)

                ###histogram updates, should be delegated to subroutine
                if agent.agent_type == AgentType.PREY:
                    prefix = 'prey'
                    max_e = self._hist_max['prey_energy']
                    shifted_speed = agent.speed + (self.PREY_SPRINT_RATIO - 1)    

                else:
                    prefix = 'predator'
                    max_e = self._hist_max['predator_energy']
                    shifted_speed = agent.speed + (self.PREDATOR_SPRINT_RATIO - 1)    

                age = self.tick - agent.born

                # Update histograms (clamp to last bin if exceeds max, also clamp to first to accommodate negative values)
                for suffix, get_val in self.hist_value_getters.items():
                    key = f'{prefix}_{suffix}'
                    if key in self._hist_max:
                        val = get_val(agent)
                        max_val = 1.0 if self._hist_max[key] == -1 else self._hist_max[key]
                        if max_val > 0:
                            bin_idx = int(val * self.HIST_BINS / max_val)
                            bin_idx = max(0, min(bin_idx, self.HIST_BINS - 1))
                        else:
                            bin_idx = 0
                        hist_counts[key][bin_idx] += 1
                        hist_new_max[key] = max(hist_new_max[key], val)
                ###end of histogram update sequence

            if not valid_prey and not valid_predators:
                continue
            
            # Find max vision radius needed in this cell
            max_r = 0
            for agent in valid_prey + valid_predators:
                if agent.vision > 0:
                    if agent.agent_type == AgentType.PREY:
                        max_vision = self.PREY_VISION_RANGE
                    else:
                        max_vision = self.PREDATOR_VISION_RANGE
                    r = int(min(agent.vision, max_vision)) + 1
                    if r > max_r:
                        max_r = r
            
            # Build neighborhood to max radius
            neighbor_ids = set(self.spatial_grid.agents_at(cx, cy))
            for r in range(1, max_r + 1):
                self._expand_ring(neighbor_ids, cx, cy, r)
            
            # Pre-compute predators in neighborhood (for prey threat detection)
            # This replaces N prey each scanning full neighbor_ids for predators
            predators_nearby = []
            for aid in neighbor_ids:
                if aid in self.agents and self.agents[aid].agent_type == AgentType.PREDATOR:
                    predators_nearby.append(self.agents[aid])
            
            # Pre-compute fertile prey in THIS cell (for quick mate check)
            fertile_prey_males_in_cell = []
            fertile_prey_females_in_cell = []
            for prey in valid_prey:
                investment = self.PREY_MALE_INVESTMENT if prey.sex == 0 else (1 - self.PREY_MALE_INVESTMENT)
                threshold = self.PREY_REPRODUCTION_THRESHOLD * investment + (prey.brutpflege if prey.sex else 0)
                if prey.energy >= threshold:
                    if prey.sex == 0:
                        fertile_prey_males_in_cell.append(prey)
                    else:
                        fertile_prey_females_in_cell.append(prey)
            
            # Pre-compute fertile predators in THIS cell
            fertile_pred_males_in_cell = []
            fertile_pred_females_in_cell = []
            for pred in valid_predators:
                investment = self.PREDATOR_MALE_INVESTMENT if pred.sex == 0 else (1 - self.PREDATOR_MALE_INVESTMENT)
                threshold = self.PREDATOR_REPRODUCTION_THRESHOLD * investment + (pred.brutpflege if pred.sex else 0)
                if pred.energy >= threshold:
                    if pred.sex == 0:
                        fertile_pred_males_in_cell.append(pred)
                    else:
                        fertile_pred_females_in_cell.append(pred)
            
            # Build cell context for optimized processing
            cell_context = {
                'neighbor_ids': neighbor_ids,
                'predators_nearby': predators_nearby,
                'fertile_prey_males': fertile_prey_males_in_cell,
                'fertile_prey_females': fertile_prey_females_in_cell,
                'fertile_pred_males': fertile_pred_males_in_cell,
                'fertile_pred_females': fertile_pred_females_in_cell,
            }
            
            # Process all agents with the shared context
            for prey in valid_prey:
                self._process_prey_turn_optimized(prey, cell_context)
            for predator in valid_predators:
                self._process_predator_turn_optimized(predator, cell_context)
            #end of cell loop

        #update histogram maxima for next turn: build into temp, then atomic assign:
        self._hist_max.update({k: max(v, 0.001) for k, v in hist_new_max.items() if self._hist_max[k] != -1})
        new_histograms = {k: [self._hist_max[k]] + counts for k, counts in hist_counts.items()}
        self._histograms = new_histograms  # atomic rebind
        self._hist_tick = self.tick  # record when histograms were built

        # PHASE 3: Prey combat (after movement, so we know who fled)
        # Combat is cell-based - prey in same cell can attack each other
        self._process_prey_combat()
        
        # Log statistics periodically
        if self.STATS_LOG_INTERVAL > 0 and self.tick % self.STATS_LOG_INTERVAL == 0:
            self._log_stats()
        
        self._check_termination()
    
    def _process_deferred_commands(self):
        """
        Process commands that were deferred during step().
        
        Called automatically at the end of step() after _stepping is set to False.
        Commands are processed in FIFO order.
        
        Note: Results from deferred commands are not returned to the original caller.
        The caller received {'status': 'deferred', ...} at call time.
        """
        while self._deferred_commands:
            name, value = self._deferred_commands.pop(0)
            # Re-call set_param - now _stepping is False so it will execute
            result = self.set_param(name, value)
            # Log if command failed (for debugging)
            if result.get('status') == 'error':
                print(f"Deferred command '{name}' failed: {result.get('message')}")
    
    def _process_prey_combat(self):
        """
        Process intra- and inter-species aggression among prey.
        
        Each prey can attack once per tick.
        Prey who fled this tick don't participate in combat.
        
        Priority:
        1. Attack other-species prey if present (prob = aggression)
        2. Else attack same-species same-sex prey (prob = 0.5 * aggression * (1 - dist²/threshold²))
        
        Never attack potential mates (opposite sex, same species).
        """
        # Build cell -> all non-fleeing prey
        cell_prey = defaultdict(list)
        for agent in self.agents.values():
            if agent.agent_type != AgentType.PREY:
                continue
            if agent.id in self._fled_this_tick:
                continue
            cell_prey[agent.cell].append(agent)
        
        # Track who has attacked this tick (one attack per agent)
        attacked_this_tick = set()
        
        for cell, prey_list in cell_prey.items():
            if len(prey_list) < 2:
                continue  # Need at least 2 prey for combat
            
            for attacker in prey_list:
                # Filter for valid attacker
                if attacker.aggression <= 0:
                    continue
                if attacker.energy < 2 * self.PREY_ATTACK_COST:
                    continue
                if attacker.id in attacked_this_tick:
                    continue
                if attacker.id not in self.agents:
                    continue  # May have died from another attack
                
                # Find potential targets
                other_species = []
                same_species_same_sex = []
                
                for target in prey_list:
                    if target.id == attacker.id:
                        continue
                    if target.id not in self.agents:
                        continue  # May have died
                    
                    genetic_dist_sq = sum((attacker.genes[i] - target.genes[i])**2 for i in range(3))
                    
                    if genetic_dist_sq > self.PREY_SPECIATION_THRESHOLD_SQ:
                        # Different species
                        other_species.append(target)
                    elif target.sex == attacker.sex:
                        # Same species, same sex (valid target)
                        same_species_same_sex.append((target, genetic_dist_sq))
                    # else: same species, opposite sex = potential mate, don't attack
                
                # Priority 1: Attack other species
                if other_species:
                    target = other_species[0]  # Just pick first one
                    p_attack = attacker.aggression
                    if attacker.sex == 1:  # Female
                        p_attack *= self.PREY_FEMALE_AGGRESSION_MULT
                    
                    if random.random() < p_attack:
                        attacker.energy -= self.PREY_ATTACK_COST
                        target.energy -= self.PREY_ATTACK_DAMAGE
                        attacked_this_tick.add(attacker.id)
                        self._stats['attacks_interspecies'] += 1
                        
                        # Repel attacker away from target (creates species separation)
                        dx = attacker.x - target.x
                        dy = attacker.y - target.y
                        dist = math.sqrt(dx * dx + dy * dy)
                        if dist > 0.001:
                            repel_x = attacker.x + self.PREY_ATTACK_REPULSION * dx / dist
                            repel_y = attacker.y + self.PREY_ATTACK_REPULSION * dy / dist
                            self._move_agent(attacker, repel_x, repel_y)
                        
                        if target.energy <= 0:
                            self.remove_agent(target.id)
                            self._stats['attack_kills'] += 1
                    continue
                
                # Priority 2: Attack same species, same sex
                if same_species_same_sex:
                    target, genetic_dist_sq = same_species_same_sex[0]  # Just pick first one
                    # Probability decreases with genetic closeness
                    p_attack = self.PREY_INTRASPECIES_MODIFIER * attacker.aggression * (1 - genetic_dist_sq / self.PREY_SPECIATION_THRESHOLD_SQ)
                    if attacker.sex == 1:  # Female
                        p_attack *= self.PREY_FEMALE_AGGRESSION_MULT
                    
                    if p_attack > 0 and random.random() < p_attack:
                        attacker.energy -= self.PREY_ATTACK_COST*self.PREY_INTRASPECIES_MODIFIER
                        target.energy -= self.PREY_ATTACK_DAMAGE*self.PREY_INTRASPECIES_MODIFIER
                        attacked_this_tick.add(attacker.id)
                        self._stats['attacks_intraspecies'] += 1
                        
                        if target.energy <= 0:
                            self.remove_agent(target.id)
                            self._stats['attack_kills'] += 1
    
    def _process_cell_interactions(self, cell: Tuple[int, int]):
        """
        Process local interactions within a cell: predation and food distribution.
        
        Order of operations:
        1. Predators pay metabolism and hunt (by ID order, eat closest catchable prey first)
        2. Surviving prey pay metabolism and share food fairly
        
        With float positions, predator catches prey if Euclidean distance < CATCH_DISTANCE.
        Metabolism costs are modified by biome cost_mult.
        """
        x, y = cell
        agent_ids = list(self.spatial_grid.agents_at(x, y))
        
        # Get biome cost multiplier for this cell
        _, cost_mult, _ = self.BIOME_MODIFIERS[self.biome[cell]]
        
        prey_in_cell = []
        predators_in_cell = []
        
        for aid in agent_ids:
            if aid not in self.agents:
                continue
            agent = self.agents[aid]
            if agent.agent_type == AgentType.PREY:
                prey_in_cell.append(agent)
            else:
                predators_in_cell.append(agent)
        
        # 1. Predators pay metabolism and hunt
        # Process by ID order (rare to have multiple predators in one cell)
        
        predators_in_cell.sort(key=lambda a: a.id)
        for predator in predators_in_cell:
            predator.energy -= (self.PREDATOR_METABOLISM_COST * cost_mult + self.PREDATOR_VISION_COST * predator.vision)
            if (self.tick - predator.born) > self.PREDATOR_SENESCENCE:
                predator.energy -= (self.tick - predator.born - self.PREDATOR_SENESCENCE) * self.PREDATOR_AGE_PENALTY
            
            if prey_in_cell:
                # Find catchable prey (within CATCH_DIST_SQ), eat closest one
                catchable = []
                for prey in prey_in_cell:
                    dx = prey.x - predator.x
                    dy = prey.y - predator.y
                    dist_sq = dx * dx + dy * dy
                    if dist_sq < self.CATCH_DIST_SQ:
                        catchable.append((dist_sq, prey))
                
                if catchable:
                    # Eat closest catchable prey
                    catchable.sort(key=lambda x: x[0])
                    target = catchable[0][1]
                    predator.energy += self.PREDATOR_ENERGY_FROM_PREY
                    if target.energy > self.PREDATOR_ENERGY_FROM_PREY:
                        predator.energy += (target.energy - self.PREDATOR_ENERGY_FROM_PREY) * 0.2 # bonus for fat prey
                    prey_in_cell.remove(target)
                    self.remove_agent(target.id)
        
        # 2. Surviving prey pay metabolism and share food fairly
        if prey_in_cell:
            # All prey pay metabolism
            for prey in prey_in_cell:
                prey.energy -= (self.PREY_METABOLISM_COST * cost_mult + self.PREY_VISION_COST * prey.vision)
                if (self.tick - prey.born) > self.PREY_SENESCENCE:
                    prey.energy -= (self.tick - prey.born - self.PREY_SENESCENCE) * self.PREY_AGE_PENALTY
            
            # Only hungry prey eat (energy < PREY_MAX_ENERGY)
            hungry_prey = [p for p in prey_in_cell if p.energy < self.PREY_MAX_ENERGY]
            
            if hungry_prey:
                food_here = self.food[(x, y)]
                n_hungry = len(hungry_prey)
                total_demand = n_hungry * self.PREY_EAT_RATE
                
                # Fair distribution: if not enough food, split equally
                if food_here >= total_demand:
                    per_prey = self.PREY_EAT_RATE
                    consumed = total_demand
                else:
                    per_prey = food_here / n_hungry
                    consumed = food_here
                
                self.food[(x, y)] -= consumed
                
                for prey in hungry_prey:
                    if prey.id in self._fled_this_tick: # value from last tick, do not eat during/right after fleeing
                        continue
                    prey.energy += per_prey
        
        # Note: Mating is handled during agent turn processing, not here
    
    # =========================================================================
    # AGENT TURN PROCESSING (called from step with pre-built neighborhood)
    # =========================================================================
    
    def _process_prey_turn_optimized(self, prey: Agent, ctx: dict):
        """
        Optimized prey turn processing using pre-computed cell context.
        
        Optimizations:
        - Skip threat scan if no predators in neighborhood
        - Check cell first for mates before wider scan
        """
        if prey.energy <= 0:
            self.remove_agent(prey.id)
            return
        
        _, cost_mult, _ = self.get_biome_modifier(prey.x, prey.y)
        
        # Calculate sprint parameters
        effective_sprint_ratio = self.PREY_SPRINT_RATIO + prey.speed
        sprint_speed = self.PREY_MOVE_SPEED * effective_sprint_ratio
        ratio_factor = effective_sprint_ratio / self.PREY_SPRINT_RATIO
        base_sprint_cost = max(self.PREY_MOVE_COST, 
                              self.PREY_MOVE_COST * self.PREY_SPRINT_COST * ratio_factor * ratio_factor)
        sprint_cost = base_sprint_cost * cost_mult
        
        # 1. Check for threats - ONLY if predators exist in neighborhood
        if prey.energy >= sprint_cost and prey.vision > 0 and ctx['predators_nearby']:
            # Scan pre-filtered predator list instead of full neighbor_ids
            threat = self._find_closest_agent_in_list(prey, ctx['predators_nearby'])
            if threat:
                prey.bearing = math.atan2(prey.y - threat.y, prey.x - threat.x)
                self._move_away_from(prey, threat.x, threat.y, speed=sprint_speed)
                prey.energy -= sprint_cost
                self._fled_this_tick.add(prey.id)  # Mark as fled (skip combat)
                if prey.energy <= 0:
                    self.remove_agent(prey.id)
                return
        
        # 2. No threat - decide between forage and seek mate
        is_satiated = prey.energy >= self.PREY_MAX_ENERGY
        if prey.sex:
            fertility_threshold = self.PREY_REPRODUCTION_THRESHOLD*(1 - self.PREY_MALE_INVESTMENT) +  prey.brutpflege
        else:
            fertility_threshold = self.PREY_REPRODUCTION_THRESHOLD * self.PREY_MALE_INVESTMENT
        is_fertile = prey.energy >= fertility_threshold
        
        if not is_satiated:
            # Forage
            self._forage(prey, cost_mult)
        elif is_fertile and prey.vision > 0:
            # Satiated and fertile - seek mate
            # Check cell first (optimization)
            opposite_sex_in_cell = ctx['fertile_prey_males'] if prey.sex == 1 else ctx['fertile_prey_females']
            mate = self._find_closest_agent_in_list(prey, opposite_sex_in_cell, exclude_self=True)
            
            if not mate:
                # No mate in cell, scan wider neighborhood
                mate = self._find_nearest_mate(prey, ctx['neighbor_ids'])
            
            if mate:
                prey.bearing = math.atan2(mate.y - prey.y, mate.x - prey.x)
                self._move_toward(prey, mate.x, mate.y, speed=self.PREY_MOVE_SPEED)
                prey.energy -= self.PREY_MOVE_COST * cost_mult
            else:
                # Satiated, fertile, but no mate visible - drift bearing and nudge forward (no energy cost)
                prey.bearing = _normalize_bearing(prey.bearing + random.gauss(0, self.PREY_BEARING_DRIFT_STDDEV))
                nudge_speed = 0.5 * self.PREY_MOVE_SPEED
                new_x = prey.x + math.cos(prey.bearing) * nudge_speed
                new_y = prey.y + math.sin(prey.bearing) * nudge_speed
                self._move_agent(prey, new_x, new_y)
        else:
            # Satiated but not fertile - drift bearing and nudge forward (no energy cost)
            prey.bearing = _normalize_bearing(prey.bearing + random.gauss(0, self.PREY_BEARING_DRIFT_STDDEV))
            nudge_speed = 0.5 * self.PREY_MOVE_SPEED
            new_x = prey.x + math.cos(prey.bearing) * nudge_speed
            new_y = prey.y + math.sin(prey.bearing) * nudge_speed
            self._move_agent(prey, new_x, new_y)
        
        # Check for death after movement
        if prey.energy <= 0:
            self.remove_agent(prey.id)
            return
        
        # 3. After movement: try to mate if fertile
        is_fertile = prey.energy >= fertility_threshold
        if is_fertile:
            # Check cell first for mating
            opposite_sex_in_cell = ctx['fertile_prey_males'] if prey.sex == 1 else ctx['fertile_prey_females']
            mate = self._find_closest_agent_in_list(
                prey, opposite_sex_in_cell, 
                max_dist_sq=self.MATING_DISTANCE * self.MATING_DISTANCE,
                exclude_self=True
            )
            if not mate:
                # Check wider neighborhood
                mate = self._find_fertile_mate_within_distance(prey, ctx['neighbor_ids'])
            if mate:
                self._mate(prey, mate)
    
    def _process_predator_turn_optimized(self, predator: Agent, ctx: dict):
        """
        Optimized predator turn processing using pre-computed cell context.
        
        Predators benefit less from optimization (fewer of them), but we still:
        - Check cell first for mates before wider scan
        
        Behavior priority when fertile but not satiated is probabilistic:
        - At fertility threshold: almost always hunt
        - At satiation: always seek mates
        - Linear interpolation in between
        """
        if predator.energy <= 0:
            self.remove_agent(predator.id)
            return
        
        _, cost_mult, _ = self.get_biome_modifier(predator.x, predator.y)
        
        # Calculate sprint parameters
        effective_sprint_ratio = self.PREDATOR_SPRINT_RATIO + predator.speed
        sprint_speed = self.PREDATOR_MOVE_SPEED * effective_sprint_ratio
        ratio_factor = effective_sprint_ratio / self.PREDATOR_SPRINT_RATIO
        base_sprint_cost = max(self.PREDATOR_MOVE_COST,
                              self.PREDATOR_MOVE_COST * self.PREDATOR_SPRINT_COST * ratio_factor * ratio_factor)
        sprint_cost = base_sprint_cost * cost_mult
        
        # Check states
        is_satiated = predator.energy >= self.PREDATOR_MAX_ENERGY
        if predator.sex:
            fertility_threshold = self.PREDATOR_REPRODUCTION_THRESHOLD*(1 - self.PREDATOR_MALE_INVESTMENT) +  predator.brutpflege
        else:
            fertility_threshold = self.PREDATOR_REPRODUCTION_THRESHOLD * self.PREDATOR_MALE_INVESTMENT
        is_fertile = predator.energy >= fertility_threshold
        
        # Decide whether to prioritize hunting or mating
        # - Not fertile: always hunt (can't mate anyway)
        # - Satiated: always seek mate
        # - Fertile but not satiated: probabilistic based on energy level and heat gene
        #   P(hunt) = 1 - energy_ratio^heat
        #   Lower heat = energy_ratio^heat is larger = lower P(hunt) = more mating focus
        prioritize_hunting = True
        if is_satiated:
            prioritize_hunting = False
        elif is_fertile:
            energy_range = self.PREDATOR_MAX_ENERGY - fertility_threshold
            if energy_range > 0:
                energy_ratio = (predator.energy - fertility_threshold) / energy_range
                # Power transformation: heat < 1 favors mating, heat > 1 favors hunting
                prioritize_hunting = random.random() > energy_ratio ** (predator.heat/100)
            else:
                prioritize_hunting = False
        
        # 1. Hunt if prioritizing hunting and prey visible
        if prioritize_hunting:
            target = self._find_nearest_in_set(predator, AgentType.PREY, ctx['neighbor_ids'])
            if target:
                predator.bearing = math.atan2(target.y - predator.y, target.x - predator.x)
                self._move_toward(predator, target.x, target.y, speed=sprint_speed)
                predator.energy -= sprint_cost
                if predator.energy <= 0:
                    self.remove_agent(predator.id)
                return
        
        # 2. Not hunting (or no prey found) - check for mating or patrol
        if is_fertile:
            # 2a. Try to mate - check cell first
            opposite_sex_in_cell = ctx['fertile_pred_males'] if predator.sex == 1 else ctx['fertile_pred_females']
            mate = self._find_closest_agent_in_list(
                predator, opposite_sex_in_cell,
                max_dist_sq=self.MATING_DISTANCE * self.MATING_DISTANCE,
                exclude_self=True
            )
            if not mate:
                mate = self._find_fertile_mate_within_distance(predator, ctx['neighbor_ids'])
            
            if mate:
                self._mate(predator, mate)
                return
            
            # 2b. Seek visible fertile mate - check cell first
            if predator.vision > 0:
                distant_mate = self._find_closest_agent_in_list(predator, opposite_sex_in_cell, exclude_self=True)
                if not distant_mate:
                    distant_mate = self._find_nearest_mate(predator, ctx['neighbor_ids'])
                
                if distant_mate:
                    predator.bearing = math.atan2(distant_mate.y - predator.y, distant_mate.x - predator.x)
                    self._move_toward(predator, distant_mate.x, distant_mate.y, speed=self.PREDATOR_MOVE_SPEED)
                    predator.energy -= self.PREDATOR_MOVE_COST * cost_mult
                    if predator.energy <= 0:
                        self.remove_agent(predator.id)
                    return
        
        # 3. Patrol (not fertile, or no mate found)
        # If we chose not to hunt but didn't find a mate, try hunting as fallback
        if not prioritize_hunting and not is_satiated:
            target = self._find_nearest_in_set(predator, AgentType.PREY, ctx['neighbor_ids'])
            if target:
                predator.bearing = math.atan2(target.y - predator.y, target.x - predator.x)
                self._move_toward(predator, target.x, target.y, speed=sprint_speed)
                predator.energy -= sprint_cost
                if predator.energy <= 0:
                    self.remove_agent(predator.id)
                return
        
        # 4. Patrol: drift bearing slightly, then move
        predator.bearing += random.gauss(0, 0.1)
        predator.bearing = predator.bearing % TWOPI
        
        dx = math.cos(predator.bearing) * self.PREDATOR_MOVE_SPEED
        dy = math.sin(predator.bearing) * self.PREDATOR_MOVE_SPEED
        next_x = predator.x + dx
        next_y = predator.y + dy
        
        if not self.is_passable(next_x, next_y):
            predator.bearing = self._random_inward_bearing(predator.x, predator.y)
            dx = math.cos(predator.bearing) * self.PREDATOR_MOVE_SPEED
            dy = math.sin(predator.bearing) * self.PREDATOR_MOVE_SPEED
            next_x = predator.x + dx
            next_y = predator.y + dy
            
            if not self.is_passable(next_x, next_y):
                next_x, next_y = predator.x, predator.y
        
        self._move_agent(predator, next_x, next_y)
        predator.energy -= self.PREDATOR_MOVE_COST * cost_mult
        
        if predator.energy <= 0:
            self.remove_agent(predator.id)
    
    def _find_closest_agent_in_list(self, searcher: Agent, candidates: List[Agent], 
                                     max_dist_sq: float = None, exclude_self: bool = False) -> Optional[Agent]:
        """
        Find closest agent from a pre-filtered list.
        
        Args:
            searcher: The agent searching
            candidates: Pre-filtered list of potential targets
            max_dist_sq: Maximum distance squared (None = use searcher's vision)
            exclude_self: If True, skip searcher in candidates list
        
        Returns:
            Closest agent within range, or None
        """
        if not candidates:
            return None
        
        if max_dist_sq is None:
            # Use searcher's vision
            if searcher.agent_type == AgentType.PREY:
                max_vision = self.PREY_VISION_RANGE
            else:
                max_vision = self.PREDATOR_VISION_RANGE
            effective_vision = min(searcher.vision, max_vision)
            if effective_vision <= 0:
                return None
            max_dist_sq = effective_vision * effective_vision
        
        best_agent = None
        best_dist_sq = max_dist_sq
        
        for agent in candidates:
            if exclude_self and agent.id == searcher.id:
                continue
            if agent.id not in self.agents:
                continue
            
            dx = agent.x - searcher.x
            dy = agent.y - searcher.y
            dist_sq = dx * dx + dy * dy
            
            if dist_sq < best_dist_sq:
                best_dist_sq = dist_sq
                best_agent = agent
        
        return best_agent
    
    def _find_fertile_mate_within_distance(self, agent: Agent, neighbor_ids: Set[int]) -> Optional[Agent]:
        """Find nearest fertile mate (opposite sex, same type, same species) within MATING_DISTANCE."""
        mating_dist_sq = self.MATING_DISTANCE * self.MATING_DISTANCE
        
        if agent.agent_type == AgentType.PREY:
            threshold = self.PREY_REPRODUCTION_THRESHOLD
            male_inv = self.PREY_MALE_INVESTMENT
            speciation_threshold_sq = self.PREY_SPECIATION_THRESHOLD_SQ
        else:
            threshold = self.PREDATOR_REPRODUCTION_THRESHOLD
            male_inv = self.PREDATOR_MALE_INVESTMENT
            speciation_threshold_sq = None  # Predators don't use speciation yet
        
        opposite_sex = 1 - agent.sex
        
        best_mate = None
        best_dist_sq = mating_dist_sq
        
        for agent_id in neighbor_ids:
            if agent_id not in self.agents or agent_id == agent.id:
                continue
            candidate = self.agents[agent_id]
            
            # Must be same type, opposite sex
            if candidate.agent_type != agent.agent_type or candidate.sex != opposite_sex:
                continue
            
            # Check speciation (prey only) - must be same species
            if speciation_threshold_sq is not None:
                genetic_dist_sq = sum((agent.genes[i] - candidate.genes[i])**2 for i in range(3))
                if genetic_dist_sq > speciation_threshold_sq:
                    continue  # Different species, can't mate
            
            # Check if within mating distance
            dx = candidate.x - agent.x
            dy = candidate.y - agent.y
            dist_sq = dx * dx + dy * dy
            if dist_sq >= mating_dist_sq:
                continue
            
            # Check if candidate is fertile
            candidate_investment = male_inv if candidate.sex == 0 else (1 - male_inv)
            candidate_threshold = threshold * candidate_investment + (candidate.brutpflege if candidate.sex else 0)
            if candidate.energy < candidate_threshold:
                continue
            
            if dist_sq < best_dist_sq:
                best_dist_sq = dist_sq
                best_mate = candidate
        
        return best_mate
    
    def _mate(self, initiator: Agent, partner: Agent):
        """
        Execute mating between two agents.
        Both must be fertile (caller should verify).
        Each agent can only mate once per tick.
        """
        # Check if either has already mated this tick
        if initiator.id in self._mated_this_tick or partner.id in self._mated_this_tick:
            return
        
        if initiator.agent_type == AgentType.PREY:
            threshold = self.PREY_REPRODUCTION_THRESHOLD
            male_inv = self.PREY_MALE_INVESTMENT
            initial_e = self.PREY_INITIAL_ENERGY
            repro_cost = self.PREY_REPRODUCTION_COST
        else:
            threshold = self.PREDATOR_REPRODUCTION_THRESHOLD
            male_inv = self.PREDATOR_MALE_INVESTMENT
            initial_e = self.PREDATOR_INITIAL_ENERGY
            repro_cost = self.PREDATOR_REPRODUCTION_COST
        
        # Determine male and female
        if initiator.sex == 0:
            male, female = initiator, partner
        else:
            male, female = partner, initiator
        
        # Split cost by investment ratio
        male.energy -= repro_cost * male_inv
        female.energy -= (repro_cost * (1 - male_inv) + female.brutpflege)
        
        # Mark both as having mated this tick
        self._mated_this_tick.add(initiator.id)
        self._mated_this_tick.add(partner.id)
        self._stats['matings'] += 1
        
        # Spawn near female
        spawn_x = female.x + random.uniform(-0.3, 0.3)
        spawn_y = female.y + random.uniform(-0.3, 0.3)
        # Clamp to world bounds
        spawn_x = max(0, min(self.width - 0.01, spawn_x))
        spawn_y = max(0, min(self.height - 0.01, spawn_y))
        
        offspring_energy = initial_e + female.brutpflege
        
        self.spawn_agent(initiator.agent_type, spawn_x, spawn_y, 
                        offspring_energy, parent_ids=(male.id, female.id))
    
    def _find_nearest_in_set(self, searcher: Agent, target_type: AgentType, candidate_ids: Set[int]) -> Optional[Agent]:
        """
        Find nearest agent of target_type within searcher's vision from candidate set.
        """
        if searcher.agent_type == AgentType.PREY:
            max_vision = self.PREY_VISION_RANGE
        else:
            max_vision = self.PREDATOR_VISION_RANGE
        
        effective_vision = min(searcher.vision, max_vision)
        if effective_vision <= 0:
            return None
        
        vision_sq = effective_vision * effective_vision
        
        best_target = None
        best_dist_sq = float('inf')
        
        for agent_id in candidate_ids:
            if agent_id not in self.agents or agent_id == searcher.id:
                continue
            agent = self.agents[agent_id]
            if agent.agent_type != target_type:
                continue
            
            dx = agent.x - searcher.x
            dy = agent.y - searcher.y
            dist_sq = dx * dx + dy * dy
            
            if dist_sq <= vision_sq and dist_sq < best_dist_sq:
                best_dist_sq = dist_sq
                best_target = agent
        
        return best_target
    
    def _find_nearest_mate(self, searcher: Agent, candidate_ids: Set[int]) -> Optional[Agent]:
        """
        Find nearest fertile mate (same type, opposite sex, same species) within searcher's vision.
        A potential mate is fertile if they can afford their share of reproduction cost.
        """
        if searcher.agent_type == AgentType.PREY:
            max_vision = self.PREY_VISION_RANGE
            threshold = self.PREY_REPRODUCTION_THRESHOLD
            male_inv = self.PREY_MALE_INVESTMENT
            speciation_threshold_sq = self.PREY_SPECIATION_THRESHOLD_SQ
        else:
            max_vision = self.PREDATOR_VISION_RANGE
            threshold = self.PREDATOR_REPRODUCTION_THRESHOLD
            male_inv = self.PREDATOR_MALE_INVESTMENT
            speciation_threshold_sq = None  # Predators don't use speciation yet
        
        effective_vision = min(searcher.vision, max_vision)
        if effective_vision <= 0:
            return None
        
        vision_sq = effective_vision * effective_vision
        opposite_sex = 1 - searcher.sex  # 0 -> 1, 1 -> 0
        
        best_mate = None
        best_dist_sq = float('inf')
        
        for agent_id in candidate_ids:
            if agent_id not in self.agents or agent_id == searcher.id:
                continue
            agent = self.agents[agent_id]
            
            # Must be same type, opposite sex
            if agent.agent_type != searcher.agent_type or agent.sex != opposite_sex:
                continue
            
            # Check speciation (prey only) - must be same species
            if speciation_threshold_sq is not None:
                genetic_dist_sq = sum((searcher.genes[i] - agent.genes[i])**2 for i in range(3))
                if genetic_dist_sq > speciation_threshold_sq:
                    continue  # Different species, can't mate
            
            # Check if potential mate is fertile
            mate_investment = male_inv if agent.sex == 0 else (1 - male_inv)
            mate_threshold = threshold * mate_investment + (agent.brutpflege if agent.sex else 0)
            if agent.energy < mate_threshold:
                continue
            
            dx = agent.x - searcher.x
            dy = agent.y - searcher.y
            dist_sq = dx * dx + dy * dy
            
            if dist_sq <= vision_sq and dist_sq < best_dist_sq:
                best_dist_sq = dist_sq
                best_mate = agent
        
        return best_mate
    
    def _random_inward_bearing(self, x: float, y: float) -> float:
        """
        Pick a random bearing that points toward the interior of the world.
        Used when a predator hits a wall or edge during patrol.
        """
        center_x = self.width / 2
        center_y = self.height / 2
        
        # Base bearing toward center
        to_center = math.atan2(center_y - y, center_x - x)
        
        # Add some randomness (±45°)
        return to_center + random.uniform(-math.pi/4, math.pi/4)
    
    # =========================================================================
    # FORAGING BEHAVIOR
    # =========================================================================
    
    def _forage(self, prey, cost_mult):
        """
        Three-tier foraging logic:
        1. If satisfied with current cell's food, drift bearing and stay
        2. If moving, check only cells in current quadrant
        3. If nothing good ahead, do full neighbor review
        """
        cell_x, cell_y = int(prey.x), int(prey.y)
        cell = (cell_x, cell_y)
        current_food = self.food[cell]

        # --- Tier 1: Satisfied here? Drift bearing, don't move. ---
        if current_food > 2 * self.PREY_EAT_RATE:
            prey.bearing = _normalize_bearing(prey.bearing + random.gauss(0, self.PREY_BEARING_DRIFT_STDDEV))
            return False  # Signal: did not move

        # --- Tier 2: Moving - check quadrant only ---
        # Map bearing [-π, π] to quadrant 0-3: Q0=East, Q1=North, Q2=West, Q3=South
        quadrant = int((prey.bearing * 2 / math.pi + 4.5) % 4)

        best_food = current_food
        best_offsets = []  # Collect all cells with best food (to randomize if equal)

        for dx, dy in self._QUADRANT_CELLS[quadrant]:
            nc = (cell_x + dx, cell_y + dy)
            if nc in self.food:
                nfood = self.food[nc]
                if nfood > best_food:
                    best_food = nfood
                    best_offsets = [(dx, dy)]
                elif nfood == best_food and best_offsets:
                    best_offsets.append((dx, dy))

        if best_offsets:
            best_offset = random.choice(best_offsets)
            target = _get_bearing(*best_offset)
            diff = _normalize_bearing(target - prey.bearing)
            prey.bearing = _normalize_bearing(prey.bearing + max(-self.BEARING_ADJUST, min(self.BEARING_ADJUST, diff)))
            self._move_in_direction(prey, prey.bearing, cost_mult)
            return True

        # --- Tier 3: Nothing good in quadrant - full neighbor review ---
        neighbors = self._get_neighbors(cell_x, cell_y)
        if not neighbors:
            prey.bearing = _normalize_bearing(prey.bearing + random.gauss(0, self.PREY_BEARING_DRIFT_STDDEV))
            return False

        best_food = current_food
        best_offsets = []  # Collect all cells with best food (to randomize if equal)

        for nx, ny in neighbors:
            nfood = self.food[(nx, ny)]
            if nfood > best_food:
                best_food = nfood
                best_offsets = [(nx - cell_x, ny - cell_y)]
            elif nfood == best_food and best_offsets:
                best_offsets.append((nx - cell_x, ny - cell_y))

        if best_offsets:
            best_offset = random.choice(best_offsets)
            target = _get_bearing(*best_offset)
            diff = _normalize_bearing(target - prey.bearing)
            prey.bearing = _normalize_bearing(prey.bearing + max(-self.BEARING_ADJUST, min(self.BEARING_ADJUST, diff)))
            self._move_in_direction(prey, prey.bearing, cost_mult)
            return True
        else:
            prey.bearing = _normalize_bearing(prey.bearing + random.gauss(0, self.PREY_BEARING_DRIFT_STDDEV))
            return False

    def _herd_consensus(self, prey, prey_in_cell):
        if prey_in_cell:
            sum_cos, sum_sin, total_weight = 0.0, 0.0, 0.0  # Weighted circular mean of same-species bearings
            for other in prey_in_cell:  # Already have this list from cell_agents
                if other.id == prey.id:
                    continue
                gdist_sq = sum((prey.genes[i] - other.genes[i])**2 for i in range(3))
                if gdist_sq < self.PREY_SPECIATION_THRESHOLD_SQ:
                    weight = 1.0 - gdist_sq / self.PREY_SPECIATION_THRESHOLD_SQ  # Closer kin → stronger pull
                    sum_cos += weight * math.cos(other.bearing)
                    sum_sin += weight * math.sin(other.bearing)
                    total_weight += weight

            if total_weight > 0:
                R = math.sqrt(sum_cos**2 + sum_sin**2) / total_weight
                consensus = math.atan2(sum_sin, sum_cos)
                diff = R*_normalize_bearing(consensus - prey.bearing)
                return diff
            else:
                return 0
        else:
            return 0


    def _move_in_direction(self, agent, bearing, cost_mult):
        """Move agent one tick in given direction, with all safety checks."""
        new_x = agent.x + math.cos(bearing) * self.PREY_MOVE_SPEED
        new_y = agent.y + math.sin(bearing) * self.PREY_MOVE_SPEED
        self._move_agent(agent, new_x, new_y)  # Handles grid, bounds, walls, marks dirty
        agent.energy -= self.PREY_MOVE_COST * cost_mult
    
    # =========================================================================
    # TERMINATION AND UTILITIES
    # =========================================================================
    
    def _log_stats(self):
        """Log simulation statistics to console."""
        prey_count, pred_count = self.count_agents()
        
        # Collect genetic diversity stats for prey
        prey_agents = [a for a in self.agents.values() if a.agent_type == AgentType.PREY]
        if prey_agents:
            aggression_vals = [a.aggression for a in prey_agents]
            avg_aggression = sum(aggression_vals) / len(aggression_vals)
        else:
            avg_aggression = 0
        
        print(f"[Tick {self.tick}] Prey:{prey_count} Pred:{pred_count} | "
              f"Matings:{self._stats['matings']} | "
              f"Attacks(inter/intra):{self._stats['attacks_interspecies']}/{self._stats['attacks_intraspecies']} "
              f"Kills:{self._stats['attack_kills']} | "
              f"AvgAggr:{avg_aggression:.3f}")
    
    def _check_termination(self):
        """Check if simulation should halt."""
        prey_count, predator_count = self.count_agents()
        
        if predator_count == 0:
            self.halted = False #no longer terminate on predator extinction
            self.halt_reason = "Predators extinct"
        elif prey_count == 0 and predator_count == 0:
            self.halted = True
            self.halt_reason = "All extinct"
        elif (prey_count + predator_count) > 100000:
            self.halted = True
            self.halt_reason = "Population exploded"
    
    def count_agents(self) -> Tuple[int, int]:
        """Return (prey_count, predator_count)."""
        prey = sum(1 for a in self.agents.values() if a.agent_type == AgentType.PREY)
        predators = sum(1 for a in self.agents.values() if a.agent_type == AgentType.PREDATOR)
        return prey, predators
   
    def report_statistics(self) -> dict:
        histograms = self._histograms
    
        def sanitize(v):
            max_val = 1.0 if v[0] == -1 else v[0]
            return [round(max_val, 2)] + v[1:]
    
        return {
            'prey': {k.replace('prey_', ''): sanitize(v) for k, v in histograms.items() if k.startswith('prey')},
            'predator': {k.replace('predator_', ''): sanitize(v) for k, v in histograms.items() if k.startswith('predator')},
        }

    @classmethod
    def report_params(cls) -> dict:
        """
        Report simulation parameters.
        
        This is the simulation's own decision about what parameters to expose.
        Server should call this rather than accessing constants directly.
        """
        return {
            'sim': {
                'food_per_turn': cls.FOOD_REGEN_PER_TURN,
                'food_ceiling': cls.FOOD_CEILING,
                'season_length': cls.SEASON_LENGTH,
                'season_strength': cls.SEASON_STRENGTH,
                'map_width': cls.DEFAULT_WIDTH,
                'map_height': cls.DEFAULT_HEIGHT,
                'num_prey': cls.DEFAULT_INITIAL_PREY,    
                'num_pred_pairs': cls.DEFAULT_INITIAL_PREDATORS,    
                'histo_bins': cls.HIST_BINS,    
            },
            'prey': {
                'initial_energy': cls.PREY_INITIAL_ENERGY,
                'move_cost': cls.PREY_MOVE_COST,
                'move_speed': cls.PREY_MOVE_SPEED,
                'sprint_ratio': cls.PREY_SPRINT_RATIO,
                'sprint_cost': cls.PREY_SPRINT_COST,
                'metabolism_cost': cls.PREY_METABOLISM_COST,
                'reproduction_threshold': cls.PREY_REPRODUCTION_THRESHOLD,
                'reproduction_cost': cls.PREY_REPRODUCTION_COST,
                'eat_rate': cls.PREY_EAT_RATE,
                'vision_range': cls.PREY_VISION_RANGE,
                'vision_cost': cls.PREY_VISION_COST,
                'default_vision': cls.PREY_DEFAULT_VISION,
                'senescence': cls.PREY_SENESCENCE,
                'age_penalty': cls.PREY_AGE_PENALTY,
                'max_energy': cls.PREY_MAX_ENERGY,            
                'foraging_consensus': cls.CONSENSUS_STRENGTH,
                'male_investment': cls.PREY_MALE_INVESTMENT,
                'speciation_threshold': cls.PREY_SPECIATION_THRESHOLD_SQ,
                'speciation_mutation_rate': cls.PREY_GENE_MUTATION_SD,
                'aggression_mutation_rate': cls.PREY_AGGRESSION_MUTATION_SD,
                'attack_cost': cls.PREY_ATTACK_COST,
                'attack_damage': cls.PREY_ATTACK_DAMAGE,
                'default_aggression': cls.PREY_DEFAULT_AGGRESSION,
                'intraspecies_aggression': cls.PREY_INTRASPECIES_MODIFIER,
                'female_aggression': cls.PREY_FEMALE_AGGRESSION_MULT,
            },
            'predator': {
                'initial_energy': cls.PREDATOR_INITIAL_ENERGY,
                'move_cost': cls.PREDATOR_MOVE_COST,
                'move_speed': cls.PREDATOR_MOVE_SPEED,
                'sprint_ratio': cls.PREDATOR_SPRINT_RATIO,
                'sprint_cost': cls.PREDATOR_SPRINT_COST,
                'metabolism_cost': cls.PREDATOR_METABOLISM_COST,
                'reproduction_threshold': cls.PREDATOR_REPRODUCTION_THRESHOLD,
                'reproduction_cost': cls.PREDATOR_REPRODUCTION_COST,
                'energy_from_prey': cls.PREDATOR_ENERGY_FROM_PREY,
                'vision_range': cls.PREDATOR_VISION_RANGE,
                'default_vision': cls.PREDATOR_DEFAULT_VISION,
                'senescence': cls.PREDATOR_SENESCENCE,
                'age_penalty': cls.PREDATOR_AGE_PENALTY,
                'max_energy': cls.PREDATOR_MAX_ENERGY,
                'male_investment': cls.PREDATOR_MALE_INVESTMENT,
                'mating_focus': cls.PREDATOR_HEAT,
                'mating_distance': cls.MATING_DISTANCE,
                'catch_distance': cls.CATCH_DIST_SQ,
            }
        }
    
    # Mapping from JSON param names to class attribute names
    _PARAM_MAP = {
        'sim.food_per_turn': 'FOOD_REGEN_PER_TURN',
        'sim.food_ceiling': 'FOOD_CEILING',
        'sim.season_length': 'SEASON_LENGTH',
        'sim.season_strength': 'SEASON_STRENGTH',
        'sim.map_width': 'DEFAULT_WIDTH',
        'sim.map_height': 'DEFAULT_HEIGHT',
        'sim.num_prey': 'DEFAULT_INITIAL_PREY',    
        'sim.num_pred_pairs': 'DEFAULT_INITIAL_PREDATORS',    
        'sim.histo_bins': 'HIST_BINS',    
        'prey.initial_energy': 'PREY_INITIAL_ENERGY',
        'prey.move_cost': 'PREY_MOVE_COST',
        'prey.move_speed': 'PREY_MOVE_SPEED',
        'prey.sprint_ratio': 'PREY_SPRINT_RATIO',
        'prey.sprint_cost': 'PREY_SPRINT_COST',
        'prey.metabolism_cost': 'PREY_METABOLISM_COST',
        'prey.reproduction_threshold': 'PREY_REPRODUCTION_THRESHOLD',
        'prey.reproduction_cost': 'PREY_REPRODUCTION_COST',
        'prey.eat_rate': 'PREY_EAT_RATE',
        'prey.vision_range': 'PREY_VISION_RANGE',
        'prey.vision_cost': 'PREY_VISION_COST',
        'prey.default_vision': 'PREY_DEFAULT_VISION',
        'prey.senescence': 'PREY_SENESCENCE',
        'prey.age_penalty': 'PREY_AGE_PENALTY',
        'prey.max_energy': 'PREY_MAX_ENERGY',
        'prey.foraging_consensus': 'CONSENSUS_STRENGTH',
        'prey.male_investment': 'PREY_MALE_INVESTMENT',
        'prey.speciation_threshold': 'PREY_SPECIATION_THRESHOLD_SQ',
        'prey.speciation_mutation_rate': 'PREY_GENE_MUTATION_SD',
        'prey.aggression_mutation_rate': 'PREY_AGGRESSION_MUTATION_SD',
        'prey.attack_cost': 'PREY_ATTACK_COST',
        'prey.attack_damage': 'PREY_ATTACK_DAMAGE',
        'prey.default_aggression': 'PREY_DEFAULT_AGGRESSION',
        'prey.intraspecies_aggression': 'PREY_INTRASPECIES_MODIFIER',
        'prey.female_aggression': 'PREY_FEMALE_AGGRESSION_MULT',
        'predator.initial_energy': 'PREDATOR_INITIAL_ENERGY',
        'predator.move_cost': 'PREDATOR_MOVE_COST',
        'predator.move_speed': 'PREDATOR_MOVE_SPEED',
        'predator.sprint_ratio': 'PREDATOR_SPRINT_RATIO',
        'predator.sprint_cost': 'PREDATOR_SPRINT_COST',
        'predator.metabolism_cost': 'PREDATOR_METABOLISM_COST',
        'predator.reproduction_threshold': 'PREDATOR_REPRODUCTION_THRESHOLD',
        'predator.reproduction_cost': 'PREDATOR_REPRODUCTION_COST',
        'predator.energy_from_prey': 'PREDATOR_ENERGY_FROM_PREY',
        'predator.vision_range': 'PREDATOR_VISION_RANGE',
        'predator.vision_cost': 'PREDATOR_VISION_COST',
        'predator.default_vision': 'PREDATOR_DEFAULT_VISION',
        'predator.senescence': 'PREDATOR_SENESCENCE',
        'predator.age_penalty': 'PREDATOR_AGE_PENALTY',
        'predator.max_energy': 'PREDATOR_MAX_ENERGY',
        'predator.male_investment': 'PREDATOR_MALE_INVESTMENT',
        'predator.mating_focus': 'PREDATOR_HEAT',
        'predator.mating_distance': 'MATING_DISTANCE',
        'predator.catch_distance': 'CATCH_DIST_SQ',
    }
    
    # Commands that operate on the world instance
    _COMMANDS = {'kill_agent', 'clone_agent', 'boost_agent', 'breeding_pair'}
    
    @classmethod
    def _set_class_param(cls, name: str, value) -> dict:
        """
        Internal classmethod to set class-level simulation parameters.
        Called by set_param() when name is in _PARAM_MAP.
        """
        attr_name = cls._PARAM_MAP[name]
        
        # Get current value to determine type
        try:
            current_value = getattr(cls, attr_name)
        except AttributeError:
            return {
                'status': 'error',
                'message': f'Internal error: attribute {attr_name} not found'
            }
        
        # Convert value to appropriate type
        try:
            if isinstance(current_value, int):
                new_value = int(value)
            elif isinstance(current_value, float):
                new_value = float(value)
            else:
                new_value = value
        except (ValueError, TypeError) as e:
            return {
                'status': 'error',
                'message': f'Invalid value for {name}: {e}'
            }
        
        # Validate value (basic sanity checks)
        if isinstance(new_value, (int, float)):
            if new_value < 0 and 'cost' not in name and 'strength' not in name:
                return {
                    'status': 'error',
                    'message': f'Value for {name} cannot be negative: {new_value}'
                }
        
        # Set the new value
        setattr(cls, attr_name, new_value)
        
        return {
            'status': 'ok',
            'name': name,
            'old_value': current_value,
            'new_value': new_value
        }
    
    def set_param(self, name: str, value) -> dict:
        """
        Set a simulation parameter or execute a command.
        
        Args:
            name: Parameter name in dotted format (e.g., 'prey.move_cost') or command name
            value: New value for parameters, or command-specific data
        
        Commands:
            kill_agent: value=agent_id - despawn the agent with this id
            clone_agent: value=agent_id - spawn new agent near existing agent with same properties
            boost_agent: value=agent_id - double the agent's energy
            breeding_pair: value=[x, y, type] - spawn male and female of type (0=prey, 1=predator) at (x,y)
        
        Returns:
            {'status': 'ok', ...} on success
            {'status': 'deferred', ...} for deferred execution
            {'status': 'error', 'message': '...'} on failure
        """
        # Check if this is a class-level parameter
        if name in self._PARAM_MAP:
            return self._set_class_param(name, value)
        
        # Check if this is a command
        if name in self._COMMANDS:
            
            # Commands that modify the agents dict must be deferred during step()
            # to avoid "dictionary changed size during iteration" errors
            if name in ('kill_agent', 'clone_agent', 'breeding_pair') and self._stepping:
                self._deferred_commands.append((name, value))
                return {
                    'status': 'deferred',
                    'command': name,
                    'message': 'Command queued for execution after current step completes'
                }
            
            if name == 'kill_agent':
                try:
                    agent_id = int(value)
                except (ValueError, TypeError):
                    return {'status': 'error', 'message': f'kill_agent requires integer agent id, got: {value}'}
                
                if agent_id not in self.agents:
                    return {'status': 'error', 'message': f'Agent {agent_id} not found'}
                
                self.remove_agent(agent_id)
                return {'status': 'ok', 'command': name, 'agent_id': agent_id}
            
            elif name == 'clone_agent':
                try:
                    agent_id = int(value)
                except (ValueError, TypeError):
                    return {'status': 'error', 'message': f'clone_agent requires integer agent id, got: {value}'}
                
                if agent_id not in self.agents:
                    return {'status': 'error', 'message': f'Agent {agent_id} not found'}
                
                original = self.agents[agent_id]
                
                # Spawn new agent in proximity (random offset within 1 unit)
                offset_x = (random.random() - 0.5) * 2  # [-1, 1]
                offset_y = (random.random() - 0.5) * 2  # [-1, 1]
                new_x = max(0, min(self.width - 0.01, original.x + offset_x))
                new_y = max(0, min(self.height - 0.01, original.y + offset_y))
                
                # Create new agent with same properties
                new_agent = self.spawn_agent(
                    original.agent_type, 
                    new_x, 
                    new_y, 
                    original.energy,
                    parent_ids=(-1, -1)  # Clone is treated as founder
                )
                
                # Copy genetic properties from original
                new_agent.sex = original.sex
                new_agent.vision = original.vision
                new_agent.brutpflege = original.brutpflege
                new_agent.speed = original.speed
                new_agent.heat = original.heat
                new_agent.aggression = original.aggression
                new_agent.genes = original.genes.copy()
                new_agent.bearing = original.bearing
                
                return {
                    'status': 'ok', 
                    'command': name, 
                    'original_id': agent_id, 
                    'new_id': new_agent.id,
                    'position': (round(new_x, 2), round(new_y, 2))
                }
            
            elif name == 'boost_agent':
                try:
                    agent_id = int(value)
                except (ValueError, TypeError):
                    return {'status': 'error', 'message': f'boost_agent requires integer agent id, got: {value}'}
                
                if agent_id not in self.agents:
                    return {'status': 'error', 'message': f'Agent {agent_id} not found'}
                
                agent = self.agents[agent_id]
                old_energy = agent.energy
                agent.energy *= 2
                self._mark_dirty(agent_id)
                
                return {
                    'status': 'ok', 
                    'command': name, 
                    'agent_id': agent_id,
                    'old_energy': round(old_energy, 2),
                    'new_energy': round(agent.energy, 2)
                }
            
            elif name == 'breeding_pair':
                # value should be [x, y, type]
                try:
                    if not isinstance(value, (list, tuple)) or len(value) != 3:
                        raise ValueError("Expected [x, y, type]")
                    x, y, agent_type_val = value
                    x = float(x)
                    y = float(y)
                    agent_type_val = int(agent_type_val)
                    if agent_type_val not in (0, 1):
                        raise ValueError("type must be 0 (prey) or 1 (predator)")
                except (ValueError, TypeError) as e:
                    return {'status': 'error', 'message': f'breeding_pair requires [x, y, type]: {e}'}
                
                # Validate position
                if not (0 <= x < self.width and 0 <= y < self.height):
                    return {'status': 'error', 'message': f'Position ({x}, {y}) out of bounds'}
                
                if not self.is_passable(x, y):
                    return {'status': 'error', 'message': f'Position ({x}, {y}) is not passable (wall)'}
                
                agent_type = AgentType.PREY if agent_type_val == 0 else AgentType.PREDATOR
                initial_energy = self.PREY_INITIAL_ENERGY if agent_type == AgentType.PREY else self.PREDATOR_INITIAL_ENERGY
                
                # Spawn male at slight offset
                cell_x, cell_y = int(x), int(y)
                male_x = cell_x + random.random()
                male_y = cell_y + random.random()
                male = self.spawn_agent(agent_type, male_x, male_y, initial_energy)
                male.sex = 0  # Force male
                
                # Spawn female at different position in same cell
                female_x = cell_x + random.random()
                female_y = cell_y + random.random()
                female = self.spawn_agent(agent_type, female_x, female_y, initial_energy)
                female.sex = 1  # Force female
                
                return {
                    'status': 'ok',
                    'command': name,
                    'position': (cell_x, cell_y),
                    'type': agent_type_val,
                    'male_id': male.id,
                    'female_id': female.id
                }
        
        # Unknown parameter or command
        return {
            'status': 'error',
            'message': f'Unknown parameter or command: {name}',
            'valid_params': list(self._PARAM_MAP.keys()),
            'valid_commands': list(self._COMMANDS)
        }

# Module-level helpers
def _get_bearing(dx, dy):
    return math.atan2(dy, dx)

def _normalize_bearing(b):
    """Keep bearing in [-pi, pi]."""
    return (b + math.pi) % TWOPI - math.pi

def run_simulation(
    width: int = 100,
    height: int = 100,
    initial_prey: int = 1000,
    initial_predators: int = 15,
    max_ticks: int = 5000,
    seed: Optional[int] = None
):
    """Run the simulation and print statistics."""
    
    # Use World.create() for initialization
    world = World.create(
        width=width,
        height=height,
        initial_prey=initial_prey,
        initial_predators=initial_predators,
        seed=seed
    )
    
    print(f"Starting simulation: {width}x{height} world, {initial_prey} prey, {initial_predators} predators")
    print("-" * 50)
    
    for tick in range(max_ticks):
        world.step()
        prey_count, pred_count = world.count_agents()
        
        # Print status every 10 ticks
        if tick % 10 == 0:
            print(f"Tick {world.tick:4d}: Prey={prey_count:5d}, Predators={pred_count:4d}")
        
        # Check if simulation halted itself
        if world.halted:
            print(f"\n*** {world.halt_reason} at tick {world.tick} ***")
            break
    else:
        print(f"\n*** Simulation completed {max_ticks} ticks ***")
    
    prey_count, pred_count = world.count_agents()
    print(f"Final state: Prey={prey_count}, Predators={pred_count}")
    
    # Demo: show delta tracking
    print("\n--- Delta tracking demo ---")
    world.halted = False  # Reset for demo
    world.mark_clean()  # Clear dirty state
    world.step()  # Do one more step
    delta = world.get_dirty_state()
    print(f"After 1 tick: {len(delta['spawned'])} spawned, {len(delta['updated'])} updated, {len(delta['despawned'])} despawned")


if __name__ == "__main__":
    run_simulation(seed=42)