
# ============================================================
# DESCRIPTION DU FICHIER
# Ce fichier permet de simuler le régime actuel avec des paramètres choisis
# et de comparer le résultat avec les statistiques cibles
# il est le fichier central pour la calibration du modèle simulé
# Fichiers appellés : parametres_regime_reference.py
# Utilisé par : simulations_reformes.ipynb (partie Calibration)
# Fichiers produits : fichiers csv optionels retracant les statistiques agrégées simulées mois par mois
# ============================================================

# packages
import numpy as np
from scipy.optimize import minimize
import pandas as pd
import os
os.getcwd()
os.chdir('/Codes_note_90/modele_simule')# à adapter en fonction de la localisation de vos codes
from parametres_regime_reference import * # fichier regroupant les paramètres utilisés et les statistiques cibles
from concurrent.futures import ProcessPoolExecutor, as_completed
from numpy.random import SeedSequence, default_rng
import math

# ============================================================
# CLASSE INDIVIDU
# Représente un agent dans la simulation du marché du travail.
# Chaque individu possède un statut (employé, chômeur indemnisé,
# chômeur non indemnisé), une éligibilité à l'assurance chômage,
# des taux de transition emploi/chômage, et un historique de
# mois travaillés sur les 72 derniers mois.
# ============================================================
class Individual:

    # ------------------------------------------------------------------
    # Fonction intialisant les caractéristiques de l'individu
    # ------------------------------------------------------------------
    def __init__(self, id, status='employed', eligible=False, job_finding_rate=0.1, job_destruction_rate=0.01, take_up_rate=0.5, months_worked=None, months_ui_left=0,take_up=0,jdr_coef=1.0, rngs=None):

        # Initialisation de l'historique de travail : 72 mois à 0 (= non travaillés) par défaut
        if months_worked is None:
            months_worked = [0] * 72 
        self.id = id
        self.status = status           # Statut courant : 'employed', 'ui_unemployed', 'noui_unemployed' (emploi, chomage indemnisé, chomage non indemnisé)
        self.eligible = eligible       # Éligibilité à l'assurance chômage 
        self.months_worked = months_worked  # Historique binaire d'emploi des 72 derniers mois (1=travaillé, 0=non travaillé)
        self.tot_months_worked = sum(self.months_worked[0:pra_max])  # Nombre total de mois travaillés sur la période de référence
        self.months_ui_left = months_ui_left  # Mois de droits au chômage restants
        self.take_up = take_up         # Indicateur de recours effectif à l'assurance chômage (1=recours, 0=non-recours)
        self.take_up_rate = take_up_rate  # Taux de recours à l'assurance chômage (probabilité)
        self.jdr_coef = jdr_coef       # Coefficient de sensibilité du taux de destruction d'emploi à l'ancienneté
        self.delai_months_ui_left = 0  # Délai restant avant déchéance des droits 

        # Marqueurs de flux mensuels pour calculer les stocks/flux agrégés (remis à 0 à chaque mise à jour)
        self.flow_into_ui = 0     # Indicateur d'une entrée en chômage indemnisé (ouverture de droits initiale)
        self.flow_out_ui = 0      # Indicateur d'une sortie du chômage indemnisé (retour à l'emploi)
        self.flow_recharg = 0     # Indicatuer d'un rechargement des droits au chômage
        self.flow_reprise = 0     # Indicateur d'une reprise de droits préexistants après perte d'emploi
        self.odd = 0              # Indicateur d'une ouverture de droits 
        self.months_odd = 0       # Nombre de mois travaillés ayant donné lieu à l'ouverture de droits
        self.epuise_uni = 0       # Indicateur d'épuisement des droits 
        

        # Générateurs de nombres aléatoires indépendants par type d'événement
        # (permet la reproductibilité du code et l'indépendance des chocs)
        self.rngs = rngs or { 
            'destruction': np.random.default_rng(),  # Destruction d'emploi
            'finding': np.random.default_rng(),       # Retour à l'emploi
            'system': np.random.default_rng(),        # Renouvellement de la population
            'take_up': np.random.default_rng(),       # Décision de recours à l'assurance chômage 
        }
        self.job_finding_rate = job_finding_rate        # Taux mensuel de retour à l'emploi
        self.job_destruction_rate = job_destruction_rate  # Taux mensuel de destruction d'emploi 

    # ------------------------------------------------------------------
    # Fonction metant à jour l'état de l'individu et ses caractéristiques (historique d'emploi, droits restants) chaque mois en     # fonction des probabilités de transition, qui sont elles mêmes dépendantes de ses caractéristiques (état en début de mois,     # historique d'emploi, droits restants), et des règles de l'assurance chômage (éligibilité, durée des droits etc). 
    # ------------------------------------------------------------------
    def update_status(self):
        #Met à jour le statut de l'individu pour un nouveau mois de simulation.

        # Remise à zéro des marqueurs de flux mensuels
        self.flow_into_ui = 0
        self.flow_recharg = 0
        self.flow_reprise = 0
        self.flow_out_ui = 0
        self.odd = 0
        self.months_odd = 0
        self.epuise_uni = 0

        self.tot_months_worked = sum(self.months_worked[0:pra_max])
        
        # Décompte du délai avant déchéance des droits (si non nul)
        if self.delai_months_ui_left >= 1:
            self.delai_months_ui_left -= 1

        # Renouvellement de la population : avec probabilité flow_system,
        # l'individu est "réinitialisé" (simule les entrées/sorties du marché du travail, calibré de manière approximative avec un renouvellement en 44 ans (62-18)
        if self.rngs['system'].random() < flow_system:
            self.months_ui_left = 0
            self.months_worked = [0] * 72  
            self.take_up = 0
            # Un chômeur indemnisé affecté par la réinitialisation passe au chômage non indemnisé (car un nouvel entrant sur le marché du travail n'a pas de droits au chômage restants)
            if self.status == 'ui_unemployed':
                self.status = 'noui_unemployed'

        # ----------------------------------------------------------------
        # Cas 1 : Individu en emploi
        # ----------------------------------------------------------------
        if self.status == 'employed':
            self.take_up = 0
            # Si le délai de déchéance est écoulé, les droits restants sont annulés
            if self.delai_months_ui_left == 0:
                self.months_ui_left = 0
            # Mise à jour de l'historique : marque le mois courant comme travaillé car l'individu est en emploi
            self.months_worked = [1] + self.months_worked[:-1]
            self.tot_months_worked = sum(self.months_worked[0:pra_max]) # Calcul du nombre total de mois travaillés sur la période de référence
            # Calcul de l'ancienneté continue (nombre de mois travaillés sans interruption)
            months_worked_continu = min(self.months_worked.index(0), 72) if 0 in self.months_worked else 72 #plafonne à la durée de la mémoire (= 5 ans)
            # Vérification de l'éligibilité à une indemnisation de l'assurance chômage
            if self.tot_months_worked >= threshold_UI or self.months_ui_left > 0: # Si l'individu a travaillé suffisament ou droit restant
                self.eligible = True
                # Taux de destruction d'emploi individualisé et décroissant avec l'ancienneté (jdr_coef est < 0)
                jdrate = self.job_destruction_rate*np.exp(self.jdr_coef*months_worked_continu) 
            else: # Si l'individu n'a pas de droit restant et n'a pas assez travaillé pour êtré éligible à un nouveau droit
                self.eligible = False
                # Taux de destruction d'emploi individualisé et décroissant avec l'ancienneté (jdr_coef est < 0), ce taux est calculé de la même façon que pour un individu éligible 
                jdrate = self.job_destruction_rate*np.exp(self.jdr_coef*months_worked_continu) 

            # Tirage du choc de destruction d'emploi
            if self.rngs['destruction'].random() < jdrate: # L'individu perd son emploi avec une probabilité jdrate
                if self.eligible: # Si l'individu est éligible à l'assurance chômage
                    self.status = 'noui_unemployed' 
                    if self.months_ui_left == 0 and self.rngs['take_up'].random() < self.take_up_rate:
                        # Ouverture de droits : l'individu n'avait pas de droits résiduels
                        # et décide de recourir à l'assurance chômage (probabilité self.take_up_rate)
                        self.flow_into_ui = 1
                        self.status = 'ui_unemployed'
                        self.odd = 1
                        self.take_up = 1
                        self.months_odd = self.tot_months_worked
                        months_worked_utile = self.months_worked[0:pra_max] # vecteur des périodes travaillées sur la période de référence
                        # Longeur de la PRC (durée entre le premier et dernier mois travaillé sur la période de référence) 
                        length_prc = (pra_max - months_worked_utile[::-1].index(1)) if 1 in self.months_worked else 0
                        # Calcul des droits ouverts : entre 6 mois et pbd_max, selon la formule de conversion
                        rights = min(pbd_max, max(round(min(self.tot_months_worked*1.7, length_prc)*coef_conversion_UI),6))
                        self.months_ui_left = rights
                        self.months_droits = rights
                        # Le délai de déchéance court sur 36 mois + durée des droits
                        self.delai_months_ui_left = 36 + rights
                        self.months_worked = [0]*72 
                    else:
                        # Reprise de droits résiduels si l'individu en avait encore
                        if self.months_ui_left > 0:
                            self.status = 'ui_unemployed'
                            self.flow_reprise = 1
                            self.take_up = 1
                else:
                    # Perte d'emploi sans éligibilité : chômage non indemnisé
                    self.status = 'noui_unemployed'
                    self.take_up = 0

        # ----------------------------------------------------------------
        # Cas 2 : Individu au chômage (indemnisé ou non)
        # ----------------------------------------------------------------
        else:
            self.months_worked = [0] + self.months_worked[:-1] # marque le mois courant comme non travaillé
            self.tot_months_worked = sum(self.months_worked[0:pra_max]) # calcule le nombre de mois travaillés sur la période de référence
            
            # Sous-cas 2a : Chômeur indemnisé (recours effectif à l'assurance chômage)
            if self.take_up > 0:
                self.months_ui_left = max(self.months_ui_left - 1, 0) # consommation d'un mois de droits
                if self.rngs['finding'].random() < self.job_finding_rate: #si retrouve un emploi
                    self.status = 'employed'
                    self.take_up = 0
                    self.flow_out_ui = 1
                else: # si ne retrouve pas d'emploi
                    if self.months_ui_left == 0: # si est arrivé au bout de son droit
                        # Vérification de la durée retravaillée utilisable pour un rechargement :
                        if sum(self.months_worked) == 0 or self.months_worked.index(1) >= 48: # si pas retravaillé 
                            self.tot_months_worked_recharg = 0
                        else:
                            # Calcul des mois travaillés sur la période de référence du rechargement
                            self.tot_months_worked_recharg = sum(self.months_worked[self.months_worked.index(1):self.months_worked.index(1)+pra_max]) # regarde le nombre de mois travaillés sur la PRA associée
                            
                        if self.tot_months_worked_recharg >= threshold_UI: # si a retravaillé au moins 6 mois => ouvre un nouveau droit (on considère un taux de recours de 1 ici notamment grace aux dispositifs automatiques de rappel)
                            # Rechargement des droits au chômage 
                            self.flow_recharg = 1 
                            self.odd = 1
                            self.months_odd = self.tot_months_worked_recharg
                            months_utiles = self.months_worked[self.months_worked.index(1):self.months_worked.index(1)+pra_max]
                            length_prc = (pra_max - months_utiles[::-1].index(1)) if 1 in self.months_worked else 0
                            # Calcul des nouveaux droits rechargés
                            rights = min(pbd_max, max(round(min(self.tot_months_worked_recharg*1.7, length_prc)*coef_conversion_UI),6))
                            self.months_ui_left = rights
                            self.months_droits = rights
                            self.delai_months_ui_left = 36 + rights
                            self.months_worked = [0]*72 #36 
                        else: # si n'a pas suffisament retravaillé: sortie du chômage pour épuisement de droits, passe au chômage non indemnisé
                            self.epuise_uni = 1
                            self.status = 'noui_unemployed'
                            self.take_up = 0
                            self.eligible = False
            # Sous-cas 2b : Chômeur non indemnisé
            else:
                if self.rngs['finding'].random() < self.job_finding_rate: # Si il retrouvent un emploi (probabilité suposée identique à celle d'un chômeur indemnisé dans ce modèle)
                    self.status = 'employed'
                    self.take_up = 0
                # sinon pas de changement d'état

                    


