بعد تنفيذ كود ماتلاب التالي  بالتأكيد! سأقوم بإعادة كتابة الكود لاستخدام خوارزمية تحسين أسراب الجزيئات (Particle Swarm Optimization - PSO) بدلاً من الخوارزمية الوراثية (GA) لحل مشكلة توجيه المركبات (Vehicle Routing Problem - VRP). سأضيف أيضًا تعليقات تفصيلية باللغة العربية لشرح كل جزء من الكود.

مقدمة:

خوارزمية PSO هي تقنية تحسين تعتمد على محاكاة سلوك أسراب الطيور أو الأسماك في البحث عن الطعام. في سياق مشكلة توجيه المركبات، سنقوم بتطبيق PSO للعثور على أفضل مسارات المركبات لتوصيل العملاء مع مراعاة قيود معينة مثل السعة والوقت.


% مسح نافذة الأوامر
clc;

% تنظيف جميع المتغيرات من مساحة العمل
clear;

%% **الجزء الأول: إعداد البيانات والمعلومات الأساسية**

% تحديد حجم السرب (عدد الجزيئات)
swarm_size = 100;

% تحديد عدد الأجيال (التكرارات) لخوارزمية PSO
max_iterations = 100;

% قراءة الإحداثيات من ملف CSV
coordinates = readcell('cordinates.csv');

% تحويل الإحداثيات إلى الصيغة العشرية باستخدام الدالة الجديدة
[latitudes, longitudes] = convert_coordinates(coordinates(:,1), coordinates(:,2));

% قراءة بيانات العملاء من ملف CSV
tbl_customers = readcell('c101customers35.csv');

% استبدال الأعمدة الثاني والثالث في tbl_customers بالإحداثيات المحولة
tbl_customers(:, 2) = num2cell(latitudes);
tbl_customers(:, 3) = num2cell(longitudes);

% حفظ البيانات المحدثة في نفس الملف
writecell(tbl_customers, 'c101customers35.csv');

% قراءة بيانات العملاء من CSV إلى مصفوفة
tbl_customers = readmatrix('c101customers35.csv');

% استخراج الوقت الكلي من مصفوفة البيانات
total_time = tbl_customers(1, 6);

% استخراج إحداثيات العملاء والمستودعات (إحداثيات خطوط العرض والطول)
latitudes = tbl_customers(:, 2);
longitudes = tbl_customers(:, 3);

% تخزين الإحداثيات في مصفوفة
vertexs = [latitudes, longitudes];

% حساب المسافات بين كل زوج من النقاط باستخدام دالة هافرساين
n = size(vertexs, 1); % عدد النقاط (العملاء والمستودعات)
arr_distance_matrix = zeros(n); % مصفوفة لتخزين المسافات

% حلقة مزدوجة لحساب المسافة بين كل زوج من النقاط
for i = 1:n
    for j = 1:n
        arr_distance_matrix(i,j) = haversine_distance(vertexs(i,1), vertexs(i,2), vertexs(j,1), vertexs(j,2));
    end
end

% تحديد الأعمدة التي سيتم حذفها (الأعمدة 2 و3)
columns_to_remove = [2, 3];

% التحقق من وجود العمود الثامن وحذفه إذا كان موجودًا
if width(tbl_customers) >= 8
    columns_to_remove = [columns_to_remove, 8];
end

% التحقق من وجود العمود التاسع وحذفه إذا كان موجودًا
if width(tbl_customers) >= 9
    columns_to_remove = [columns_to_remove, 9];
end

% حذف الأعمدة من الجدول
tbl_customers(:, columns_to_remove) = [];

% تخزين الأعمدة المتبقية في مصفوفة
arr_customers = tbl_customers;

% الحصول على عدد العملاء بناءً على عدد الصفوف في arr_customers
num_customers = size(arr_customers, 1);

% تحديد السعة الكلية للمركبات
total_capacity = 200;

% إنشاء قائمة تضم جميع العملاء
all_customers_list = 1:num_customers;

% إنشاء متغير لخزن أفضل حل
global_best_solution = [];
global_best_cost = inf;

% إنشاء مصفوفة لخزن الحلول الأفضل لكل جزيء
particle_best_solutions = cell(swarm_size, 1);
particle_best_costs = inf(swarm_size, 1);

% مبادئ التحكم في PSO
w = 0.729;          % وزن القصور الذاتي
c1 = 1.49445;       % معامل التعلم الشخصي
c2 = 1.49445;       % معامل التعلم الاجتماعي

