الگوریتم جدیدی برای حل مسأله مسیریابی-موجودی با ارسال مستقیم

نوع مقاله: مقاله پژوهشی

نویسندگان

1 دانشجوی دکتری مهندسی صنایع دانشکده فنی و مهندسی دانشگاه تربیت مدرس

2 دانشیار دانشکده فنی و مهندسی دانشگاه تربیت مدرس

چکیده

این مقاله به بررسی مسأله مسیریابی-موجودی چند محصولی چند دوره­ای در یک زنجیره تأمین دو سطحی؛ شامل یک تولیدکننده و مجموعه­ای از خرده­فروشان اختصاص دارد. در مسأله مورد بررسی، علاوه بر مدیریت موجودی و برنامه­ریزی توزیع، برنامه­ریزی تولید نیز در نظر گرفته شده است. مسأله با هدف کمینه­سازی مجموع هزینه­های سیستم شامل هزینه­های راه­اندازی، توزیع و نگهداری موجودی مدلسازی شده است. محصولات توسط ناوگانی از وسایل حمل همسان با ظرفیت محدود تحت استرات‍ژی ارسال مستقیم به خرده­فروشان تحویل داده می­شوند. همچنین، ظرفیت تولید و نگهداری محدود و کمبود غیرمجاز فرض شده است. نشان داده شده است که مسایل مشابه بدون در نظر داشتن برنامه­ریزی توزیع در زمره مسایل با پیچیدگی سخت قرار دارند، بنابراین مسأله فوق نیز، مسأله­ای با پیچیدگی سخت است. از این رو، در این مقاله الگوریتم بهینه­سازی گروه ذرات بهبودیافته جدیدی برای حل آن توسعه داده شده است. الگوریتم پیشنهادی از دو بخش مجزا تشکیل شده است. نخست، مقادیر متغیرهای صفرویک با استفاده از الگوریتم پیشنهادی تعیین و سپس با حل یک مدل برنامه­ریزی خطی، مقادیر متغیرهای پیوسته محاسبه می­شود. کارایی الگوریتم پیشنهادی با استفاده از مسایل نمونه تصادفی متعددی با الگوریتم­های ژنتیک و بهینه­سازی گروه ذرات مقایسه شده است. نتایج محاسباتی بیانگر عملکرد بهتر الگوریتم پیشنهادی است.

کلیدواژه‌ها


عنوان مقاله [English]

A New Algorithm for Solving the Inventory Routing Problem with Direct Shipment

نویسندگان [English]

  • Ali Hossein Mirzaei 1
  • Isa Nakhai Kamalabadi 2
  • Seyed Hessameddin Zegordi 2
1 Ph.D. Student in Industrial Engineering, Faculty of Engineering, Tarbiat Modares University
2 Associate Professor, Faculty of Engineering, Tarbiat Modares University
چکیده [English]

In this paper a multi-commodity multi-period inventory routing problem in a two-echelon supply chain consisting of a manufacturer and a set of retailers has been studied. In addition to inventory management and distribution planning, production planning has also been considered in the above problem. The objective is to minimize total system cost that consists of production setup, inventory holding and distribution costs. The commodities are delivered to the retailers by an identical fleet of limited capacity vehicles through direct shipment strategy. Also it is assumed that production and storage capacity is limited and stockout is not allowed. Since similar problems without distribution planning are known as NP-hard, this is also an NP-hard problem. Therefore, in this paper, a new improved particle swarm optimization algorithm has been developed consisting of two distinguished phases for problem solving. First, the values of binary variables are determined using the proposed algorithm and then, the continuous variables are calculated by solving a linear programming model. Performance of the proposed algorithm has been compared with genetic and original particle swarm optimization algorithms using various samples of random problems. The findings imply significant performance of the proposed algorithm.

کلیدواژه‌ها [English]

  • Supply Chain
  • Inventory routing problem
  • Direct shipping strategy
  • Particle Swarm Optimization
  • Production-distribution planning

1- مقدمه

مسأله مسیریابی-موجودی[1] بسط مهمی از مسأله مسیریابی وسیله نقلیه[2] است که در آن تصمیمات کنترل موجودی و مسیریابی در هم ادغام می­شوند (کوردیو و همکاران[3]، 2007). مسأله مسیریابی-موجودی بیشتر در سیستم­های مدیریت موجودی توسط فروشنده[4] (VMI) کاربرد دارد. در سیستم­های مدیریت موجودی توسط فروشنده، فروشنده قادر است تا زمانبندی و اندازه تحویل محصول به خرده­فروشان را کنترل نماید. در قبال این آزادی عمل، فروشنده تضمین می­کند که مشتریان با کمبود مواجه نمی­شوند. در روابط سنتی­تر میان فروشنده و مشتری که در آن مشتریان درخواست سفارش محصولات را به فروشنده می­دادند، به دلیل زمانبندی سفارشات مشتریان، ممکن است کارایی به شدت کاهش و به نوبه آن هزینه­های موجودی و توزیع به شدت افزایش یابد. با وجود این، تحقق کاهش هزینه­های ناشی از به کارگیری سیستم­های VMI در عمل ساده نیست به ویژه با افزایش تعداد و تنوع مشتریان این امر دشوارتر نیز می­شود. با استفاده از مسأله مسیریابی-موجودی دستیابی به این هدف امکانپذیر است. در مسأله مسیریابی-موجودی با تعیین برنامه توزیع بهینه­ای که مجموع هزینه­ها را کمینه سازد، می­توان به این هدف دست یافت (کوردیو و همکاران، 2007). مسأله مسیریابی-موجودی در پژوهش­های متعددی بررسی شده است که مرور جامعی از پژوهش­های پیشین توسط اندرسون و همکاران[5] (2010) ارائه شده است. نویسندگان با بررسی ابعاد صنعتی مسأله، طبقه­بندی و مرور جامعی از پژوهش­های موجود ارائه داده­اند.

از جمله استراتژی­های توزیعی که می­توان در مسأله مسیریابی-موجودی به کار گرفت، استراتژی توزیع ارسال مستقیم[6] است. ارسال مستقیم به حالتی اشاره می­کند که در آن در هر دوره هر وسیله حمل تنها به یک خرده­فروش محصول تحویل می­دهد. به دلیل سادگی پیاده­سازی این استراتژی توزیع در عمل، طبیعی است که در مسأله مسیریابی-موجودی ابتدا چنین استراتژیی بررسی شود. از این رو، استراتژی ارسال مستقیم در مسأله مسیریابی-موجودی، بطور مکرر در سیستم­های توزیع صنعتی استفاده می­شود. کارایی استراتژی توزیع ارسال مستقیم حداقل به بزرگی جذر کوچکترین نرخ بهره­برداری از وسایل حمل است. به بیان دیگر، در صورتی که نرخ تقاضای هر خرده­فروش 90 درصد ظرفیت وسایل حمل ضرب در کران بالای فرکانس تحویل به خرده­فروشان باشد، آنگاه کارایی استراتژی ارسال مستقیم در حدود 86/94 درصد خواهد بود. کارایی ارسال مستقیم در شرایط خاص، حتی به 100 درصد نیز می­رسد. از این رو، در بسیاری از شرایط استفاده از استراتژی ارسال مستقیم توجیه منطقی دارد (لی و همکاران[7]، 2010). مقالات مختلفی با تمرکز بر استراتژی ارسال مستقیم، به بررسی مسأله مسیریابی-موجودی پرداخته­اند که از آن جمله می­توان به موارد زیر اشاره کرد. لی و همکاران (2010) به بررسی تحلیلی استراتژی ارسال مستقیم پرداخته­اند و کارایی ارسال مستقیم را در شرایط مختلف محاسبه کرده­اند. لی و همکاران (2008) مسأله مسیریابی-موجودی را در حالتی در نظر گرفته­اند که تأمین­کننده تنها یک وسیله حمل در اختیار داشته و در هر دوره تنها می­تواند برای یک مشتری موجودی ارسال کند. نویسندگان الگوریتمی ابتکاری برای یافتن توالی شدنی ارائه کرده­اند. کمبل و هاردین[8] (2005) با بررسی حداقل تعداد وسایل حمل مورد نیاز در مسأله مسیریابی-موجودی با ارسال مستقیم، برای حل مسأله، الگوریتمی حریصانه[9] ارائه داده­اند. چنگ و دوران[10] (2004) مدلی برای مسأله مسیریابی-موجودی در زنجیره تأمین جهانی نفت خام ارائه کرده­اند که در آن تقاضای مشتریان و طول زمان سفر، غیر قطعی است. همچنین، تقاضای مشتریان متغیر با زمان (پویا) است. در این مقاله یک سیستم پشتیبانی تصمیم­گیری برای طراحی و کنترل سیستم موجودی و حمل توسعه داده شده است. کلیوگت و همکاران[11] (2002) با بررسی مسأله مسیریابی-موجودی تصادفی با ارسال­های مستقیم و مدلسازی آن به صورت فرآیند تصمیم­گیری مارکوفی زمان گسسته، روشی تخمینی بر مبنای برنامه­ریزی پویا برای حل ارائه کرده­اند. بارنز شوستر و بسک[12] (1997) به ارزیابی کارایی استراتژی ارسال مستقیم در مسأله مسیریابی-موجودی با افق زمانی نامحدود پرداخته­اند که در آن تقاضای خرده­فروشان احتمالی با تابع توزیع معین و کمبود مجاز است.