# ===================================================================
# CLASSE SIMULATION
# Orchestre la simulation sur un ensemble d'individus et un nombre
# de mois donné. Agrège les statistiques mensuelles dans un DataFrame.
# ===================================================================
  

class Simulation:
    # ------------------------------------------------------------------
    # Fonction initialisant les individus
    # ------------------------------------------------------------------
    def __init__(self, num_individuals=100, months=months_sim,
                 job_finding_rate=0.1, take_up_rate=0.3, job_destruction_rate=0.01,
                 jdr_coef=1.0, base_seed=None): 

        # Utilisation d'une SeedSequence pour garantir la reproductibilité
        # et l'indépendance des flux aléatoires entre individus
        ss_master = SeedSequence(base_seed)
        individual_seeds = ss_master.spawn(num_individuals)

        self.num_individuals = num_individuals
        self.months = months
        self.individuals = []

        for i in range(num_individuals):
            # Création de 4 flux RNG indépendants par individu
            ss_ind = individual_seeds[i].spawn(4)
            rngs = {
                'destruction': default_rng(ss_ind[0]),
                'finding': default_rng(ss_ind[1]),
                'system': default_rng(ss_ind[2]),
                'take_up': default_rng(ss_ind[3]),
            }

            ind = Individual(
                id=i,
                job_finding_rate=job_finding_rate,
                take_up_rate=take_up_rate,
                job_destruction_rate=job_destruction_rate,
                jdr_coef=jdr_coef,
                rngs=rngs
            ) 
            self.individuals.append(ind)

        # DataFrame de stockage de l'historique de la simulation
        self.history = pd.DataFrame(columns=['month', 'employed', 'ui_unemployed', 'noui_unemployed'])
        
    # ------------------------------------------------------------------
    # fonction exécutant la simulation chaque mois et enregistrant les statistiques agrégées
    # ------------------------------------------------------------------
    def run(self):
        
        records = []
        
        # Execute la simulation sur tous les mois
        for month in range(self.months):
            # ----------------------------------------------------------
            # Calcul des stocks et flux agrégés en début de mois
            # ----------------------------------------------------------
            employed = sum(1 for ind in self.individuals if ind.status == 'employed')
            employed_eligible = sum(1 for ind in self.individuals if ind.status == 'employed' and ind.eligible)
            ui_unemployed = sum(1 for ind in self.individuals if ind.status == 'ui_unemployed')
            noui_unemployed = sum(1 for ind in self.individuals if ind.status == 'noui_unemployed')
            noui_eligible = sum(1 for ind in self.individuals if ind.status == 'noui_unemployed' and ind.eligible)
            # Employés ayant entre 6 et 8 mois d'ancienneté (segment d'intérêt pour le calage)
            employed_6_8_history = sum(1 for ind in self.individuals if ind.status == 'employed' and 6 <= ind.tot_months_worked < 8)
            # Agrégation des flux mensuels
            flow_out_ui = sum(ind.flow_out_ui for ind in self.individuals)
            flow_into_ui = sum(ind.flow_into_ui for ind in self.individuals)
            flow_recharg = sum(ind.flow_recharg for ind in self.individuals)
            flow_reprise = sum(ind.flow_reprise for ind in self.individuals)
            epuise_uni = sum(ind.epuise_uni for ind in self.individuals)
            odd = sum(ind.odd for ind in self.individuals)  # Total des ouvertures/rechargements de droits
            # Part des ouvertures de droits pour les individus ayant 6-8 mois travaillés
            odd_6_8 = sum(1 for ind in self.individuals if ind.odd == 1 and 6 <= ind.months_odd < 8)
            # Part des ouvertures de droits conduisant à une durée de 15-18 mois
            odd_pbd_15_18 = sum(1 for ind in self.individuals if ind.odd == 1 and 15 < ind.months_droits <= 18)
            
            duree_droits_recharg = sum(ind.months_droits for ind in self.individuals if ind.flow_recharg == 1) 
            duree_droits_init = sum(ind.months_droits for ind in self.individuals if ind.flow_into_ui == 1)
            
            tax_rev = tax * employed # cotisations
            ui_cost = ui_unemployed*b # allocations
            rsa_cost = b_rsa * noui_unemployed * rsa_eligibility_rate # minima sociaux
            net_budget = tax_rev - ui_cost - rsa_cost # dépense mensuelle nette du système

            # Enregistrement des statistiques du mois
            records.append({
                'month': month,
                'employed': employed,
                'employed_eligible': employed_eligible,
                'ui_unemployed': ui_unemployed,
                'noui_unemployed': noui_unemployed,
                'noui_eligible': noui_eligible,
                'employed_6_8_history': employed_6_8_history,
                'flow_out_ui': flow_out_ui,
                'flow_into_ui': flow_into_ui,
                'flow_recharg': flow_recharg,
                'flow_reprise': flow_reprise,
                'odd': odd,
                'odd_6_8': odd_6_8,
                'odd_pbd_15_18': odd_pbd_15_18,
                'epuise_uni': epuise_uni,
                'duree_droits_recharg': duree_droits_recharg,
                'duree_droits_init':duree_droits_init,
                'tax_rev': tax_rev,
                'ui_cost': ui_cost,
                'rsa_cost': rsa_cost,
                'net_budget': net_budget
            })
            # mise à jour des individus
            for ind in self.individuals:
                ind.update_status()
        self.history = pd.concat([self.history, pd.DataFrame(records)], ignore_index=True)

    def get_history(self):
        return self.history