% تهيئة السرب (الجزيئات)
particles = cell(swarm_size, 1);          % مسارات الجزيئات
particles_velocities = cell(swarm_size, 1); % سرعات الجزيئات

% تهيئة الجزيئات بحلول عشوائية
for i = 1:swarm_size
    % إنشاء حل عشوائي (مسار مركبة)
    particles{i} = generate_initial_solution(arr_customers, total_capacity, total_time, arr_distance_matrix);
    
    % تهيئة أفضل حل للجزيء
    particle_best_solutions{i} = particles{i};
    
    % حساب تكلفة الحل
    particle_best_costs(i) = calculate_total_cost(particles{i}, arr_distance_matrix, total_capacity, total_time, arr_customers);
    
    % تحديث الحل الأفضل عالميًا
    if particle_best_costs(i) < global_best_cost
        global_best_cost = particle_best_costs(i);
        global_best_solution = particles{i};
    end
    
    % تهيئة سرعة الجزيء (يمكن أن تكون صفرية أو قيمة عشوائية)
    particles_velocities{i} = initialize_particle_velocity(particles{i});
end

%% **الجزء الثاني: الحلقة الرئيسية لخوارزمية PSO**

for iteration = 1:max_iterations
    for i = 1:swarm_size
        % تحديث سرعة الجزيء
        particles_velocities{i} = update_velocity(particles_velocities{i}, particles{i}, particle_best_solutions{i}, global_best_solution, w, c1, c2);
        
        % تحديث موقع الجزيء (الحل)
        particles{i} = update_position(particles{i}, particles_velocities{i});
        
        % ضمان صلاحية الحل بعد التحديث
        particles{i} = repair_solution(particles{i}, arr_customers, total_capacity, total_time, arr_distance_matrix);
        
        % حساب تكلفة الحل الجديد
        current_cost = calculate_total_cost(particles{i}, arr_distance_matrix, total_capacity, total_time, arr_customers);
        
        % تحديث أفضل حل للجزيء إذا كان الحل الحالي أفضل
        if current_cost < particle_best_costs(i)
            particle_best_costs(i) = current_cost;
            particle_best_solutions{i} = particles{i};
        end
        
        % تحديث الحل الأفضل عالميًا
        if particle_best_costs(i) < global_best_cost
            global_best_cost = particle_best_costs(i);
            global_best_solution = particle_best_solutions{i};
        end
    end
    
    % طباعة معلومات التكرار الحالي
    fprintf('Iteration %d: Best Cost = %f\n', iteration, global_best_cost);
end

%% **الجزء الثالث: عرض الحل الأفضل**

% استخراج المسارات من الحل الأفضل
best_routes = global_best_solution;

% حساب المسافات لكل مسار
best_distances = zeros(length(best_routes), 1);
for i = 1:length(best_routes)
    best_distances(i) = calculate_route_distance(best_routes{i}, arr_distance_matrix);
end

% رسم الخريطة
drawMap(best_routes, vertexs, arr_distance_matrix);

% عرض معلومات المسارات
for i = 1:length(best_routes)
    fprintf('Route %d: Distance = %f\n', i, best_distances(i));
    fprintf('Route %d Sequence: %s\n', i, mat2str(best_routes{i}));
end

%% **الدوال المستخدمة في الكود**

%% دالة توليد حل ابتدائي عشوائي
function initial_solution = generate_initial_solution(arr_customers, total_capacity, total_time, arr_distance_matrix)
    num_customers = size(arr_customers, 1);
    customers = 2:num_customers; % العملاء بدون المستودع (نفترض أن المستودع هو 1)
    % خلط العملاء عشوائيًا
    customers = customers(randperm(length(customers)));
    routes = {};
    current_route = [1]; % بدء المسار من المستودع
    current_capacity = 0;
    current_time = 0;
    for customer = customers
        demand = arr_customers(customer, 2);
        ready_time = arr_customers(customer, 3);
        due_time = arr_customers(customer, 4);
        service_time = arr_customers(customer, 5);
        travel_time = arr_distance_matrix(current_route(end), customer);
        arrival_time = current_time + travel_time;
        wait_time = max(0, ready_time - arrival_time);
        finish_time = arrival_time + wait_time + service_time;
        % التحقق من القيود
        if current_capacity + demand <= total_capacity && finish_time <= due_time && finish_time + arr_distance_matrix(customer, 1) <= total_time
            % إضافة العميل إلى المسار الحالي
            current_route = [current_route, customer];
            current_capacity = current_capacity + demand;
            current_time = finish_time;
        else
            % إنهاء المسار الحالي والبدء بمسار جديد
            current_route = [current_route, 1];
            routes{end+1} = current_route;
            current_route = [1, customer];
            current_capacity = demand;
            current_time = travel_time + service_time;
        end
    end
    % إنهاء المسار الأخير
    if current_route(end) ~= 1
        current_route = [current_route, 1];
        routes{end+1} = current_route;
    end
    initial_solution = routes;