مسأله مسیریابی وسایل نقلیه در زمره مسایل با پیچیدگی سخت[13] قرار دارد (لنسترا و رینوی[14]، 1981)؛ در نتیجه، مسأله مسیریابی-موجودی که توسعه­ای از مسأله مسیریابی وسایل نقلیه است، یک مسأله با پیچیدگی سخت است. بنابراین، در مقالات بسیاری، با استفاده از الگوریتم­های فراابتکاری، روش­های حل تقریبی برای مسأله توسعه داده شده است. عزیز و معین[15] (2007) مسأله را در حالت چند محصولی، چند دوره­ای با چند تأمین­کننده و یک کارخانه مونتاژ و با هدف کمینه­سازی مجموع هزینه‌های حمل و نگهداری موجودی بررسی کرده­اند. آنان با بررسی دو نحوه نمایش مختلف جواب، یک الگوریتم ژنتیک ترکیبی بر اساس رویکرد نخست تخصیص و سپس مسیریابی ارائه کرده­اند. معین و همکاران[16] (2010) در تحقیق مشابه­ای، الگوریتم ژنتیک ترکیبی بهبودیافته­ای ارائه داده­اند. بودیا و همکاران[17] (2007) در تحقیق خود با بررسی یک مسأله مسیریابی-موجودی که در آن برنامه­ریزی تولید نیز لحاظ شده است، الگوریتم­های حل مبتنی بر رویه جستجوی انطباقی حریصانه تصادفی[18] توسعه داده­اند. بودیا و پرینز[19] (2009) در تحقیق مشابه­ای، ساختاری مبتنی بر الگوریتم ممتیک[20] با مدیریت جمعیت جواب­ها توسعه داده­اند. ژائو و همکاران[21] (2008) مدل یکپارچه­ای برای مسأله مسیریابی-موجودی در زنجیره تأمین سه رده­ای ارائه کرده و برای حل آن الگوریتم جستجوی همسایگی بزرگ متغیری[22] توسعه داده­اند. ژائو و همکاران (2007) رویکرد حل جدیدی بر اساس الگوریتم فراابتکاری جستجوی ممنوع[23] برای حل مسأله مسیریابی-موجودی در یک زنجیره تأمین دو رده­ای توسعه داده­اند. اسپارچی-الکازار و همکاران[24] (2007) با پیاده­سازی الگوریتم ژنتیک برای حل مسأله مسیریابی-موجودی با چند محصول، تأثیر مقادیر مختلف پارامترهای ورودی الگوریتم ژنتیک را ارزیابی کرده­اند تا بهترین مجموع مقادیر پارامترهای الگوریتم بدست آید. راسدیانسیا و سائو[25] (2005) با توسعه مدلی برای مسأله مسیریابی-موجودی ماشین­های سکه­ای، به ارائه الگوریتمی دومرحله­ای مبتنی بر الگوریتم­های الحاق[26] و جستجوی ممنوع پرداخته­اند.

این مقاله به بررسی مسأله مسیریابی-موجودی در حالت چند دوره­ای چند محصولی و با هدف کمینه­سازی مجموع هزینه­های سیستم می­پردازد. هزینه­های سیستم شامل هزینه­های راه­اندازی، توزیع و نگهداری موجودی است. علاوه بر مدیریت موجودی و برنامه­ریزی توزیع، برنامه­ریزی تولید نیز در مسأله، ادغام شده است. در این مقاله، مسأله مسیریابی-موجودی در یک زنجیره تأمین دوبخشی متشکل از یک تولیدکننده و مجموعه­ای از خرده‌فروشان بررسی شده است که در آن، در هر دوره مقدار مشخصی از هریک از محصولات با استفاده از ناوگان همسانی از وسایل حمل با ظرفیت محدود و تحت استراتژی ارسال مستقیم، بین خرده­فروشان توزیع می­گردد. خدمت­رسانی به هریک از خرده­فروشان توسط وسایل حمل، یک دوره زمانی طول می­کشد. بنابراین، هر وسیله حمل در هر دوره تنها می­تواند به یک خرده­فروش محصول تحویل دهد. همچنین فرض شده است که محصول ارسالی به هریک از خرده­فروشان در هر دوره تنها باید طی یک ارسال به خرده­فروش تحویل داده شود. ظرفیت انباشت و تولید محدود است و کمبود مجاز نیست. مسأله مسیریابی-موجودی چند دوره­ای چند محصولی با ارسال مستقیم تحت مفروضات فوق، حالت توسعه یافته­ای از مسأله تعیین اندازه انباشته با ظرفیت محدود[27] در حالت چند محصولی، با ظرفیت تولید و هزینه­های نگهداری ثابت است. در پژوهش‌های پیشین نشان داده شده است که مسأله تعیین اندازه انباشته در شرایط فوق، یک مسأله با پیچیدگی سخت است (بیترن و یاناس[28]، 1981). بنابراین، می­توان نتیجه گرفت که مسأله بررسی شده در این مقاله نیز، با پیچیدگی سخت است. در نتیجه، حل مسأله مذکور با استفاده از روش­های حل دقیق[29] به ویژه برای مسایل با ابعاد بزرگ، در زمان محاسباتی معقول امکانپذیر نیست. در این صورت، می­توان با استفاده از الگوریتم­های فراابتکاری به جواب­های تقریبی مناسب در زمان محاسباتی معقول دست یافت. بنا به دانش نویسندگان، در پژوهش­های پیشین الگوریتم ابتکاری یا فراابتکاری برای حل مسأله مسیریابی-موجودی با مفروضات فوق، توسعه داده نشده است. از این رو، برای حل مسأله، الگوریتم بهینه­سازی گروه ذرات بهبودیافته[30] (IPSO) کارایی توسعه داده شده است. کارایی الگوریتم پیشنهادی با استفاده از مسایل نمونه تصادفی متعددی با الگوریتم­های ژنتیک و بهینه­سازی گروه ذرات[31] (PSO) مقایسه شده است.

ساختار ارائه مطالب بدین صورت خواهد بود که پس از مقدمه ارائه شده در این بخش، در بخش دوم، مدل ریاضی مسأله ارائه خواهد شد، در بخش سوم ساختار الگوریتم حل پیشنهادی تشریح خواهد گردید، بخش چهارم به بیان نتایج محاسباتی اختصاص یافته است و سرانجام در بخش پنجم نتیجه­گیری ارائه شده است.

 

2- مدل ریاضی

مسأله مسیریابی-موجودی چند دوره­­ای چند محصولی با ارسال مستقیم را می­توان در قالب مدل برنامه­ریزی عدد صحیح آمیخته[32] فرموله کرد. پیش از بیان مدل، ابتدا نمادهای استفاده شده برای توصیف مدل ارائه می­شود:

p : نشانگر محصول

j : نشانگر گره­ها (گره 0 نشانگر تولیدکننده و گره­های 1 تا N نشانگر خرده­فروشان)

t: نشانگر دوره‌های برنامه‌ریزی

N: تعداد خرده­فروشان

: هزینه ثابت راه­اندازی تولید برای تولیدکننده در دوره t.

: میزان موجودی محصول نوع p در پایان دوره t در مرکز j.

: هزینه حمل محصول از تولیدکننده به خرده­فروشی j .

: ضریب مصرف ظرفیت انباشت توسط محصول نوع p از فضای انبار یا وسیله حمل.

: ضریب مصرف ظرفیت تولید توسط محصول نوع p.

: ظرفیت انباشت در مرکز j .

: مقدار تقاضای محصول نوع p در خرده‌فروشی j در دوره t.

: ظرفیت تولیدکننده برای تولید محصول.

: ظرفیت وسیله حمل و نقل.

: هزینه نگهداری یک واحد محصول نوع p در مرکز j.

: میزان محصول نوع p ارسال شده از تولیدکننده به مرکز توزیع j در طی دوره t.

: میزان تولید محصول نوع p در طی دوره t توسط تولیدکننده.

: متغیر صفر و یک نشانگر ارسال محصول به خرده­فروش j در دوره t.

: متغیر صفر و یک نشانگر تولید محصول در دوره t.

مدل ریاضی مسأله، با استفاده از معادله‌های (1) الی (8) ارائه شده است:

 

 

(1) 

(2)  

(3)           

(4)                          

(5)               

(6)                        

(7)                               

(8)                        

معادله (1) بیانگر تابع هدف مدل پیشنهادی است که شامل هزینه­های ثابت راه­اندازی، هزینه ارسال محصول به خرده­فروشان، هزینه­های نگهداری محصول توسط تولیدکننده و خرده­فروشان است. معادله (2) بیانگر توزان موجودی بین تقاضا برای محصول نوع p در خرده‌فروشی j در دوره t و مجموع محصول نوع p ارسال شده از تولیدکننده به خرده‌فروشی l در دوره t است. توازن میان مجموع محصول نوع p ارسال شده از تولیدکننده به تمامی خرده‌فروشی‌ها در دوره t و مجموع محصول نوع p که در دوره t ام تولید شده است، با استفاده از معادله (3) بیان می­شود. معادله (4) محدودیت ظرفیت تولید را تضمین می­کند. محدودیت ظرفیت وسایل حمل و نقل با استفاده از معادله (5) ارائه شده است. معادله‌های (4) و (5) بطور ضمنی نشان می‌دهند که آیا در دوره t به ترتیب، توسط تولیدکننده تولیدی صورت گرفته است یا محصولی برای خرده­فروشی j ام ارسال شده است. معادله (6) تضمین­کننده محدودیت ظرفیت نگهداری برای تولیدکننده و خرده‌فروشی­ها است. محدودیت‌های فنی روی متغیرهای تصمیم از طریق معادله‌های (7) و (8) برقرار شده است.

 

3- الگوریتم بهینه­سازی گروه ذرات پیشنهادی

الگوریتم فراابتکاری بهینه­سازی گروه ذرات روش محاسباتی تکاملی مبتنی بر جمعیت جواب­ها است. مانند سایر الگوریتم­های فراابتکاری، الگوریتم مذکور ابزار بهینه­سازیی است که می­تواند برای حل انواع مختلفی از مسایل بهینه­سازی به­کار گرفته شود. این الگوریتم از جدیدترین روش­های فراابتکاری است که با الهام­گیری از رفتار اجتماعی گروهی از پرندگان مهاجر که در تلاش برای دستیابی به مقصد ناشناخته­ای هستند، توسط کندی و ابرهارت[33] (1995) در سال 1995 میلادی توسعه داده شده است. در الگوریتم PSO، جمعیت جواب­ها، گروه[34] نامیده می­شود و هر جواب مانند یک پرنده در گروهی از پرندگان است و ذره[35] نام دارد و شبیه کرموزوم در الگوریتم ژنتیک است. تمامی ذرات دارای مقدار شایستگی[36] هستند که با استفاده از تابع شایستگی[37] محاسبه می­گردند و تابع شایستگی ذرات باید بهینه گردد. جهت حرکت هر ذره توسط بردار سرعت[38] آن ذره معین می­شود. برخلاف الگوریتم ژنتیک، در فرآیند تکاملی الگوریتم مذکور، پرندگان جدیدی از نسل قبل (جواب­های جدید از جواب‌های قبلی) ایجاد نمی­گردد، بلکه هر پرنده رفتار اجتماعی خود را با توجه به تجربیاتش و رفتار سایر پرندگان گروه تکامل بخشیده و مطابق آن حرکت خود را به سوی مقصد بهبود می­دهد. به عبارت دیگر، در این الگوریتم عملگرهای تکاملی چون تقاطع[39] و جهش[40] وجود ندارد (هو و همکاران[41]، 2004؛ جائو و همکاران[42]، 2006). ساختار کلی الگوریتم بهینه­سازی گروه ذرات در شکل 1 ارائه شده است.

 

