"""
Macau Casino-Tourism Integrated Simulation Model
Final Version 3.0 | Python 3.10+ | Updated: 2025-02-23
"""

import simpy
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import random
from mesa import Agent, Model
from mesa.space import MultiGrid, PropertyLayer
from mesa.time import RandomActivation
from deap import base, creator, tools, algorithms
from sklearn.multioutput import MultiOutputRegressor
from sklearn.ensemble import GradientBoostingRegressor
from datetime import datetime
import warnings

# Suppress warnings
warnings.filterwarnings("ignore", category=FutureWarning)

# Setting up random seed
random.seed(41)
np.random.seed(41)

# ======================
# IMPROVED CASINO MODEL
# ======================

class CasinoModel:
    """Optimized casino operations with anti-corruption factors"""
    
    def __init__(self, peak_hours):
        self.env = simpy.Environment()
        self.tables = simpy.Resource(self.env, capacity=50)
        self.peak_hours = peak_hours
        self.wait_times = []
        self.revenue = 0
        self.staff_cost = 0
        self.hourly_revenue = []
        
        # Anti-corruption impact factors
        self.corruption_factor = 0.71  # From research [2]
        
        self.arrival_rates = {
            'morning': 80 * self.corruption_factor,
            'afternoon': 120 * self.corruption_factor,
            'peak': 300 * self.corruption_factor,
            'night': 200 * self.corruption_factor
        }

    def player_generator(self):
        """Batch processing for efficiency"""
        player_id = 0
        while True:
            current_hour = self.env.now % 24
            rate = self._get_arrival_rate(current_hour)
            
            # Batch processing
            batch_size = 100
            yield self.env.timeout(np.mean(
                [np.random.exponential(60/rate) for _ in range(batch_size)]
            ))
            
            player_id += batch_size
            self.env.process(self.play_table(player_id))

    def _get_arrival_rate(self, hour):
        """Time-dependent arrival rates"""
        if hour in self.peak_hours:
            return self.arrival_rates['peak']
        elif 6 <= hour < 12:
            return self.arrival_rates['morning']
        elif 12 <= hour < 18:
            return self.arrival_rates['afternoon']
        else:
            return self.arrival_rates['night']

    def play_table(self, player_id):
        """Table gameplay process"""
        arr_time = self.env.now
        with self.tables.request() as request:
            yield request
            wait_time = self.env.now - arr_time
            self.wait_times.append(wait_time)
            
            game_type = np.random.choice(
                ['Baccarat', 'Blackjack', 'Roulette'],
                p=[0.6, 0.3, 0.1]
            )
            
            service_time = self._get_service_time(game_type)
            yield self.env.timeout(service_time)
            
            revenue = self._calculate_revenue(game_type)
            self.revenue += revenue
            self.hourly_revenue.append(revenue)

    def _get_service_time(self, game_type):
        """Service time distributions"""
        return {
            'Baccarat': np.random.triangular(15, 30, 45),
            'Blackjack': np.random.triangular(10, 25, 40),
            'Roulette': np.random.triangular(20, 35, 50)
        }[game_type]

    def _calculate_revenue(self, game_type):
        """Revenue calculation with anti-corruption impact"""
        return {
            'Baccarat': 1500 * (1 - 0.3 * self.corruption_factor),
            'Blackjack': 1200 * (1 - 0.25 * self.corruption_factor),
            'Roulette': 800 * (1 - 0.2 * self.corruption_factor)
        }[game_type]

    def run(self, sim_hours):
        """Execute simulation"""
        self.env.process(self.player_generator())
        self.env.run(until=sim_hours)
        
        self.staff_cost = 120 * (sim_hours//24) * 50
        return {
            'total_players': len(self.wait_times),
            'total_revenue': self.revenue,
            'staff_cost': self.staff_cost,
            'profit': self.revenue - self.staff_cost,
            'hourly_revenue': self.hourly_revenue
        }

# ======================
# TOURISM FLOW MODEL
# ======================

class TouristAgent(Agent):
    """Enhanced tourist agent with spending behaviors"""
    
    def __init__(self, unique_id, model):
        super().__init__(unique_id, model)
        self.type = np.random.choice(
            ['Gambler', 'Cultural', 'Family', 'CraftEnthusiast'],
            p=[0.4, 0.3, 0.2, 0.1]
        )
        self.initial_budget = {
            'Gambler': np.random.lognormal(8, 0.5),
            'Cultural': np.random.lognormal(7, 0.3),
            'Family': np.random.lognormal(7.5, 0.4),
            'CraftEnthusiast': np.random.lognormal(9, 0.6)
        }[self.type]
        self.current_budget = self.initial_budget
        self.total_spending = 0
        self.dwell_time = 0
        self.craft_interest = np.random.beta(2,5) if self.type == 'CraftEnthusiast' else 0

    def move(self):
        """Movement with craft cluster preference"""
        possible_steps = self.model.grid.get_neighborhood(
            self.pos,
            moore=True,
            include_center=False
        )
        
        if self.type == 'CraftEnthusiast':
            craft_clusters = []
            for pos in possible_steps:
                x, y = pos
                cell_data = self.model.cluster_layer.data[x][y]
                if cell_data['cluster_type'] == 'craft':
                    craft_clusters.append(pos)
            
            if craft_clusters:
                possible_steps = craft_clusters
                
        new_position = self.random.choice(possible_steps)
        self.model.grid.move_agent(self, new_position)

    def spend(self):
        """Spending behavior with budget constraints"""
        if self.type == 'Gambler':
            expenditure = min(np.random.exponential(500), self.current_budget)
        elif self.type == 'Cultural':
            expenditure = min(np.random.normal(300, 50), self.current_budget)
        elif self.type == 'Family':
            expenditure = min(np.random.gamma(2, 150), self.current_budget)
        else:
            expenditure = min(1500 + np.random.pareto(2.5), self.current_budget)
            
        self.current_budget -= expenditure
        self.total_spending += expenditure
        self.model.total_spending += expenditure

    def satisfaction_score(self):
        """Calculate tourist satisfaction"""
        return min(1.0, (self.total_spending / self.initial_budget) * 
                  (1 + 0.2 * self.dwell_time/24))

    def step(self):
        self.move()
        self.spend()
        self.dwell_time += 1

class TourismModel(Model):
    """Enhanced model with craft clusters"""
    
    def __init__(self, N, width, height):
        super().__init__()
        self.num_agents = N
        self.width = width
        self.height = height
        self.grid = MultiGrid(width, height, torus=False)
        self.schedule = RandomActivation(self)
        self.total_spending = 0
        
        # Initialize property layer
        self.cluster_layer = PropertyLayer(
            name="craft_clusters",
            width=width,
            height=height,
            default_value={'cluster_type': None, 'visitors': 0, 'revenue': 0},
            dtype=object
        )
        self.grid.add_property_layer(self.cluster_layer)
        
        self._init_craft_clusters()
        
        # Create agents
        for i in range(self.num_agents):
            agent = TouristAgent(i, self)
            self.schedule.add(agent)
            x = self.random.randrange(self.grid.width)
            y = self.random.randrange(self.grid.height)
            self.grid.place_agent(agent, (x, y))

    def _init_craft_clusters(self):
        """Initialize craft clusters"""
        cluster_positions = [
            (10, 10), (35, 10), (10, 40), (35, 40),
            (20, 25), (45, 30), (5, 35)
        ]
        
        for pos in cluster_positions:
            self.cluster_layer.set_cell(
                position=pos,
                value={
                    'cluster_type': 'craft',
                    'visitors': 0,
                    'revenue': 0
                }
            )

    def step(self):
        self.schedule.step()
        self._update_cluster_stats()

    def _update_cluster_stats(self):
        """Update cluster statistics"""
        for x in range(self.width):
            for y in range(self.height):
                cell_data = self.cluster_layer.data[x][y]
                if cell_data['cluster_type'] == 'craft':
                    agents_here = self.grid.get_cell_list_contents((x, y))
                    visitors = len([a for a in agents_here 
                                  if a.type == 'CraftEnthusiast'])
                    revenue = sum(a.total_spending * a.craft_interest 
                                for a in agents_here 
                                if a.type == 'CraftEnthusiast')
                    
                    self.cluster_layer.set_cell(
                        position=(x, y),
                        value={
                            'cluster_type': 'craft',
                            'visitors': visitors,
                            'revenue': revenue
                        }
                    )

# ======================
# OPTIMIZATION FRAMEWORK
# ======================

class Optimizer:
    """Multi-objective optimization with surrogate modeling"""
    
    def __init__(self):
        self.toolbox = base.Toolbox()
        self._setup_evolution()
        self.surrogate = MultiOutputRegressor(GradientBoostingRegressor())
        self.gen = 0
        # Define operational costs and revenue parameters
        self.cost_per_table = 5000  # Daily fixed cost per table
        self.cost_per_staff = 1200  # Daily cost per staff
        self.revenue_per_table = 15000  # Average daily revenue per table
        self.vip_multiplier = 1.5  # VIP revenue multiplier

    def _setup_evolution(self):
        creator.create("FitnessMulti", base.Fitness, weights=(-1.0, -1.0, -1.0))
        creator.create("Individual", list, fitness=creator.FitnessMulti)
        
        self.toolbox.register("attr_int", np.random.randint, 20, 100)
        self.toolbox.register("attr_float", np.random.uniform, 0.5, 2.0)
        self.toolbox.register("attr_marketing", np.random.uniform, 1.0, 5.0)
        # Custom crossover operator
        def custom_cx(ind1, ind2):
            # Integer crossover for tables and staff
            c1, c2 = tools.cxTwoPoint([ind1[0]], [ind2[0]])
            ind1[0], ind2[0] = int(c1[0]), int(c2[0])
            
            c1, c2 = tools.cxTwoPoint([ind1[1]], [ind2[1]])
            ind1[1], ind2[1] = int(c1[0]), int(c2[0])
            
            # Blend crossover for marketing budget
            alpha = 0.5
            gamma = (1. + 2. * alpha) * random.random() - alpha
            ind1[2] = (1. - gamma) * ind1[2] + gamma * ind2[2]
            ind2[2] = gamma * ind1[2] + (1. - gamma) * ind2[2]
            
            return ind1, ind2
        self.toolbox.register("individual", tools.initCycle, creator.Individual,
                             (self.toolbox.attr_int, self.toolbox.attr_int, 
                              self.toolbox.attr_float), n=1)
        self.toolbox.register("population", tools.initRepeat, list, 
                            self.toolbox.individual)
        self.toolbox.register("attr_table", np.random.randint, 50, 150)
        self.toolbox.register("attr_staff", np.random.randint, 70, 300)
        self.toolbox.register("mate", tools.cxBlend, alpha=0.5)
        self.toolbox.register("mutate", tools.mutGaussian, mu=0, sigma=1, indpb=0.2)
        self.toolbox.register("select", tools.selNSGA2)

    def _calculate_profit(self, tables, staff, marketing):
        """Realistic profit calculation"""
        # Revenue calculation
        base_revenue = tables * self.revenue_per_table * 7  # Weekly revenue
        vip_revenue = base_revenue * (marketing/5.0) * self.vip_multiplier
        total_revenue = base_revenue + vip_revenue
        
        # Cost calculation
        table_costs = tables * self.cost_per_table * 7
        staff_costs = staff * self.cost_per_staff * 7
        marketing_costs = marketing * 1e6  # Convert to actual cost
        
        return total_revenue - (table_costs + staff_costs + marketing_costs)
    
    def _population_similarity(self, population):
        """Calculate population similarity based on fitness values"""
        if not population:
            return 0
        
        # Get fitness values for valid individuals
        fitnesses = [ind.fitness.values[0] for ind in population 
                    if ind.fitness.valid]
        
        if not fitnesses:
            return 0
            
        # Calculate similarity based on fitness variance
        mean_fitness = sum(fitnesses) / len(fitnesses)
        variance = sum((f - mean_fitness) ** 2 for f in fitnesses) / len(fitnesses)
        
        # Return normalized similarity score (0-1)
        similarity = 1.0 / (1.0 + variance)
        return similarity

    def _update_population(self, population):
        """Maintain population diversity"""
        if self._population_similarity(population) > 0.9:  # Too similar
            num_random = len(population) // 4
            population[-num_random:] = self.toolbox.population(n=num_random)

    def run_optimization(self, ngen=40):
        pop = self.toolbox.population(n=50)
        hof = tools.ParetoFront()

        # Initial evaluation
        fitnesses = list(map(self._true_evaluation, pop))
        for ind, fit in zip(pop, fitnesses):
            ind.fitness.values = fit

        # Train surrogate
        X_train = [list(ind) for ind in pop]
        y_train = [ind.fitness.values for ind in pop]
        self.surrogate.fit(X_train, y_train)

        for self.gen in range(1, ngen+1):
            offspring = tools.selTournament(pop, len(pop), tournsize=3)
            offspring = list(map(self.toolbox.clone, offspring))

            # Apply variation
            for child1, child2 in zip(offspring[::2], offspring[1::2]):
                if random.random() < 0.8:
                    self.toolbox.mate(child1, child2)
                    del child1.fitness.values
                    del child2.fitness.values

            for mutant in offspring:
                if random.random() < 0.2:
                    self.toolbox.mutate(mutant)
                    del mutant.fitness.values

            # Evaluate
            invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
            if self.gen % 10 == 0:
                fitnesses = list(map(self._true_evaluation, invalid_ind))
                X_train.extend([list(ind) for ind in invalid_ind])
                y_train.extend(fitnesses)
                self.surrogate.fit(X_train, y_train)
            else:
                fitnesses = self.surrogate.predict([list(ind) for ind in invalid_ind])
                for ind, fit in zip(invalid_ind, fitnesses):
                    ind.fitness.values = tuple(fit)
            # Add diversity maintenance
            self._update_population(pop)
            pop[:] = offspring
            hof.update(pop)
            print(f"Generation {self.gen} completed")
            # Set the staff and tables to integer values
            for ind in pop:
                ind[0], ind[1] = map(int, [ind[0], ind[1]])
            for ind in hof:
                ind[0], ind[1] = map(int, [ind[0], ind[1]])
            print(f"Best individual: {hof[0]}")
        return pop, hof

    def _true_evaluation(self, individual):
        tables, staff, marketing = individual
        penalty = 0

         # Staff-to-table ratio constraints
        min_staff = tables * 1.5  # Minimum 1.5 staff per table
        max_staff = tables * 3    # Maximum 3 staff per table
        
        if staff < min_staff or staff > max_staff:
            penalty += 1e6 * (min(abs(staff - min_staff), abs(staff - max_staff)))
        
        # Calculate profit with noise
        profit = self._calculate_profit(tables, staff, marketing)
        noise = np.random.uniform(0.95, 1.05)

        casino = CasinoModel(peak_hours=[13,14,21,22,23])
        casino_results = casino.run(24*7)
        
        tourism = TourismModel(N=int(marketing*1000), width=50, height=50)
        for _ in range(24*7):
            tourism.step()
            
        return (
            -(profit * noise - penalty),
            -tourism.total_spending * noise,
            staff * self.cost_per_staff*7
        )

# ======================
# VISUALIZATION & ANALYSIS
# ======================

def plot_results(casino_data, tourism_data):
    """Generate analysis visualizations"""
    fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(15, 10))
    
    # Casino performance
    ax1.plot(casino_data['hour'], casino_data['hourly_revenue'], label='Revenue')
    ax1.plot(casino_data['hour'], casino_data['hourly_cost'], label='Cost')
    ax1.set_title('Hourly Casino Performance')
    ax1.set_xlabel('Hour')
    ax1.set_ylabel('HKD')
    ax1.legend()
    
    # Wait times
    sns.histplot(casino_data['wait_times'], bins=30, ax=ax2)
    ax2.set_title('Wait Time Distribution')
    ax2.set_xlabel('Minutes')
    
    # Tourist types
    tourism_data['type'].value_counts().plot(kind='bar', ax=ax3)
    ax3.set_title('Tourist Type Distribution')
    
    # Spending patterns
    sns.boxplot(x='type', y='initial_budget', data=tourism_data, ax=ax4)
    ax4.set_title('Initial Budget by Tourist Type')
    
    plt.tight_layout()
    plt.savefig('casino_tourism_analysis.png')
    plt.show()