end

%% دالة حساب تكلفة الحل (يمكن تعديلها لتشمل عوامل أخرى)
function total_cost = calculate_total_cost(solution, arr_distance_matrix, total_capacity, total_time, arr_customers)
    total_distance = 0;
    for i = 1:length(solution)
        route = solution{i};
        if ~validate_route(route, arr_customers, total_capacity, total_time, arr_distance_matrix)
            total_cost = inf;
            return;
        end
        total_distance = total_distance + calculate_route_distance(route, arr_distance_matrix);
    end
    % يمكن تعديل دالة التكلفة لتشمل عدد المركبات (المسارات)
    total_cost = total_distance; % + alpha * number_of_vehicles (إذا رغبت في ذلك)
end

%% دالة حساب مسافة المسار
function route_distance = calculate_route_distance(route, arr_distance_matrix)
    route_distance = 0;
    for i = 1:length(route)-1
        route_distance = route_distance + arr_distance_matrix(route(i), route(i+1));
    end
end

%% دالة التحقق من صلاحية المسار
function is_valid = validate_route(route, arr_customers, total_capacity, total_time, arr_distance_matrix)
    current_capacity = 0;
    current_time = 0;
    is_valid = true;
    for i = 2:length(route)-1
        customer = route(i);
        demand = arr_customers(customer, 2);
        ready_time = arr_customers(customer, 3);
        due_time = arr_customers(customer, 4);
        service_time = arr_customers(customer, 5);
        travel_time = arr_distance_matrix(route(i-1), customer);
        arrival_time = current_time + travel_time;
        wait_time = max(0, ready_time - arrival_time);
        finish_time = arrival_time + wait_time + service_time;
        if current_capacity + demand > total_capacity || finish_time > due_time || finish_time + arr_distance_matrix(customer, 1) > total_time
            is_valid = false;
            return;
        end
        current_capacity = current_capacity + demand;
        current_time = finish_time;
    end
    % التحقق من العودة إلى المستودع في الوقت المحدد
    if current_time + arr_distance_matrix(route(end-1), 1) > total_time
        is_valid = false;
    end
end

%% دالة تهيئة سرعة الجزيء
function velocity = initialize_particle_velocity(solution)
    % في VRP، السرعة يمكن أن تمثل تبادلات أو تحولات في المسارات
    % هنا، يمكننا تهيئة السرعة كقائمة من العمليات (مثل التبادلات العشوائية)
    velocity = {};
end

%% دالة تحديث سرعة الجزيء
function new_velocity = update_velocity(current_velocity, current_solution, personal_best_solution, global_best_solution, w, c1, c2)
    % تحديث السرعة بناءً على المعادلة القياسية لـ PSO
    % نظرًا لأن الحلول هنا هي مسارات، سنحتاج إلى تمثيل السرعة بعمليات على المسارات (مثل التبادلات)
    % لتبسيط الأمر، سنستخدم التبادلات العشوائية
    new_velocity = current_velocity;
    % إضافة تأثير التعلم الشخصي
    new_velocity = [new_velocity, generate_exchange_operations(current_solution, personal_best_solution, c1)];
    % إضافة تأثير التعلم الاجتماعي
    new_velocity = [new_velocity, generate_exchange_operations(current_solution, global_best_solution, c2)];
    % تطبيق وزن القصور الذاتي
    % في هذا السياق، يمكن أن يكون ذلك بتقليل عدد العمليات في السرعة
    new_velocity = new_velocity(1:round(w * length(new_velocity)));
end

%% دالة تحديث موقع الجزيء (الحل)
function new_solution = update_position(current_solution, velocity)
    new_solution = current_solution;
    % تطبيق العمليات المحددة في السرعة على الحل الحالي
    for i = 1:length(velocity)
        operation = velocity{i};
        new_solution = apply_exchange_operation(new_solution, operation);
    end
end