شکل 1. ساختار الگوریتم بهینه­سازی گروه ذرات

 

الگوریتم PSO با گروهی از جواب­های (ذرات) تصادفی آغاز و سپس با به هنگام­سازی ذرات در هر تکرار به دنبال جواب بهینه می­گردد (هو و همکاران، 2004). اگر متغیرهای تصمیم و به نوبه آن موقعیت ذرات، از نوع صفر و یک باشند؛ بردارهای سرعت و موقعیت هریک از ذرات در هر تکرار الگوریتم، طبق روابط (9) الی (12) محاسبه می­شوند (انگلبرکت[43]، 2005):

(9)(10)                              

(11)                                   

(12)                          

طبق رابطه (9)بردار سرعت جدید هر ذره بر اساس سرعت قبلی خود ذره (Vi(t-1))، بهترین موقعیتی که ذره تاکنون به آن دست یافته است (pBesti) و موقعیت بهترین ذره در همسایگی ذره که تابحال بدست آمده است (nBesti)، محاسبه می‌گردد. در صورتی که همسایگی هر ذره شامل تمام ذرات گروه باشد، آنگاه nBesti بیانگر موقعیت بهترین ذره در میان گروه است که با gBest به آن اشاره می­شود. r1 و r2 دو عدد تصادفی (با توزیع یکنواخت بین [0,1]) هستند که مستقل از یکدیگر تولید می­شوند. c1 و c2 که با نام ضرایب یادگیری به آنان اشاره شده است، تأثیر pBest و nBest را بر فرآیند جستجو کنترل می­نمایند. w بیانگر ضریب وزنی اینرسی است. بردار سرعت ذرات با مقدار Vmax محدود شده است. Vmax به عنوان محدودیتی است که قابلیت جستجوی جهانی گروه ذرات را کنترل می­کند. با استفاده از رابطه (11) بردار سرعت هریک از ذرات به بردار احتمال تغییر، تبدیل می­شود. در رابطه فوق، si بیانگر احتمال آن است که xit  برابر با 1 شود. سپس با استفاده از رابطه (12) بردار موقعیت هریک از ذرات بهنگام می­گردد. در رابطه فوق،عددی تصادفی با توزیع یکنواخت بین صفر و یک است.

مسأله مسیریابی-موجودی چند دوره­ای چند محصولی با ارسال مستقیم تحت مفروضات مورد بررسی در این مقاله، از جمله مسایل با درجه پیچیدگی سخت است. در نتیجه در این مقاله، الگوریتم بهینه­سازی گروه ذرات بهبودیافته جدیدی برای حل مسأله فوق توسعه داده شده است.

با توجه به معادلات (1) الی (8) مسأله مسیریابی-موجودی، یک مسأله برنامه­ریزی عدد صحیح آمیخته است که در آن متغیرهای تصمیم به دو دسته متغیرهای پیوسته (wpjt، Ipjt و Ppt) و متغیرهای صفر و یک (Zt و Xjt) قابل تفکیک هستند. از آنجایی که دشواری حل مسأله مذکور ناشی از متغیرهای صفر و یک است؛ بنابراین، با تفکیک روی متغیرهای تصمیم، می­توان راه­حلی دو بخشی برای حل مسأله توسعه داد. در راه­حل پیشنهادی ابتدا با استفاده از الگوریتم بهینه­سازی گروه ذرات بهبودیافته، مقادیر متغیرهای صفر و یک (Zt و Xjt) تعیین می­شود، سپس با توجه به مشخص بودن مقادیر متغیرهای صفر و یک، مسأله مسیریابی-موجودی، به یک مسأله برنامه­ریزی خطی برای تعیین مقادیر متغیرهای پیوسته (wpjt، Ipjt و Ppt) با هدف کمینه­سازی مجموع هزینه­ها تبدیل می­شود که به سادگی می­توان مسأله برنامه­ریزی خطی حاصل را با روش­های موجود برای حل مسایل برنامه­ریزی خطی حل کرد و سرانجام، به جواب نهایی مسأله اصلی دست یافت.

ساختار کلی الگوریتم پیشنهادی در شکل 2 ارائه شده است. همانگونه که در شکل نشان داده شده است، ساختار پیشنهادی مبتنی بر الگوریتم بهینه­سازی گروه ذرات، طراحی شده است. مطابق با شکل 2، الگوریتم پیشنهادی از دو مرحله تشکیل شده است که مرحله اول خود متشکل از دو بخش مجزاست. در بخش نخست، ابتدا به تعداد SSize ذره با استفاده از روشی شبه تصادفی به عنوان ذرات اولیه، ایجاد می­شوند. پس از آن، همسایگی تصادفی به اندازه NSize ذره، حول هر جواب اولیه، شکل می­گیرد. در طی گام­های الگوریتم، کیفیت و پراکندگی جواب­ها بهبود می­یابد و حداکثر به تعداد SSize ذره از جواب­های خوب و غیرتکراری انتخاب و در مجموعه­ای با نام مجموعه مرجع -برای استفاده در مراحل بعدی الگوریتم- قرار داده می­شود. در پایان بخش نخست، جواب­های موجود در مجموعه مرجع به عنوان ذرات اولیه برای شروع بخش دوم در نظر گرفته می­شوند. در صورتی که تعداد ذرات موجود در مجموعه مرجع کمتر از SSize ذره باشد، ذرات باقیمانده موردنیاز به صورت تصادفی ایجاد و افزوده می­شوند. در بخش دوم نیز فرآیند بهبود جواب­ها با مکانیزمی متفاوت ادامه می­یابد و در هر تکرار الگوریتم جواب­های خوب و غیرتکراری در مجموعه مرجع قرار داده می­شوند. در پایان بخش دوم، جواب­های موجود در مجموعه مرجع به عنوان ذرات اولیه برای شروع بخش نخست در نظر گرفته می­شوند. در این بخش نیز، در صورتی که تعداد ذرات موجود در مجموعه مرجع کمتر از SSize ذره باشد، ذرات باقیمانده موردنیاز به صورت تصادفی تولید و اضافه می­شوند. بخش اول و دوم، هر کدام به ترتیب پس از P1S1 و P1S2 مرتبه تکرار یا در صورتی که به ترتیب پس از طی CR11 و CR12 مرتبه تکرار در مقدار شایستگی بهترین جواب تغییری حاصل نشود، خاتمه می­یابند. همچنین، مرحله اول الگوریتم -رویه تکرار پیاپی بخش اول و دوم- پس از طی P1 مرتبه تکرار پایان می­پذیرد. در پایان مرحله نخست الگوریتم، مجموعه مرجع تشکیل شده در طی این مرحله، به عنوان جواب­های آغازین مرحله دوم الگوریتم، در نظر گرفته و الگوریتم وارد مرحله دوم می­شود. در مرحله دوم، ابتدا با استفاده از مدل برنامه­ریزی خطی، مقادیر بهینه متغیرهای تصمیم پیوسته با توجه به مقادیر متغیرهای صفر و یک (Zt و Xjt) تعیین و از بین جواب­ها، بهترین جواب انتخاب می­شود. پس از آن با استفاده از رویه جستجوی محلی کیفیت جواب منتخب بهبود می­یابد. مرحله دوم الگوریتم، بعد از P2 مرتبه تکرار یا در صورتی که پس از CR2 مرتبه تکرار تغییری در مقدار بهترین ذره حاصل نگردد، خاتمه می­یابد.

هرچند در تمام مراحل الگوریتم پیشنهادی می‌توان مقادیر متغیرهای تصمیم پیوسته (wpjt، Ipjt و Ppt) را با استفاده از حل یک مدل برنامه­ریزی خطی بدست آورد؛ با وجود این، به دلیل افزایش سرعت الگوریتم، در مرحله اول، روشی تقریبی برای تعیین مقادیر متغیرهای تصمیم پیوسته توسعه و استفاده می­شود.

در مرحله اول الگوریتم، بردارهای سرعت و موقعیت هریک از ذرات در هر تکرار الگوریتم، طبق روابط (9) الی (12) محاسبه می­شوند. در بخش دوم از مرحله نخست الگوریتم، در رابطه (9) مقدار nBesti با مقدار بهترین ذره در میان تمام ذرات (gBest) جایگزین می­شود.


 

شکل 2. ساختار الگوریتم بهینه­سازی گروه ذرات پیشنهادی

 

 

 

 

 

3-1- نحوه نمایش ذرات

یکی از مهمترین تصمیماتی که در زمان طراحی یک الگوریتم فراابتکاری می­بایستی اتخاذ شود، نحوه نمایش جواب­ها[44] و نیز چگونگی برقراری ارتباطی مؤثر، یکتا و قابل شناسایی میان جواب­ها و فضای جستجوی مسأله مورد بررسی است. در این مقاله برای نمایش ذرات، با توسعه روش عزیز و معین (2007)، از شیوه "نمایش زمان ارسال" برای نمایش بردار موقعیت ذرات استفاده شده است. در این شیوه نمایش از مقادیر صفر و یک برای نمایش زمان ارسال محصولات به خرده­فروشان مختلف و زمان تولید استفاده شده است (مقادیر متغیرهای Zt و Xjt). هر 1 نشاندهنده تولید (1=Zt) یا ارسال محصول به خرده­فروش در آن دوره (1=Xjt) و مقدار صفر بیانگر عدم تولید (0=Zt) یا ارسال محصول در آن دوره (0=Xjt) است. نمایی کلی از بردار موقعیت ذرات در نحوه نمایش زمان ارسال برای حل مسأله مسیریابی-موجودی با ارسال مستقیم برای حالتی با 6 دوره برنامه­ریزی و 3 خرده­فروش، در شکل 3 ارائه شده است. موجودی خرده­فروشان و نیز تولیدکننده در ابتدای دوره برنامه­ریزی صفر فرض شده است، بنابراین، در نخستین دوره برنامه­ریزی تمام مقادیر مربوط به محصولات و خرده­فروشان برابر با 1 است.

 