# -------------------------------------------------------------------
# Fonction exécutant la simulation pour une population d'individus donnée
# -------------------------------------------------------------------

def run_one(sim_id, num_individuals=10000, months=months_sim, seed=None,
            job_finding_rate=0.1, take_up_rate=0.9, job_destruction_rate=0.01,
            jdr_coef=0): 
    #Lance une simulation unique identifiée par sim_id et retourne son historique.

    sim = Simulation(num_individuals=num_individuals, months=months,
                     job_finding_rate=job_finding_rate, take_up_rate=take_up_rate,
                     job_destruction_rate=job_destruction_rate,
                     jdr_coef=jdr_coef, base_seed=seed) 
    sim.run()
    hist = sim.get_history().copy()
    hist["run_id"] = sim_id  # Ajout de l'identifiant du thread pour le suivi
    return hist

# -------------------------------------------------------------------
# Fonction de simulation parallélisée 
# Permet de paralléliser la simulation sur plusieurs cœurs
# pour réduire le temps de calcul sur de grandes populations.
# -------------------------------------------------------------------
def run_all(job_finding_rate=0.1, job_destruction_rate=0.01, take_up_rate=0.9, months=months_sim,
            jdr_coef=0): 
    #Lance nbr_threads simulations en parallèle et agrège les résultats par mois.

    # Répartition de la population totale entre les threads
    tot_pop = tot_pop_guess
    X = nbr_threads
    num_individuals_thread = int(tot_pop / X)
    # Graines différentes par thread pour garantir l'indépendance des simulations
    base_seed = 12345
    seeds = [base_seed + i * 1000 for i in range(X)]

    futures = []
    results = []
    # Soumission des tâches au pool de processus
    with ProcessPoolExecutor(max_workers=X) as ex:
        for i in range(X):
            futures.append(ex.submit(run_one, i, num_individuals_thread, months, seeds[i],
                                     job_finding_rate, take_up_rate, job_destruction_rate,
                                     jdr_coef)) 
        # Collecte des résultats au fur et à mesure de leur complétion
        for f in as_completed(futures):
            results.append(f.result())

    # Agrégation par mois : somme des stocks/flux sur tous les threads
    avg_history = pd.concat(results).groupby('month').sum().reset_index()
    return avg_history