%% دالة لإصلاح الحل في حالة وجود تعارضات
function repaired_solution = repair_solution(solution, arr_customers, total_capacity, total_time, arr_distance_matrix)
    % يمكننا التحقق من كل مسار وإجراء تعديلات إذا لزم الأمر
    repaired_solution = solution;
    % في هذا المثال، لن نقوم بأي شيء محدد، ولكن يمكن تنفيذ خوارزميات إصلاح أكثر تقدمًا
end

%% دالة لتوليد عمليات التبادل بين الحلين
function operations = generate_exchange_operations(current_solution, target_solution, coefficient)
    operations = {};
    % نظرًا لأن الحلول هي مسارات، يمكننا تحديد العملاء المختلفين بين الحلين
    % ثم إنشاء عمليات تبادل لإحضار الحل الحالي نحو الحل المستهدف
    % هنا، سنختار عددًا من العمليات بناءً على معامل التعلم
    num_operations = round(coefficient * 5); % مثال على استخدام المعامل
    for i = 1:num_operations
        % اختيار عشوائي لمسار وموضع في الحلين
        route_idx = randi([1, min(length(current_solution), length(target_solution))]);
        route_current = current_solution{route_idx};
        route_target = target_solution{route_idx};
        % اختيار عميل عشوائي لتبديله
        if length(route_current) > 2 && length(route_target) > 2
            pos_current = randi([2, length(route_current)-1]);
            pos_target = randi([2, length(route_target)-1]);
            customer_current = route_current(pos_current);
            customer_target = route_target(pos_target);
            % إنشاء عملية التبادل
            operations{end+1} = [route_idx, customer_current, customer_target];
        end
    end
end

%% دالة لتطبيق عملية التبادل على الحل
function new_solution = apply_exchange_operation(solution, operation)
    new_solution = solution;
    route_idx = operation(1);
    customer_current = operation(2);
    customer_target = operation(3);
    route = new_solution{route_idx};
    % إيجاد مواقع العملاء في المسار
    idx_current = find(route == customer_current, 1);
    idx_target = find(route == customer_target, 1);
    if ~isempty(idx_current) && ~isempty(idx_target)
        % تبديل العملاء
        route([idx_current, idx_target]) = route([idx_target, idx_current]);
        new_solution{route_idx} = route;
    end
end

%% دالة لتحويل الإحداثيات من الدرجات والدقائق والثواني إلى صيغة عشرية
function [lat_dec, lon_dec] = convert_coordinates(lat, lon)
    n = length(lat);
    lat_dec = zeros(n, 1);
    lon_dec = zeros(n, 1);
    for i = 1:n
        lat_parts = sscanf(lat{i}, '%d° %d''%d"N');
        lon_parts = sscanf(lon{i}, '%d° %d''%d"E');
        lat_dec(i) = lat_parts(1) + lat_parts(2)/60 + lat_parts(3)/3600;
        lon_dec(i) = lon_parts(1) + lon_parts(2)/60 + lon_parts(3)/3600;
    end
end

%% دالة لحساب المسافة باستخدام صيغة هافرساين
function d = haversine_distance(lat1, lon1, lat2, lon2)
    R = 6371; % نصف قطر الأرض بالكيلومترات
    % تحويل الإحداثيات من درجات إلى راديان
    lat1 = deg2rad(lat1);
    lon1 = deg2rad(lon1);
    lat2 = deg2rad(lat2);
    lon2 = deg2rad(lon2);
    % تطبيق صيغة هافرساين
    dlat = lat2 - lat1;
    dlon = lon2 - lon1;
    a = sin(dlat/2)^2 + cos(lat1) * cos(lat2) * sin(dlon/2)^2;
    c = 2 * atan2(sqrt(a), sqrt(1-a));
    d = R * c;
end