شکل 3. نمونه­ای از نحوه نمایش ذرات

 

3-2- تولید جواب­های اولیه

به طور کلی کیفیت جواب­های آغازین بر عملکرد الگوریتم­های فراابتکاری تأثیر بسزایی دارد. بنابراین، طراحی روشی مؤثر برای تولید جواب­های اولیه از اهمیت بالایی برخوردار است. از این رو، در الگوریتم پیشنهادی، برای تولید ذرات اولیه، روشی شبه تصادفی طراحی شده است که با به­کارگیری آن ذرات بین دو حالت کمترین و بیشترین تعداد دفعات ارسال ایجاد شوند. گام­های رویه پیشنهادی برای تولید ذرات آغازین در زیر تشریح شده است:

1-    تمام مؤلفه­های ذره نخست را برابر 1 قرار دهید. این ذره، ذره­ای با بیشترین تعداد دفعات ارسال و همواره شدنی است. به عبارت دیگر، ذره نخست بیانگر حالتی است که در آن، در تمام دوره­ها تولید و در تمام دوره­ها و برای تمام خرده­فروشان ارسال داشته باشیم.

2-    با استفاده از رویه ارائه شده در شکل 4، ذره دوم را ایجاد کنید. ذره دوم، ذره­ای با کمترین تعداد دفعات ارسال است. برای تولید ذره دوم، محدودیت انباشت که با استفاده از رابطه (6) بیان می­شود، آزاد شده است.

 

شکل 4. شبه کد تولید ذره با کمترین تعداد دفعات تولید و ارسال

3-    ذره اول و دوم را به عنوان والد اول و دوم قرار دهید. پس از آن، با استفاده از عملگر تقاطع پراکنده[45] سایر ذرات موردنیاز را ایجاد کنید. برای انجام عملگر تقاطع پراکنده، ابتدا یک آرایه تصادفی صفر و یک با اندازه برابر با اندازه ذرات ایجاد می­گردد. مؤلفه i ام ذره فرزند (ذره جدید) اگر مؤلفه i ام آرایه تصادفی 1 بود، برابر با مؤلفه i ام والد 1 (ذره اول) و در غیر اینصورت برابر با مؤلفه i ام والد 2 (ذره دوم) خواهد بود. نحوه انجام عملگر تقاطع پراکنده در شکل 5 ارائه شده است.

 

شکل 5. نحوه انجام عملگر تقاطع پراکنده

 

3-3. محاسبه مقادیر شایستگی

به طور معمول، در الگوریتم­های فراابتکاری هر جواب در هر مرحله با استفاده از تابع شایستگی ارزیابی و مقدار شایستگی آن محاسبه می­شود. در این مقاله نیز مانند سایر مسایل تک هدفه از تابع هدف مسأله که توسط رابطه (1) بیان شده است، به عنوان تابع شایستگی ذرات استفاده می­شود. در ساختار الگوریتم پیشنهادی، موقعیت هر ذره، تنها بیانگر متغیرهای Zt و Xjtاست. بنابراین، برای محاسبه تابع هدف لازم است، wpjt، Ipjt و Ppt و در پی آن مقدار تابع هدف بر اساس Zt و Xjtمحاسبه شود. به ازای مقادیر ثابتی از متغیرهای Zt و Xjt، می‌توان مقادیر مختلفی برای متغیرهای wpjt، Ipjt و Ppt بدست آورد که از میان این مقادیر، wpjt، Ipjt و Ppt را باید به گونه­ای تعیین کرد که ضمن برقراری محدودیت­های مسأله، کمترین هزینه­ها را نیز دربرداشته باشد. در الگوریتم پیشنهادی برای یافتن مقادیر بهینه مذکور، از یک مدل برنامه­ریزی خطی استفاده می­شود. بدین منظور، با توجه به روابط (2)و (3) و با توجه به اینکهخواهیم داشت:

(13)  

(14)          

با جایگذاریاز روابط (13)و (14) در رابطه (1) و با مشخص بودن زمان­های تولید و ارسال، پس از ساده­سازی خواهیم داشت:

(15)

همچنین با جایگذاریاز روابط (13)و (14) در رابطه (6) خواهیم داشت:

(16)      

(17)  

 

همچنین با توجه به روابط (13)، (14) و چون داریم ، می­توان محدودیت­های (2) و (3) را با محدودیت­های (18) و (19) جایگزین کرد:

(18)           

(19)                 

بنابراین از روابط فوق و با توجه به معادلات (1)الی(8)، مدل ساده شده به صورت زیر است:

(20)

(21)           

(22)         

(23)                                  

 

(24)    

(25)  

 

(26)                 

(27)                                  

و معادلات (16) و (19).

که در رابطه (22)،
  است. مدل ساده شده فوق به تعداد P(N+1)T متغیر و حداکثرN(P+1) محدودیت کمتر از مدل اصلی دارد. حال در هر مرحله از الگوریتم، با مشخص بودن زمان­های تولید و ارسال (و) و با استفاده از مدل مذکور، می­توان متغیرهای wpjt و Ppt و نیز مقادیر شایستگی هریک از ذرات را محاسبه و با توجه به مقادیر حاصل از حل مدل برنامه­ریزی خطی برای متغیرهای و، در صورت نیاز، اصلاح لازم در زمان­های تولید و ارسال را نیز اعمال کرد -زیرا ممکن است با حل مدل برنامه­ریزی خطی مقدار بهینه برخی از مقادیر تولید یا ارسال صفر بدست آید که در این صورت زمان ارسال متناظر با آن نیز صفر خواهد شد.

محاسبه مقدار شایستگی ذرات با استفاده از حل مدل برنامه­ریزی خطی تا حدی به زمان محاسباتی بالایی نیاز دارد، بنابراین برای افزایش سرعت الگوریتم، فقط در مرحله دوم الگوریتم از این روش استفاده می­شود و برای محاسبه مقدار شایستگی در مرحله نخست الگوریتم، روشی تقریبی که در ادامه تشریح خواهد شد، توسعه داده شده است:

1-    با استفاده از رابطه (28) مقدار اولیه wpjt را محاسبه کنید:

(28)               

که در رابطه مذکور، r نشانگر دوره بعدی است که ارسال بعدی در آن صورت می­گیرد.

2-    با استفاده از روابط (29) تا (31) مقدار بدست آمده برای wpjt را تعدیل کنید به گونه­ای که تمام مقادیر wpjt شدنی شود ():

(29)   

(30)                

 

(31)

 

3-    مقدار اولیه Ppt را با استفاده از رابطه (32) محاسبه کنید:

(32)                

که در رابطه مذکور، r نشانگر دوره بعدی است که تولید بعدی در آن صورت می­گیرد.

4-    با استفاده از روابط (33) تا (35) مقدار بدست آمده برای Ppt را تعدیل کنید به گونه­ای که تمام مقادیر Ppt شدنی شود ():

(33)               

(34)   

(35)

 

5-    با توجه به روابط (36) و (37) مقدار Ip0t و Ipjt را محاسبه کنید:

(36)        

(37)           

6-    با به­کارگیری روابط (38) تا (40) میزان نشدنی بودن جواب­های حاصل را بدست آورید:

(38)                          

(39)                    

(40)

که در روابط بالا،  ضریب جریمه برای تخطی از محدودیت k ام است.

7-    با استفاده از رابطه (41) مقدار شایستگی ذره را محاسبه کنید:

(41)                

که در رابطه فوق، Z(X) بیانگر تابع هدف حاصل از رابطه (1) است.

در پایان رویه تقریبی نیز، مانند حل مدل برنامه­ریزی خطی، با توجه به مقادیر حاصل برای متغیرهای  و ، در صورت نیاز اصلاح لازم در زمان­های تولید و ارسال اعمال می­شود. در مرحله دوم الگوریتم، در صورتی که ذره­ای نشدنی باشد، حل مدل برنامه­ریزی خطی امکانپذیر نیست. از این رو، در مرحله دوم الگوریتم، در صورتی که ذره­ای نشدنی باشد نیز، از رویه تقریبی فوق برای محاسبه مقدار شایستگی آن ذره استفاده می­شود.

تفاوت سرعت الگوریتم پیشنهادی در دو حالت محاسبه مقدار شایستگی ذرات با استفاده از حل مدل برنامه­ریزی خطی در هر دو مرحله (LP-IPSO) و محاسبه مقدار شایستگی ذرات با استفاده از رویه تقریبی فوق در مرحله نخست و استفاده از مدل برنامه­ریزی خطی فقط در مرحله دوم الگوریتم (IPSO) را می­توان با استفاده از چند مثال بررسی کرد. میانگین نتایج حاصل از 10 بار اجرای الگوریتم پیشنهادی در دو حالت مذکور برای حل چند مسأله نمونه (مسایل نمونه در بخش 4 تشریح خواهند شد)، در جدول 1 ارائه شده است.

جدول 1. نتایج حاصل از حل مسایل نمونه با استفاده از LP-IPSO و IPSO

شماره مسأله

تابع هدف

نشدنی بودن

زمان (ثانیه)

IPSO

LP-IPSO

IPSO

LP-IPSO

IPSO

LP-IPSO

1

2358

2358

0

0

4.5

160.4

2

3487

3482

0

0

5.8

233

3

3533.6

3496.3

0

0

7.1

278.6

4

5127.6

5048

0

0

8.7

418.4

 

 

شکل 6. میزان اختلاف زمان محاسباتی موردنیاز برای حل مسایل نمونه با استفاده از LP-IPSO و IPSO

 