# ------------------------------------------------------------------
# Fonction de calcul des statistiques agréees :
# Calcule les statistiques agrégées issues la simulation sur les 12 derniers mois.
# Ces statistiques servent à caler les paramètres du modèle sur des données cibles.
# Exportation optionnelle vers des fichiers CSV.
# ------------------------------------------------------------------
def moments_computation(avg_history, file1=None, file2=None):

    # Restriction aux 12 derniers mois (état stationnaire supposé atteint)
    avg_history_last_12 = avg_history[avg_history['month'] >= (months_sim - 12)].copy()
    # stocks/flux (moyennes sur les 12 derniers mois)
    employed_avg = avg_history_last_12['employed'].mean()
    employed_eligible_avg = avg_history_last_12['employed_eligible'].mean()
    ui_unemployed_avg = avg_history_last_12['ui_unemployed'].mean()
    noui_unemployed_avg = avg_history_last_12['noui_unemployed'].mean()
    noui_eligible_avg = avg_history_last_12['noui_eligible'].mean()
    employed_6_8_history = avg_history_last_12['employed_6_8_history'].mean()
    flow_out_ui_avg = avg_history_last_12['flow_out_ui'].mean()
    flow_into_ui_avg = avg_history_last_12['flow_into_ui'].mean()
    flow_recharg_avg = avg_history_last_12['flow_recharg'].mean()
    flow_reprise_avg = avg_history_last_12['flow_reprise'].mean()
    odd_6_8_hist = avg_history_last_12['odd_6_8'].mean()
    odd_pbd_15_18_hist = avg_history_last_12['odd_pbd_15_18'].mean()
    odd_hist = avg_history_last_12['odd'].mean()
    epuise_uni_avg = avg_history_last_12['epuise_uni'].mean()
    
    duree_droits_recharg_avg = avg_history_last_12['duree_droits_recharg'].mean()/flow_recharg_avg
    duree_droits_init_avg = avg_history_last_12['duree_droits_init'].mean()/flow_into_ui_avg
    # Liste ordonnée des moments retournés (utilisée dans objective_function)
    stocks = [
        employed_avg, ui_unemployed_avg, noui_unemployed_avg, noui_eligible_avg, employed_eligible_avg,
        employed_6_8_history, flow_out_ui_avg, flow_into_ui_avg, flow_recharg_avg, flow_reprise_avg,
        odd_6_8_hist, odd_pbd_15_18_hist, odd_hist, epuise_uni_avg, duree_droits_recharg_avg, duree_droits_init_avg
    ]

    # exports optionels vers des fichiers csv externes
    if file1 is not None:
        # Exportation de la moyenne des moments à l'échelle de l'économie réelle
        s_last = (avg_history_last_12.mean(numeric_only=True) * economy_scale)
        
        s_last.to_csv(file1, index=False)

    if file2 is not None:
        # Exportation de la série mensuelle complète à l'échelle de l'économie réelle
        s_last_full = (avg_history_last_12 * economy_scale)
        
        s_last_full.to_csv(file2, index=False)
    return stocks