# ======================
# MAIN EXECUTION
# ======================

def main():
    """Main simulation workflow"""
    print("Macau Casino-Tourism Simulation")
    print("===============================")
    
    try:
        # Configuration
        days = 7
        num_tourists = 10000
        grid_size = 100
        generations = 10
        
        # Initialize models
        casino = CasinoModel(peak_hours=[13,14,21,22,23])
        optimizer = Optimizer()
        
        # Run casino simulation
        print("\nRunning casino simulation...")
        casino_results = casino.run(24 * days)
        
        # Prepare casino data with length verification
        hours = range(24 * days)
        hourly_revenue = casino_results['hourly_revenue']
        wait_times = casino.wait_times
        hourly_cost = [casino_results['staff_cost']/(24*days)] * len(hours)

        # Ensure all arrays have the same length by truncating to minimum length
        min_length = min(len(hours), len(hourly_revenue), len(wait_times))

        casino_hourly = pd.DataFrame({
            'hour': list(hours)[:min_length],
            'hourly_revenue': hourly_revenue[:min_length],
            'hourly_cost': hourly_cost[:min_length],
            'wait_times': wait_times[:min_length]
        })

        
        # Run tourism simulation
        print("\nRunning tourism simulation...")
        tourism = TourismModel(N=num_tourists, width=grid_size, height=grid_size)
        tourism_data = []
        for agent in tourism.schedule.agents:
            tourism_data.append({
                'type': agent.type,
                'initial_budget': agent.initial_budget,
                'total_spending': agent.total_spending,
                'dwell_time': agent.dwell_time,
                'satisfaction': agent.satisfaction_score()
            })
        tourism_df = pd.DataFrame(tourism_data)
        
        # Run optimization
        print(f"\nOptimizing over {generations} generations...")
        pop, hof = optimizer.run_optimization(ngen=generations)
        
        # Generate outputs
        plot_results(casino_hourly, tourism_df)
        
        print("\nOptimization Results:")
        for idx, sol in enumerate(hof[:3]):
            print(f"Solution {idx+1}:")
            print(f"  Tables: {sol[0]}")
            print(f"  Staff: {sol[1]}")
            print(f"  Marketing: {sol[2]:.2f}M HKD")
            print(f"  Expected Profit: {-sol.fitness.values[0]/1e6:.2f}M HKD\n")
            
    except Exception as e:
        print(f"\nSimulation error: {str(e)}")
        raise

if __name__ == "__main__":
    main()