همانگونه که در جدول 1 و شکل 6 ارائه شده است، هرچند LP-IPSO به مقادیر تابع هدف بهتری نسبت به IPSO دست یافته است، اما به زمان محاسباتی به مراتب بیشتری نیاز دارد. با افزایش ابعاد مسأله، تفاوت سرعت الگوریتم در دو حالت افزایش می­یابد.

تعیین مناسب ضریب جریمه از اهمیت بسزایی برخوردار است. زیرا، اگر ضریب جریمه بسیار کوچک انتخاب شود، نیروی کافی برای جلوگیری از تخطی از محدودیت­ها وارد نمی­شود. از سوی دیگر، اگر ضریب جریمه بسیار بزرگ انتخاب شود، تمام جواب­های نشدنی حتی در گام­های ابتدایی الگوریتم حذف می­شوند و توانایی جستجوی الگوریتم کاهش می­یابد (انگلبرکت، 2005). از این رو، در ساختار الگوریتم پیشنهادی، به منظور حفظ توانایی جستجوی الگوریتم در گام­های ابتدایی و تاکید بر حذف جواب­های نشدنی در گام­های نهایی، طی گام­های مرحله نخست الگوریتم، با استفاده از رابطه (42) ضریب جریمه با روندی صعودی، افزایش می­یابد.

(42)

در رابطه بالا،ضریب جریمه تخطی از محدودیت k ام در تکرار t ام،میزان افزایش در ضریب جریمه و infeasibility(X) میزان نشدنی بودن بهترین ذره است. ضریب جریمه در طی گام­های مرحله دوم الگوریتم ثابت باقی می­ماند.

 

 

 

3-4. ایجاد همسایگی تصادفی ذرات

در مرحله نخست الگوریتم، به منظور افزایش پراکندگی ذرات و جستجوی مناطق بیشتری از محدوده جواب، با بهره­گیری از عملگر جهش، حول هر ذره همسایگی تصادفی تشکیل می­شود. بدین منظور:

1-    تعداد NSize-1 کپی از ذره ایجاد می­شود.

2-    به تعداد NSize-1 آرایه تصادفی صفر و یک، با تعداد مؤلفه­های برابر با مؤلفه­های یک ذره، ایجاد می­گردد. HDRate درصد از مؤلفه­های آرایه تولیدشده، به صورت تصادفی برابر 1 قرار داده می­شود.

3-    مطابق با مؤلفه­های با مقدار 1 در آرایه تصادفی، مؤلفه متناظر در ذره کپی معکوس می­گردد؛ به عبارت دیگر، مقدار 0 به 1 و بالعکس مقدار 1 به 0 تبدیل می­شود.

نمایی کلی از نحوه انجام عملگر جهش در شکل 7 ارائه شده است.

 

شکل 7. نحوه انجام عملگر جهش

 

3-5. بهبود همسایگی ذرات

پس از ایجاد همسایگی ذرات، طی گام­های الگوریتم، جواب­های موجود در هر همسایگی مستقل از ذرات سایر همسایگی­ها بهبود می­یابند. اما به منظور افزایش میزان بهبود در کیفیت جواب­ها طی گام­های الگوریتم، لازم است میان همسایگی­های مختلف ذرات، به گونه­ای تبادل اطلاعات صورت گیرد. بدین منظور، پس از بهنگام­سازی مقادیر nBest در هر تکرار الگوریتم، با به­کارگیری روش انتخاب مسابقه[46] به تعداد 2SSize ذره از میان ذرات nBest به عنوان والدین[47] انتخاب می­شوند. پس از آن با بهره­گیری از عملگر تقاطع پراکنده که در بخش 3-2 تشریح شد، SSize ذره جدید (فرزند[48]) ایجاد می­شود. ذرات جدید، جایگزین بدترین ذره در هر همسایگی می­شوند.

 

3-6. تشکیل و بهنگام­سازی مجموعه مرجع (RSet)

در ساختار الگوریتم پیشنهادی، در مرحله نخست الگوریتم، به منظور حفظ جواب­های غیرتکراری خوب در طی گام­های الگوریتم و همچنین به منظور فراهم­سازی جواب­های اولیه با تنوع و کیفیت بالا برای گام­های آتی الگوریتم، مجموعه مرجع تشکیل می­گردد. بدین منظور از میان تمامی ذرات غیرتکراری موجود در گروه ذرات، حداکثر به تعداد SSize ذره با بهترین مقادیر شایستگی انتخاب و در مجموعه مرجع قرار داده می­شوند. در تکرارهای بعدی الگوریتم، تنها ذراتی در مجموعه مرجع قرار می­گیرند که تکراری نبوده و به مقادیر شایستگی بهتری دست یافته باشند. رویه تشکیل مجموعه مرجع در دو بخش اول و دوم از مرحله نخست الگوریتم پیشنهادی، یکسان است.

  

3-7. متنوع­سازی ذرات

در ساختار الگوریتم پیشنهادی، به منظور جلوگیری از همگرایی زودهنگام و ایجاد جواب­های اولیه با پراکندگی و کیفیت خوب برای مراحل بعدی الگوریتم، مکانیزم متنوع­سازی ذرات انجام می­گیرد. به همین منظور، در بخش نخست مرحله اول الگوریتم، بهترین ذره در میان مجموعه مرجع انتخاب و به تعداد MRate درصد از ذرات گروه، از ذره منتخب کپی ایجاد می­شود. ذرات کپی با به­کارگیری عملگر جهش که در بخش 3-4 تشریح شد، تغییر می­یابند و جایگزین بدترین ذرات گروه می­شوند. با استفاده از این روش، با اعمال تغییر تصادفی محدود در بهترین ذره گروه، با افزایش مؤلفه­های خوب در میان جواب­ها، فرآیند جستجو بهبود می­یابد.

 

3-8. جستجوی محلی

در مرحله دوم از ساختار الگوریتم پیشنهادی، پس از انتخاب بهترین ذره، کیفیت ذره منتخب با استفاده از جستجوی محلی بهبود می­یابد. برای این هدف، یکی از مؤلفه­های ذره به صورت تصادفی انتخاب و مقدار آن مؤلفه معکوس می­شود؛ به عبارت دیگر، اگر مؤلفه برابر 1 باشد به 0 و بالعکس اگر 0 باشد به 1 تبدیل می­شود. بدین ترتیب ذره جدیدی تنها با یک مؤلفه متفاوت نسبت به ذره اولیه بدست می­آید. با استفاده از مدل برنامه­ریزی خطی که در بخش 3-3 ارائه شد، مقدار شایستگی ذره جدید محاسبه می­گردد. اگر مقدار شایستگی ذره جدید از ذره اولیه بهتر باشد، ذره جدید جایگزین ذره اولیه می­شود و در غیر اینصورت ذره جدید حذف و ذره اولیه بدون تغییر باقی می­ماند. رویه مذکور حداکثر پس از P2 مرتبه تکرار یا در صورتی که در بهترین مقدار شایستگی که ذره بدان دست یافته است، طی CR2 تکرار متوالی بهبودی حاصل نشود، الگوریتم متوقف می­گردد.

 

4- نتایج محاسباتی

با استفاده از مسایل نمونه تصادفی متعددی، عملکرد الگوریتم پیشنهادی در مقایسه با الگوریتم بهینه­سازی گروه ذرات، الگوریتم ژنتیک و نرم­افزار Lingo (تنها برای مسایل با ابعاد کوچک) ارزیابی می­شود. با مقایسه عملکرد الگوریتم پیشنهادی با الگوریتم بهینه­سازی گروه ذرات، می­توان به میزان بهبود حاصل از ساختار جدید طراحی شده برای الگوریتم پیشنهادی نسبت به ساختار پایه­ای الگوریتم بهینه­سازی گروه ذرات پی برد. الگوریتم ژنتیک نیز، از جمله معروفترین و پرکاربردترین الگوریتم­های فراابتکاری است که توسط پ‍ژوهش­های متعددی در حل مسایل بهینه­سازی مختلفی به­کار گرفته شده است. به همین دلیل، برای ارزیابی عملکرد الگوریتم پیشنهادی، از دو الگوریتم فوق به عنوان الگوریتم­های معیار استفاده شده است.

الگوریتم­های پیشنهادی و بهینه­سازی گروه ذرات، در محیط Matlab 7.5 برنامه­نویسی و الگوریتم ژنتیک با استفاده از جعبه ابزار نرم­افزار Matlab 7.5 پیاده­سازی شده است. تمامی الگوریتم­های مذکور در کامپیوتری با مشخصات پردازنده AMD، 64 بیتی، 3 گیگاهرتز، ویندوز XP و 2 گیگابایت رم، اجرا شده­اند.

 

 

4-1- چگونگی ایجاد مسایل نمونه

مسایل نمونه در دو گروه با ابعاد کوچک و بزرگ به صورت تصادفی ایجاد شده­اند. ابعاد مسایل، تعداد دوره­های برنامه­ریزی و تعداد محصولات در جدول 2 ارائه شده است. همانگونه که در جدول مذکور نشان داده شده است، در تولید مسایل نمونه، هزینه حمل محصولات از مسایل معیار مسأله مسیریابی وسیله نقلیه با ظرفیت محدود[49] استخراج شده است (مجموعه داده­های مسیریابی وسیله نقلیه[50]). اطلاعات لازم برای تولید مسایل نمونه در جدول 3 ارائه شده است. در این مقاله، همانند بودیا و پرینز (2009) ظرفیت تولید و وسایل حمل متناسب با تقاضا، هزینه ثابت راه­اندازی تولید و همچنین ظرفیت انبار تولیدکننده متناسب با ظرفیت تولید در نظر گرفته شده است. ظرفیت انبار خرده‌فروشان نیز با ظرفیت وسایل حمل متناسب است. بدین منظور، ظرفیت وسایل حمل 2 برابر () میانگین تقاضای یک خرده­فروش و ظرفیت تولید 5/3 برابر میانگین تقاضای روزانه تمام خرده­فروشان () تعیین شده است. هزینه ثابت راه­اندازی تولید 5/1 برابر () و ظرفیت انبار تولیدکننده برابر با () ظرفیت تولید در نظر گرفته شده است. ظرفیت انبار تمام خرده­فروشان یکسان و برابر با ظرفیت وسیله حمل () تعیین شده است. در جدول 3، dp بیانگر میانگین تقاضای محصول p ام در یک دوره است. در تولید مسایل نمونه با 5 محصول، 5 سطح متفاوت تقاضا شامل تقاضای خیلی کم (d1jt)، کم (d2jt)، متوسط (d3jt)، زیاد (d4jt) و خیلی زیاد (d5jt) مدنظر قرار گرفته است که مقدار تقاضای هر محصول در هر دوره به صورت یکنواخت در بازه ارائه شده در جدول 3، ایجاد می­شود. برای مسایل با 3 محصول، 3 سطح تقاضای خیلی کم (d1jt)، متوسط (d3jt) و خیلی زیاد (d5jt) لحاظ شده است.

 

 

