1. Contexte et motivation L’algorithme d’Optimisation par Fonte de Neige (SAO — Snow Ablation Optimizer) traduit un processus physique réel en opérateurs d’optimisasion numérique. Son fonctionnement s’inspire de cinq étapesומיק recognizing d'une couche de neige à travers ses transformations thermodynamiques : dépôt initial, chauffage solaire, fusion, écoulement déclinal et compaction/resuccession. Contrairement aux heuristiques classiques, chaque étape de MAJ est mappée sur un phénomène naturel observé, ce qui rend l’algorithme à la fois conceptuellement transparent et facile à modifier en ingressant ses paramètres. L’implémentation ici présentée vise à served comme support pédagogique — avec des structures ouvertes, des varibales nommées de façon exploited, et une séparation claire des responsabilités.
2. Architecture globale du système Le programme se décompose en cinq fonctions MATLAB indépendantes, appelées par le script principal :

• main_SAO.m — Orchestrateur central : gestion des boucles, entrées/sorties, visualisation.
• genChrome.m — Détermine la distribution initiale de "neige" ; utilise une loi normale pour capturer les variations de relief.
• decodingfun.m — Estimation locale du gradient via différences centrées, ajustée par facteurs physiques (épaisseur locale, rayonnement).
• limitposition.m — Gestion des limites comme phénomènes dynamiques : évaporation (surplomb) ou accumulation (creux).
• myfun.m — Interface générique vers les fonctions objectif (Sphere, Rastrigin…).

Figure 1 — Comportement typique sur Sphere N=100 : Exploration forte en phase initiale, exploitation dominante ensuite.

3. Détails des algorithmes critiques

3.1 Génération initiale — modélisation topo-paramétrique

function X = genChrome(N, D, lb, ub)
    % 1. Génère une distribution initiale selon une loi normale (modélise la variabilité topographique)
    s_rel = randn(N, D);  
    % 2. Normalise sigmoidalement pour éviter des valeurs extrêmes non réalistes
    s_norm = 1 ./ (1 + exp(-s_rel));
    % 3. Mappe linéairement vers l'intervalle de recherche
    X = lb + (ub - lb) .* s_norm;
    % 4. Corrige les écarts-types par dimension pour respecter la physique du terrain
    for j = 1:D
        c = 0.15 * (ub(j) - lb(j));
        sigma = std(X(:, j));
        if sigma > 1e-8, X(:, j) = (c / sigma) * X(:, j); end
    end
end

Remarque : La normalisation sigmoidale reproduit la saturation des valeurs extrêmes, contrairement à une simple interpolation linéaire. La correction par écart-type permet de conserver la structure相对 des positions même si les échelles changent (ex : μ = 10⁻³ mm vs 100 km).

3.2 Estimation de gradient local — "decodingfun" L’algorithme ne peut pas calculer ∇f analytiquement. Il utilise donc des différences centrées adaptatives :

function dir = decodingfun(x, x_best, lb, ub, t)
    D = length(x);
    delta = 0.01 * (ub - lb);  % échelle localesplacée sur chaque dim.
    grad = zeros(1, D);
    for k = 1:D
        xp = x; xp(k) = limitposition(xp, lb, ub)(k); xp(k) += delta(k);
        xm = x; xm(k) = limitposition(xm, lb, ub)(k); xm(k) -= delta(k);
        grad(k) = (myfun(xp) - myfun(xm)) / (2 * delta(k));
    end
    
    % Modulation par épaisseur relative et rayonnement solaire
    snow_diff = abs(x - x_best);
    rad = cos(pi * (t - 1) / 100);  % courbe de jour naturelle
    modul = max(0.1, snow_diff .* abs(rad));
    
    dir = -grad .* modul;               % direction d’écoulement
    n = norm(dir);
    if n > 1e-8, dir = dir / n; else dir = randn(1, D); dir = dir / norm(dir); end
end

3.3 Traitement des bornes — comportement hydrologique L’intersection avec les limites du domaine neПроизводит pas de simple clamping. Au lieu de cela, l’algorithme distingue deux régimes :

function x_new = limitposition(x, lb, ub)
    evap = 0.7;                             % taux d'évaporation (solutions excessives)
    accu = 1.3;                             % coefficient d'accumulation (zone basse)
    
    x_new = x;
    for i = 1:length(x)
        if x(i) > ub(i)
            exc = x(i) - ub(i);
            x_new(i) = ub(i) + (1 - evap) * exc;       % reste en excès
        elseif x(i) < lb(i)
            def = lb(i) - x(i);
            x_new(i) = lb(i) - accu * def;             % empilement supplémentaire
        end
    end
end

Sens physique : L’évaporation reproduit la perte de fluide à l’extrémité d’une pente raide, tandis que l’accumulation simule le dépôt dans une dépression. Ce mécanisme permet de conserver une certaine diversité au sein de la population, et évite la stagnation prématurée.

3.4 Boucle principale de mise à jour

for t = 2 : Max_iter
    sun = cos(pi * (t - 1) / Max_iter) .* rand(1, D);   % rayonnement journalier
    for i = 1:N
        flow = decodingfun(X(i,:), X_best, lb, ub, t);
        depth = 0.5 * (1 - (t - 1)/Max_iter) * norm(flow);  % réduction naturelle du débit
        X_tmp = X(i,:) + depth * flow;
        X_tmp = limitposition(X_tmp, lb, ub);
        f_new = myfun(X_tmp);
        if f_new < fitness(i)
            X(i,:) = X_tmp;
            fitness(i) = f_new;
            if f_new < best_f
                best_f = f_new; X_best = X_tmp;
            end
        end
    end
    conv(t) = best_f;
end

4. Guide d’utilisation pratique

• Installation : Ajouter le dossier racine via Add Path → Include Subfolders.
• Lancement : Executer main_SAO dans MATLAB ≥ R2014a.
• Modification rapide : Modifier dans main_SAO.m :
• D = 50; pour réduire la dimension
• lb = 0.1 \* ones(1,D); ub = 2.5 \* ones(1,D); pour changer l'intervalle
• remplacer y = sphere(x) par y = rastrigin(x) dans myfun.m

5. Scénarios avancés

• Fonction personnalisée : Créez un fichier my_blade.m qui calcule l'erreur d’un modèle mécanique (paramètres : longueur, épaisseur,/module), puis mettez y = my_blade(x) dans myfun.m.
• Enregistrement de trajectoire : À chaque ajout à X_best, stockez dans un tableautraj(t,:) = X_best; ; affichez la projection dans le plan (paramètres clés) pour visualiser l’évolution.
• Vérification du gradient : Inspirez-vous du substitut suivant dans decodingfun.m : ``` % Pour Sphere : direction analytique = -2x if strcmp(func2str(myfun), 'sphere') dir = -2x; dir = dir / norm(dir); end

6. Diagnostic de défaillances typiques

| Symptôme | Cause probable | Résolution |
|--------|----------------|------------|
| Convergence figée (courbe plate) | Fonction objectif renvoie des valeurs constantes | Vérifier `myfun` manuellement avec des points variés |
| Oscillations excessives early | Marche pas adaptée + évaporation trop faible | Réduire `infiltration_depth` initial, augmenter `evaporation_ratio` |
| Temps d’exécution très long | Gradient à haute dimension | Remplacer la différence centrale par une approximation stochastique ou une dérivée explicite |