# ------------------------------------------------------------------
# Fonction objectif pour le calage des paramètres structurels du modèle. 
# C'est cette fonction qui est appellée par le notebook
# Mesure l'écart quadratique entre les statistiques simulés et les statistiques cibles.
# Paramètres intervenant dans l'optimisation : taux de retour à l'emploi, taux de destruction d'emploi
# (niveau et sensibilité à l'ancienneté) et taux de recours.
# Moments cibles intervenant dans l'optimisation : Nombre d'individus au chômage indemnisé (pénalisé plus fortement,
# Nombre d'ouvertures de droits initiales, nombre de sorties vers l'emploi, nombre de reprises de droits et part des ouvertures # de droits avec 6 à 8 mois de durée d'affiliation (utilisé pour améliorer la précision du calage pour les réformes de         # l'éligibilité
# ------------------------------------------------------------------

def objective_function(params, file1=None, file2=None):

    jfr, jdr, jdr_coef, take_up_r= params 
    job_finding_rate = np.exp(jfr) # taux de retour à l'emploi
    job_destruction_rate = np.exp(jdr)# delta 0 (partie constante du taux de destruction d'emploi)
    jdr_coef =  jdr_coef # delta 1 (partie du taux de destruction d'emploi dépendante de l'ancienneté)
    take_up_rate = take_up_r # taux de recours
    
    
    # Lancement de la simulation avec les paramètres courants
    avg_history = run_all(job_finding_rate=job_finding_rate, job_destruction_rate=job_destruction_rate,
                           take_up_rate=take_up_rate, jdr_coef = jdr_coef) 
    stocks = moments_computation(avg_history, file1, file2)
    
    # Déballage des moments simulés
    (employed_avg, ui_unemployed_avg, noui_unemployed_avg, noui_eligible_avg, employed_eligible_avg,
     employed_6_8_history, flow_out_ui_avg, flow_into_ui_avg, flow_recharg_avg, flow_reprise_avg,
     odd_6_8_hist, odd_pbd_15_18_hist, odd_hist, epuise_uni_avg, duree_droits_recharg_avg, duree_droits_init_avg) = stocks
    
    # Calcul des indicateurs dérivés utilisés dans la fonction objectif
    flow_into_UI = flow_into_ui_avg
    flow_out_UI = flow_out_ui_avg
    # Taux de chômage global (indemnisé + non indemnisé)
    u_rate = (ui_unemployed_avg + noui_unemployed_avg) / max((employed_avg + ui_unemployed_avg + noui_unemployed_avg), 1e-12)
    # Part des employés avec 6-8 mois d'ancienneté parmi les employés
    share_68 = employed_6_8_history / max(employed_avg, 1e-12)
    # Part des ouvertures de droits pour les individus 6-8 mois parmi toutes les ouvertures
    share_68_odd = odd_6_8_hist / max(odd_hist, 1e-12)
    # Part des ouvertures conduisant à une durée de 15-18 mois
    share_pbd_15_18_odd = odd_pbd_15_18_hist/ max(odd_hist, 1e-12)
    # Taux d'épuisement des droits ARE
    rate_epuise_uni = epuise_uni_avg / max(ui_unemployed_avg, 1e-12)
    epuise_uni = epuise_uni_avg
    
    # Flux total d'entrée en ARE (première ouverture + reprises de droits)
    flow_ui_full = flow_into_UI + flow_reprise_avg
    flow_ui_full_target = target_UI_inflow + target_reprise_inflow

    # ------------------------------------------------------------------
    # Calcul de la fonction objectif :
    # Somme des écarts quadratiques normalisés entre moments simulés et cibles
    # ------------------------------------------------------------------
    obj_value = 0.0
       
    obj_value += ((ui_unemployed_avg - target_UI_stock)/target_UI_inflow)**2 #
    obj_value += ((flow_into_UI - target_UI_inflow)/target_UI_inflow)**2 #
    obj_value += ((flow_out_UI - target_UI_outflow)/target_UI_outflow)**2
    obj_value += ((flow_reprise_avg - target_reprise_inflow)/target_reprise_inflow)**2
    obj_value += ((share_68_odd - target_share68_inflow_odd)/target_share68_inflow_odd)**2
    
    # moments alternatifs/supplémentaires qu'il serait possible d'inclure dans une autre version: nombre d'épuisements de droits (exclut car fortement affecté par le fait que la montée en charge de la réforme de 2023 est encore en cours dans les données), nombre d'individus en emploi, flux total d'entrants (au lieu de plusieurs flux séparés), flux de rechargements (exclut car difficile à définir exactement dans les données agrégées)
    
    # Affichage du diagnostic d'optimisation à chaque évaluation
    print("--------------------------------------------------")
    print("Valeur de la fonction objectif :", obj_value)
    print(f"taux_retour_emploi={job_finding_rate}, taux_destruction_emploi={job_destruction_rate}, "
          f"taux_recours={take_up_rate}, coef_jdr={jdr_coef}") 
    print(f"stock ARE={ui_unemployed_avg} vs {target_UI_stock}, flux entrée ARE={flow_into_UI} vs {target_UI_inflow}, emploi={employed_avg} vs {nbr_individuals_emp}, "
          f"flux sortie ARE={flow_out_UI} vs {target_UI_outflow}, part flux entrants 6-8 mois={share_68_odd} vs {target_share68_inflow_odd}")
    print(f"stock ARE={ui_unemployed_avg}, stock emploi={employed_avg}, stock chômage non indemnisé={noui_unemployed_avg}")
    print(f"flux rechargement = {flow_recharg_avg} vs {target_recharg_inflow}, flux reprise = {flow_reprise_avg} vs {target_reprise_inflow}, épuisement vers non-indemnisé ={epuise_uni} vs {target_UI_epuise} , flux entrée total = {flow_ui_full} vs {flow_ui_full_target}, part durée 15-18 mois = {share_pbd_15_18_odd}", )
    
    print(f"durée moyenne des droits ouverts lors des rechargements = {duree_droits_recharg_avg}")
    print(f"durée moyenne des droits ouverts lors des odd initiales = {duree_droits_init_avg}")
    
    return obj_value