جدول 2. چگونگی ایجاد مسایل نمونه

شماره مسأله

بعد مسأله

(T×N×P)

نام مسأله

شماره مسأله

بعد مسأله

(T×N×P)

نام مسأله

1

کوچک (3×1×10)

P-n16-k8

13

بزرگ (3×101×10)

M-n101-k10

2

کوچک (5×1×10)

P-n16-k8

14

بزرگ (5×101×10)

M-n101-k10

3

کوچک (3×1×15)

P-n16-k8

15

بزرگ (3×101×15)

M-n101-k10

4

کوچک (5×1×15)

P-n16-k8

16

بزرگ (5×101×15)

M-n101-k10

5

کوچک (3×5×10)

P-n16-k8

17

بزرگ (3×121×10)

M-n121-k7

6

کوچک (5×5×10)

P-n16-k8

18

بزرگ (5×121×10)

M-n121-k7

7

کوچک (3×5×15)

P-n16-k8

19

بزرگ (3×121×15)

M-n121-k7

8

کوچک (5×5×15)

P-n16-k8

20

بزرگ (5×121×15)

M-n121-k7

9

کوچک (3×10×10)

P-n16-k8

21

بزرگ (3×151×10)

M-n151-k12

10

کوچک (5×10×10)

P-n16-k8

22

بزرگ (5×151×10)

M-n151-k12

11

کوچک (3×10×15)

P-n16-k8

23

بزرگ (3×151×15)

M-n151-k12

12

کوچک (5×10×15)

P-n16-k8

24

بزرگ (5×151×15)

M-n151-k12

 

جدول 3. نحوه تولید پارامترهای مسایل نمونه

نام پارامتر

مقدار پارامتر

نام پارامتر

مقدار پارامتر

Q

 

d1jt

U[1,3]

Pmax

 

d2jt

U[7,10]

ft

 

d3jt

U[15,20]

Imax0

 

d4jt

U[25,35]

Imaxj

 

d5jt

U[45,60]

hpj

1

ap, kp

1

 

 


4-2- مفروضات و پارامترهای الگوریتم­ها

برای تنظیم الگوریتم­ها، آزمایش­های زیادی با مجموعه مقادیر مختلف پارامترها انجام گرفته است که در پایان با استفاده از مجموعه مقادیر زیر بهترین نتایج حاصل شده است:

الگوریتم ژنتیک: تعداد تکرارها، اندازه جمعیت، تعداد جواب نخبه[51] و ضریب جریمه برای مسایل نمونه کوچک به ترتیب برابر با 200، 100، 10 و 100 و برای مسایل نمونه بزرگ به ترتیب برابر با 1500، 350، 50 و 5^10 قرار داده شده است. سایر پارامترها همانند گزینه­های پیش فرض جعبه­ابزار Matlab در نظر گرفته شده است.

الگوریتم بهینه­سازی گروه ذرات: تعداد تکرارها، تعداد ذرات، ضرایب یادگیری (c1 و c2Vmax و ضریب جریمه برای مسایل نمونه کوچک به ترتیب برابر با 200، 100، 2، 2، 6 و 100 و برای مسایل نمونه بزرگ به ترتیب برابر با 1500، 350، 2، 2، 6 و 5^10 قرار داده شده است. ضریب اینرسی به صورت خطی از مقدار 1.4 به 0.9 در طول اجرای الگوریتم کاهش می­یابد.

الگوریتم بهینه­سازی گروه ذرات بهبودیافته: تعداد تکرارهای بخش اول (P1S2) و دوم (P1S1)، مرحله اول (P1) و دوم (P2)، حداکثر تکرار بدون بهبود در بهترین جواب برای بخش اول (C11) و دوم (C12) و مرحله دوم (C2)، تعداد ذرات (SSize)، اندازه همسایگی ذرات (NSize)، ضرایب یادگیری (c1 و c2Vmax در بخش اول و دوم، HDRate، MRate، ضریب جریمه اولیه برای بخش اول () و دوم () و مرحله دوم () و میزان افزایش در ضریب جریمه در هر گام () برای مسایل نمونه کوچک به ترتیب برابر است با 50، 50، 10، 250، 25، 25، 50، 20، 10، 2، 2، 3، 6، 0.07، 0.1، 10، 10، 75 و 0.1 و برای مسایل نمونه بزرگ به ترتیب برابر با 75، 75، 10، 350، 35، 35، 75، 30، 10، 2، 2، 6، 6، 0.1، 0.1، 100، 100، 750 و 0.5 قرار داده شده است. ضریب اینرسی به صورت خطی از مقدار 1.4 به 0.9 در طول اجرای الگوریتم کاهش می­یابد.

 

4-3- نتایج عددی

هریک از الگوریتم­ها برای حل هر مسأله نمونه، 10 بار اجرا شده­اند که با تحلیل نتایج حاصل، از دیدگاه 3 معیار عملکردی[52] به ارزیابی عملکرد الگوریتم پیشنهادی پرداخته شده است: 1) میانگین مقادیر تابع هدف بدست آمده، 2) میانگین زمان محاسباتی و 3) میزان استواری[53] نتایج.

هر اندازه تغییرات مقادیر تابع هدف در اجراهای مختلف الگوریتمی مشخص، کمتر باشد؛ میزان استواری الگوریتم بیشتر خواهد بود. معیار عملکردی میزان استواری نتایج با استفاده از رابطه (43) محاسبه می­گردد (انگلبرکت، 2005):

(43)             

در رابطه فوق، میانگین و واریانس مقادیر تابع هدف بدست آمده در اجراهای مختلف الگوریتم است. هر اندازه واریانس نتایج و به نوبه آن بازه ارائه شده در رابطه (43) کوچکتر باشد، میزان استواری الگوریتم بیشتر است.

در جداول 4 و 5، میانگین مقادیر تابع هدف بدست آمده و زمان­های محاسباتی صرف شده توسط سه الگوریتم برای حل هر مسأله، ارائه شده است. در برخی از موارد، در حل یک مسأله مشخص، پس از اختتام الگوریتم، جوابی نشدنی حاصل شده است. در این صورت، میزان تخطی از محدودیت­ها در هر اجرای الگوریتم برای آن مسأله تعیین و ذخیره و در پایان میانگین میزان تخطی از محدودیت­ها در تمام اجراهای الگوریتم برای آن مسأله محاسبه می­شود. در جداول 4 و 5، میانگین مقدار تخطی از محدودیت­ها در حل یک مسأله، در ستون نشدنی بودن ارائه شده است. با توجه به جداول مذکور، الگوریتم بهینه­سازی گروه ذرات بهبودیافته پیشنهادی در تمام مسایل نمونه نسبت به الگوریتم­های معیار، به نتایج بهتری دست یافته است. در مسایل نمونه با ابعاد کوچک، الگوریتم پیشنهادی توانسته است به مقادیری بسیار نزدیک به بهینه و در مسأله نمونه 1 به مقدار بهینه دست یابد. همچنین، الگوریتم پیشنهادی در تمام مسایل نمونه به جوابی شدنی رسیده است، در حالی که الگوریتم بهینه­سازی گروه ذرات در برخی از موارد به نتایج نشدنی منتج شده است.

 

 

جدول 4. میانگین مقدار تابع هدف و زمان محاسباتی حل مسایل نمونه تصادفی با ابعاد کوچک

شماره مسأله

Lingo

IPSO

GA

PSO

تابع هدف

زمان (ثانیه)

تابع هدف

نشدنی بودن

زمان (ثانیه)

تابع هدف

نشدنی بودن

زمان (ثانیه)

تابع هدف

نشدنی بودن

زمان (ثانیه)

1

2358

1

2358

0

4.5

2390.2

0

1.2

2394.8

0

0.9

2

3482

1

3487

0

5.8

3533.9

0

1.4

3497.2

0

1.3

3

3481

3

3533.6

0

7.1

3637.4

0

1.4

3656.9

0

1.4

4

5048

6

5127.6

0

8.8

5308.4

0

1.9

5409.5

0

1.8

5

12554.4

26

12581.4

0

14.5

13301.4

0

3.1

13230.2

0

3.3

6

18150.4

21

18238

0

22.1

18771.9

0

4.9

18723.7

0

4.9

7

16438.6

*10800

18396

0

23.8

18990.3

0

4.8

18888.7

0

5

8

26183

6858

26435.5

0

39.3

28584.2

0

6.9

27947.4

0.1

7.6

9

25184.6

510

25565.8

0

28.4

26992.8

0

5.9

26802

0.2

6.2

10

36411.6

7943

37163.7

0

43.4

39228.4

0

8.9

39337.8

0

9.8

11

36608

10800

37551.9

0

49.1

39929.2

0

9

40152.8

0.4

9.3

12

52418

10800

53141.2

0

73.7

57121.9

0

15

57275.5

0.3

15.6

میانگین

19859.8

3980.8

20298.3

0

26.7

21482.5

0

5.4

21443

0.1

5.6

- بهترین مقدار بدست آمده برای هر مسأله نمونه پررنگ شده است.

*  نرم­افزار به جواب بهینه نرسیده و متوقف شده است، بهترین جواب بدست آمده گزارش شده است.

 

جدول 5. میانگین مقدار تابع هدف و زمان محاسباتی حل مسایل نمونه تصادفی با ابعاد بزرگ

شماره مسأله

IPSO

GA

PSO

تابع هدف

نشدنی بودن

زمان (ثانیه)

تابع هدف

نشدنی بودن