%% دالة لرسم الخريطة
function drawMap(best_routes, vertexs, arr_distance_matrix)
    clf;
    customer = vertexs(1:end,:);
    num_vehicules = size(best_routes, 1);
    hold on;
    box on;
    title('Optimal Roadmap')
    hold on;
    C = linspecer(num_vehicules); % تعيين لون فريد لكل مسار
    % تخزين الألوان لكل مسار
    route_colors = cell(num_vehicules, 1);
    for i = 1:num_vehicules
        route_colors{i} = C(i, :);
    end
    for i = 1:num_vehicules
        part_seq = best_routes{i};
        len = length(part_seq);
        color = route_colors{i}; % الحصول على لون المسار الحالي
        for j = 1:len-1
            c1 = customer(part_seq(j), :);
            c2 = customer(part_seq(j+1), :);
            plot([c1(1), c2(1)], [c1(2), c2(2)], '-', 'color', color, 'linewidth', 1);
            % عرض المسافة
            mid_x = (c1(1) + c2(1)) / 2;
            mid_y = (c1(2) + c2(2)) / 2;
            text(mid_x, mid_y, num2str(arr_distance_matrix(part_seq(j), part_seq(j + 1))), 'Color', color, 'FontSize', 8, 'FontWeight', 'bold');
        end
    end
    % رسم نقاط العملاء
    plot(customer(:, 1), customer(:, 2), 'ro', 'linewidth', 1);
    % إضافة أرقام العملاء
    for i = 2:size(customer, 1)  % بدءًا من 2 لتخطي المستودع
        text(customer(i, 1), customer(i, 2), num2str(i-1), 'HorizontalAlignment', 'center', 'VerticalAlignment', 'bottom', 'Color', 'k', 'FontWeight', 'bold');
    end
    hold on;
    % رسم المستودع
    plot(vertexs(1, 1), vertexs(1, 2), 's', 'linewidth', 2, 'MarkerEdgeColor', 'b', 'MarkerFaceColor', 'b', 'MarkerSize', 10);
    % إضافة تسمية المستودع
    text(vertexs(1, 1), vertexs(1, 2), 'Depot', 'HorizontalAlignment', 'center', 'VerticalAlignment', 'top', 'Color', 'b', 'FontWeight', 'bold');

    % تمكين وظيفة التكبير
    zoom on;
end
شرح تفصيلي للكود:

إعداد البيانات:
نقوم بقراءة بيانات الإحداثيات والعملاء من ملفات CSV وتحويل الإحداثيات إلى صيغة عشرية.
نحسب مصفوفة المسافات بين جميع النقاط باستخدام دالة هافرساين.
تهيئة خوارزمية PSO:
نحدد حجم السرب (عدد الجزيئات) وعدد التكرارات.
نهيئ الجزيئات بحلول عشوائية ونحسب تكلفة كل حل.
نحتفظ بأفضل حل لكل جزيء (البعد الشخصي) وأفضل حل عالمي (البعد الاجتماعي).
الحلقة الرئيسية لـ PSO:
في كل تكرار، نقوم بتحديث سرعة وموقع كل جزيء بناءً على المعادلات القياسية لـ PSO.
نقوم بتطبيق عمليات الإصلاح على الحلول لضمان تلبية القيود.
نحدث الحلول الأفضل الشخصية والعالمية.
عرض النتائج:
بعد انتهاء التكرارات، نعرض الحل الأفضل من خلال طباعة المسارات والمسافات ورسم الخريطة.
الدوال المساعدة:
تم تعريف عدة دوال مساعدة لتنفيذ مختلف العمليات، مثل:
توليد الحلول العشوائية.
حساب تكلفة الحل.
تحديث السرعة والموقع للجزيئات.
التحقق من صلاحية المسارات.
رسم الخريطة.
ملاحظات:

في خوارزمية PSO، تم تمثيل السرعة كعمليات تبادل بين العملاء في المسارات.
تم استخدام معاملات قياسية لـ PSO (w، c1، c2) ولكن يمكن تعديلها لتحسين الأداء.
يمكن تحسين دالة الإصلاح لضمان جودة الحلول بشكل أفضل.
تم تبسيط بعض الأجزاء للمحافظة على وضوح الكود، ولكن يمكن تطويرها لتحسين الأداء والدقة.
الخلاصة:

قمنا بإعادة كتابة الكود لاستخدام خوارزمية PSO بدلاً من GA لحل مشكلة توجيه المركبات. الكود يتضمن تهيئة السرب، تحديث السرعات والمواقع، الحفاظ على أفضل الحلول، وأخيرًا عرض النتائج.

ظهر الخطا التالي   Error using randi
First input must be a positive scalar integer value IMAX, or two integer values [IMIN IMAX] with IMIN less than or equal to
IMAX.

Error in Full_PSO>generate_exchange_operations (line 317)
        route_idx = randi([1, min(length(current_solution), length(target_solution))]);

Error in Full_PSO>update_velocity (line 285)
    new_velocity = [new_velocity, generate_exchange_operations(current_solution, global_best_solution, c2)];

Error in Full_PSO (line 126)
        particles_velocities{i} = update_velocity(particles_velocities{i}, particles{i}, particle_best_solutions{i}, global_best_solution, w, c1, c2);
 
>>