زمان (ثانیه)

تابع هدف

نشدنی بودن

زمان (ثانیه)

13

256257

0

863

273682

0

1760

272084

0

1894

14

370912

0

1390

385176

0

2706

386357

0

2801

15

385860

0

1681

406516

0

2668

421797

0

2865

16

558963

0

2698

578444

0

4136

594274

0

4228

17

312315

0

1096

326621

0

2114

335792

0

2262

18

449007

0

1964

462582

0

3253

477178

0

3356

19

472147

0

1989

498801

0

3210

590122

1

3425

20

676639

0

3205

685537

0

4921

816348

2

5077

21

369996

0

1198

380098

0

2636

404811

0

2811

22

541093

0

2224

547732

0

4041

573921

0

4278

23

556473

0

2484

575109

0

4005

870552

35

4253

24

812303

0

4127

830505

0

6126

1276856

51

6339

میانگین

480164

0

2077

495900

0

3465

585008

7

3632

- بهترین مقدار بدست آمده برای هر مسأله نمونه پررنگ شده است.

 

 

برای آزمودن وجود تفاوتی معنی­دار میان عملکرد الگوریتم پیشنهادی با الگوریتم­های معیار از تحلیل واریانس دو طرفه[54] بر روی مقادیر تابع هدف بدست آمده توسط سه الگوریتم در حل مسایل با ابعاد مختلف استفاده می­کنیم. بدین منظور الگوریتم­ها به عنوان تیمار[55] و مسایل با ابعاد مختلف به عنوان بلوک[56] در نظر گرفته شده­اند. تحلیل واریانس با استفاده از نرم­افزار Minitab 14 انجام شده است که نتایج تحلیل واریانس مسایل با ابعاد کوچک و بزرگ، به ترتیب در شکل­های 8 و 9 ارائه شده است. با توجه به شکل­های 8 و 9، در هر دو گروه مسایل نمونه با ابعاد کوچک و بزرگ در سطح معنی­دار 01/0 () می­توان یکسان نبودن عملکرد الگوریتم­ها را نتیجه گرفت.

 

شکل 8. تحلیل واریانس دو طرفه بر روی نتایج حاصل از حل مسایل با ابعاد کوچک

 

شکل 9. تحلیل واریانس دو طرفه بر روی نتایج حاصل از حل مسایل با ابعاد بزرگ

مقایسه میان زمان­های محاسباتی صرف شده توسط هر سه الگوریتم در حل مسایل با ابعاد کوچک و بزرگ به ترتیب در شکل­های 10 و 11 نمایش داده شده است.

 

شکل 10. مقایسه زمان محاسباتی الگوریتم پیشنهادی و الگوریتم­های معیار در حل مسایل با ابعاد کوچک

 

شکل 11. مقایسه زمان محاسباتی الگوریتم پیشنهادی و الگوریتم­های معیار در حل مسایل با ابعاد بزرگ

الگوریتم پیشنهادی در حل مسایل با ابعاد کوچک به زمان محاسباتی بیشتری نسبت به الگوریتم­های معیار نیاز دارد. در حالی که در حل مسایل با ابعاد بزرگ زمان محاسباتی الگوریتم پیشنهادی در قیاس با الگوریتم­های معیار به مراتب کمتر است. این امر را اینگونه می­توان تفسیر کرد که در مسایل نمونه با ابعاد کوچک، هریک از الگوریتم­ها برای دستیابی به جوابی مناسب، می­بایست فضای جستجوی کوچکی را پوشش دهند. بنابراین، الگوریتم­های معیار که تنها از یک بخش تشکیل شده­اند در زمان محاسباتی اندکی قادرند به جوابی مناسب همگرا شوند. با وجود این، ساختار الگوریتم پیشنهادی از چند بخش تشکیل شده است و هر بخش به حداقل زمان محاسباتی نیاز دارد. به همین دلیل، الگوریتم پیشنهادی در حل مسایل با ابعاد کوچک، برای رسیدن به جوابی مناسب زمان محاسباتی بیشتری نسبت به دو الگوریتم معیار صرف می­کند. اما در مسایل با ابعاد بزرگ، ساختار چند بخشی الگوریتم پیشنهادی موجب پوشش مؤثرتر و کاراتر فضای جستجو می­شود. از این رو، الگوریتم پیشنهادی می­تواند در زمان محاسباتی کمتری به نتایج مطلوب­تری نسبت به دو الگوریتم معیار دست یابد.

میزان استواری هر سه الگوریتم در حل مسایل کوچک و بزرگ، به ترتیب در جداول 6 و 7 ارائه شده است. مطابق با جداول مذکور، به طور میانگین الگوریتم پیشنهادی در مقایسه با دو الگوریتم معیار در هر دو گروه مسایل کوچک و بزرگ از تغییرپذیری کمتر و به نوبه آن از استواری بیشتری برخوردار است.

عملکرد بهتر الگوریتم پیشنهادی را می­توان ناشی از ساختار جدید طراحی شده برای آن دانست. این مطلب از چهار بعد قابل بررسی است: 1) نحوه تولید جواب­های اولیه در الگوریتم پیشنهادی به گونه­ای است که جواب­هایی خوب برای مراحل بعدی فراهم می­گردد، 2) با تشکیل مجموعه مرجع در مرحله نخست الگوریتم، جواب­های با کیفیت و منحصربفردی که در هر تکرار تولید شده­اند، حفظ می­شوند. بدین ترتیب جواب­هایی خوب و متنوع برای گام­های بعدی بدست می­آید و جواب­های نامناسب از چرخه الگوریتم حذف می­شوند، 3) با تشکیل همسایگی حول ذرات در بخش نخست الگوریتم، همسایگی هر ذره با دقت بیشتری جستجو و امکان بهبود در ذرات بررسی می­شود، 4) استفاده از رویه تقریبی در مرحله نخست و مدل برنامه­ریزی خطی در مرحله دوم الگوریتم پیشنهادی برای محاسبه تابع هدف ذرات، به شدت بر سرعت الگوریتم می­افزاید. بنابراین، می­توان نتیجه گرفت که الگوریتم پیشنهادی با سرعت و ثبات بیشتری به نتایج بهتری دست می­یابد.

 

4-4- تحلیل حساسیت

تعداد دفعات و زمان بهینه ارسال محصولات به خرده­فروشان و تولید محصولات توسط تولیدکننده و در پی آن ساختار مسایل نمونه، به دو عامل نسبت مجموع هزینه­های راه­اندازی و ارسال به هزینه­های نگهداری و نسبت ظرفیت تولید، انبار یا وسایل حمل به تقاضا، وابسته است. با کاهش نسبت مجموع هزینه­های راه­اندازی و ارسال به هزینه­های نگهداری و نیز با کاهش نسبت ظرفیت تولید، انبار یا وسایل

 

 

 

جدول 6. میزان استواری الگوریتم­ها در حل مسایل نمونه تصادفی با ابعاد کوچک

شماره مسأله

IPSO

استواری

GA

استواری

PSO

استواری

کران بالای بازه

کران پایین بازه

طول بازه

کران بالای بازه

کران پایین بازه

طول بازه

کران بالای بازه

کران پایین بازه

طول بازه

1

2358

2358

0

2414.2

2366.2

48

2446.8

2342.8

104

2

3487

3487

0

3597.9

3469.9

128

3519.2

3475.2

44

3

3588.6

3478.6

110

3700.4

3574.4

126

3782.9

3530.9

252

4

5230.6

5024.6

206

5393.4

5223.4

170

5516.5

5302.5

214

5

12600.4

12562.4

38

13594.4

13008.4

586

13502.2

12958.2

544

6

18318

18158

160

18951.9

18591.9

360

18756.7

18690.7

66

7

18531

18261

270

19253.3

18727.3

526

19010.7

18766.7

244

8

26903.5

25967.5

936

29125.2

28043.2

1082

28664.4

27230.4

1434

9

25585.8

25545.8

40

27758.8

26226.8

1532

27529

26075

1454

10

37191.7

37135.7

56

39778.4

38678.4

1100

40181.8

38493.8

1688

11

38337.9

36765.9

1572

41481.2

38377.2

3104

40893.8

39411.8

1482

12

54130.2

52152.2

1978

58538.9

55704.9

2834

58079.5

56471.5

1608

میانگین

20521.9

20074.7

447.2

21965.7

20999.3

966.3

21823.6

21062.5

761.2

جدول 7. میزان استواری الگوریتم­ها در حل مسایل نمونه تصادفی با ابعاد بزرگ

شماره مسأله

IPSO

استواری

GA

استواری

PSO

استواری

کران بالای بازه

کران پایین بازه

طول بازه

کران بالای بازه

کران پایین بازه

طول بازه

کران بالای بازه

کران پایین بازه

طول بازه

13

256668

255846

822

282911

264453

18458

277379

266789

10590

14

371296

370528

768

401076

369276

31800

390853

381861

8992

15

391425

380295

11130

421730

391302

30428

431619

411975

19644

16

559402

558524

878

599754

557134

42620

603362

585186

18176

17

313514

311116

2398

337422

315820

21602

342242

329342

12900

18

450226

447788

2438

473314

451850

21464

494191

460165

34026

19

473220

471074

2146

518520

479082

39438

659377

520867

138510

20

678929

674349

4580

706873

664201

42672

931064

701632

229432

21

370509

369483

1026

386830

373366

13464

423408

386214

37194

22

541587

540599

988

554919

540545

14374

598245

549597

48648

23

564760

548186

16574

590488

559730

30758

944539

796565

147974

24

812893

811713

1180

840330

820680

19650

1353288

1200424

152864

میانگین

482036

478292

3744

509514

482287

27227

620797

549218

71579

 

 

 

حمل به تقاضا، ارسال و تولید محصولات در فواصل

زمانی کوتاه­تر رخ می­دهد. بنابراین، برای ارزیابی جامع­تر عملکرد الگوریتم پیشنهادی، نسبت به مقادیر پارامترهای ، ، ،  و  تحلیل حساسیت صورت می­گیرد. بدین منظور، با در نظر گرفتن ،  و ، مطابق با جدول 8، 5 مقدار متفاوت به ترتیب برای پارامترهای ، ، ،  و  تعیین می­شود.

جدول 8. مقادیر پارامترهای مسأله در حالات مختلف برای تحلیل حساسیت

حالت

پارامتر

1

2

3

4

5

 

1.5

1.75

2

2.25

2.5

 

1.5

2.5

3.5

4.5

5.5

 

0.5

1

1.5

2

2.5

 

0.5

0.75

1

1.25

1.5

 

0.5

0.75

1

1.25

1.5

 

هریک از الگوریتم­ها برای حل هر حالت، 10 بار اجرا شده­اند که میانگین مقادیر تابع هدف بدست آمده و نیز درصد اختلاف میان میانگین مقادیر تابع هدف بدست آمده توسط هریک از الگوریتم­های معیار با الگوریتم پیشنهادی، در نمودارهای 12 الی 16 نشان داده شده است. این نکته قابل ذکر است که هر سه الگوریتم در تمامی حالت­ها به جوابی شدنی دست یافته­اند.

با توجه به نمودارهای مذکور، الگوریتم پیشنهادی تنها در یک مورد نسبت به الگوریتم ژنتیک از عملکرد اندکی ضعیف­تر برخوردار است (در تحلیل حساسیت نسبت به پارامتر  در حالت نخست) و در سایر موارد نسبت به دو الگوریتم معیار از عملکرد مناسب­تری برخوردار است. بنابراین، در کل می­توان عملکرد بهتر الگوریتم پیشنهادی را در ساختارهای مختلف مسایل نمونه نتیجه گرفت.

 

شکل 12. تحلیل حساسیت نسبت به پارامتر Q

 

 

 

شکل 13. تحلیل حساسیت نسبت به پارامتر

 

 

 

 

 

 شکل 14. تحلیل حساسیت نسبت به پارامتر

 

 

شکل 15. تحلیل حساسیت نسبت به پارامتر

 

 

شکل 16. تحلیل حساسیت نسبت به پارامتر

 

 

5- نتیجه­ گیری

این مقاله به بررسی مسأله مسیریابی-موجودی در حالت چند دوره­­ای چند محصولی و با هدف کمینه­سازی مجموع هزینه­های سیستم شامل هزینه‌های راه­اندازی، توزیع و نگهداری موجودی پرداخته است. در مسأله بررسی شده در این مقاله، فرض شده است که چند محصول در یک کارخانه ساخته و با استفاده از ناوگانی از وسایل حمل همسان با ظرفیت محدود برای مجموعه­ای از خرده­فروشان با شیوه توزیع مستقیم، ارسال می­شود. همچنین، هر خرده­فروش در هر دوره تنها باید توسط یک وسیله حمل و یکبار خدمت­دهی شود. تولیدکننده و خرده‌فروشان به نگهداری موجودی با ظرفیت محدود می­پردازند، ظرفیت تولید محدود است و کمبود مجاز نیست. به دلیل پیچیدگی محاسباتی بالای مسأله، استفاده از روش­های حل دقیق به ویژه برای مسایل با ابعاد بزرگ در زمان محاسباتی معقول امکانپذیر نیست. بنابراین، در این مقاله الگوریتم حل تقریبی کارایی مبتنی بر الگوریتم بهینه­سازی گروه ذرات توسعه داده شده است. کارایی الگوریتم پیشنهادی با استفاده از مسایل نمونه متعددی که به صورت تصادفی ایجاد شده­اند با الگوریتم­های ژنتیک و بهینه­سازی گروه ذرات مقایسه شده است. نتایج عددی بیانگر عملکرد بهتر الگوریتم پیشنهادی است.

در تحقیقات آتی، توسعه مسأله مذکور برای سایر شرایطی که می­تواند در صنعت وجود داشته باشد، بسیار مفید خواهد بود. در این راستا پیشنهاد می­گردد مسأله برای حالتی که کمبود مجاز باشد یا امکان ارسال چندین وسیله حمل در هر دوره زمانی برای یک خرده­فروش وجود داشته باشد، توسعه داده شود.



[1] Inventory Routing Problem (IRP)

[2] Vehicle Routing Problem (VRP)

[3] Cordeau et al.

[4] Vendor Managed Inventory (VMI)

[5] Andersson et al.

[6] Direct shipment

[7] Li et al.

[8] Campbell and Hardin

[9] Greedy

[10] Cheng and Duran

[11] Kleywegt et al.

[12] Barnes-Schuster and Bassok

[13] NP-hard

[14] Lenstra and Rinnooy

[15] Aziz and Moin

[16] Moin et al.

[17] Boudia et al.

[18] Greedy randomized adaptive search procedure (GRASP)

[19] Boudia and Prins

[20] Memetic algorithm

[21] Zhao et al.

[22] Variable Large Neighborhood Search (VLNS)

[23] Tabu Search (TS)

[24] Esparcia-Alcazar et al.

[25] Rusdiansyah and Tsao

[26] Insertion

[27] The capacitated lot size problem

[28] Bitran and Yanasse

[29] Exact methods

[30] Improved Particle Swarm Optimization (IPSO)

[31] Particle Swarm Optimization (PSO)

[32] Mixed integer programming

[33] Kennedy and Eberhart

[34] Swarm

[35] Particle

[36] Fitness value

[37] Fitness function

[38] Velocity

[39] Crossover

[40] Mutation

[41] Hu et al.

[42] Jiao et al.

[43] Engelbrecht

[44] Solution representation

[45] Scattered crossover

[46] Tournament selection

[47] Parent

[48] Child

[49] Capacitated vehicle routing problem (CVRP)

[50] Vehicle routing data sets

[51] Elite count

[52] Performance measure

[53] Robustness

[54] Two-way analysis of variance (ANOVA)

[55] Treatment

[56] Block

Andersson, H., Hoff, A., Christiansen, M., Hasle, G., & Løkketangen, A. (2010). Industrial Aspects and Literature Survey: Combined Inventory Management and Routing. Computers & Operations Research, 37(9), 1515-1536.

Aziz, N. A. B., & Moin, N. H. (2007). Genetic Algorithm Based Approach for the Multi Product Multi Period Inventory Routing Problem. Paper presented at the Proceedings of the 2007 IEEE IEEM, Singapour.

Barnes-Schuster, D., & Bassok, Y. (1997). Direct Shipping and the Dynamic Single-Depot/Multi-Retailer Inventory System. European Journal of Operational Research, 101(3), 509-518.

Bitran, G. R., & Yanasse, H. H. (1981). Computational Complexity of the Capacitated Lot Size Problem. Cambridge, MA: Massachusetts Institute of Technology.

Boudia, M., Louly, M. A. O., & Prins, C. (2007). A Reactive Grasp and Path Relinking for a Combined Production–Distribution Problem. Computers & Operations Research, 34(11), 3402-3419.

Boudia, M., & Prins, C. (2009). A Memetic Algorithm with Dynamic Population Management for an Integrated Production-Distribution Problem. European Journal of Operational Research, 195(3), 703-715.

Campbell, A. M., & Hardin, J. R. (2005). Vehicle Minimization for Periodic Deliveries. European Journal of Operational Research, 165(3), 668-684.

Cheng, L., & Duran, M. A. (2004). Logistics for World-Wide Crude Oil Transportation Using Discrete Event Simulation and Optimal Control. Computers & Chemical Engineering, 28(6-7), 897-911.

Cordeau, J.-F., Laporte, G., Savelsbergh, M. W. P., & Vigo, D. (2007). Vehicle Routing. In C. Barnhart & G. Laporte (Eds.), Hanbook in Operations Research and Management Science (Vol. 14, Transportation, (367-428).

Engelbrecht, A. P. (2005). Fundamentals of Computational Swarm Intelligence. West Sussex, England: John wiley & Sons, Ltd.

Esparcia-Alcazar, A. I., Cardos, M., & Merelo, J. J. (2007). Configuring an Evolutionary Tool for the Inventory and Transportation Problem. Paper presented at the GECCO’07, London, England, United Kingdom.

Hu, X., Shi, Y., & Eberhart, R. (2004). Recent Advances in Particle Swarm. Paper presented at the Congress on Evolutionary Computation, CEC2004.

Jiao, B., Lian, Z., & Gu, X. (2006). A Dynamic Inertia Weight Particle Swarm Optimization Algorithm. Chaos, Solitons and Fractals, 37(3), 698-705.

Kennedy, J., & Eberhart, R. C. (1995). Particle Swarm Optimization. Paper presented at the IEEE Conference on Neural Networks, Piscataway, NJ.

Kleywegt, A. J., Nori, V. S., & Savelsbergh, M. W. P. (2002 ). The Stochastic Inventory Routing Problem with Direct Deliveries. Transportation Science, 36(1), 94-118.

Lenstra, J. K., & Rinnooy, K. A. H. G. (1981). Complexity of vehicle routing and scheduling problems. Networks, 11, 221–227.

Li, J., Chen, H., & Chu, F. (2010). Performance Evaluation of Distribution Strategies for the Inventory Routing Problem. European Journal of Operational Research, 202 (2), 412-419.

Li, J.-A., Wu, Y., Lai, K. K., & Liu, K. (2008). Replenishment Routing Problems between a Single Supplier and Multiple Retailers with Direct Delivery. European Journal of Operational Research, 190(2), 412-420.

Moin, N. H., Salhi, S., & Aziz, N. A. B. (2010). An Efficient Hybrid Genetic Algorithm for the Multi-Product Multi-Period Inventory Routing Problem. International Journal of Production Economics, In Press, Corrected Proof.

Rusdiansyah, A., & Tsao, D.-b. (2005). An Integrated Model of the Periodic Delivery Problems for Vending-Machine Supply Chains. Journal of Food Engineering, 70(3), 421-434.

Vehicle routing data sets. from: http://www.coin-or.org/SYMPHONY/branchandcut/VRP/data.

Zhao, Q.-H., Chen, S., & Zang, C.-X. (2008). Model and Algorithm for Inventory/Routing Decision in a Three-Echelon Logistics System. European Journal of Operational Research, 191(3), 623-635.

Zhao, Q.-H., Wang, S.-Y., & Lai, K. K. (2007). A Partition Approach to the Inventory/Routing Problem. European Journal of Operational Research, 177(2), 786-802.