جوليا وسرب الجسيمات


نواصل دراسة طرق التحسين متعدد الأبعاد ، والخط التالي هو طريقة سرب الجسيمات ، التي تبحث عن حد أدنى عالمي.


النظرية



الخوارزمية بسيطة جدًا:



يتم حساب موضع كل جسيم في لحظة معينة بواسطة المعادلة:


Vi،t+1=AcVi،t+Cprp(pixi،t)+Cgrg(gixi،t)


xi،t+1=xi،t+Vi،t+1


أين pi- تنسيق أفضل حل لجسيم معين ، gi- تنسيق أفضل حل لجميع الجسيمات لهذا العصر ، Cpو Cg- عوامل الترجيح (مختارة لنموذج معين) ، Acهو معامل القصور الذاتي ، ويمكن جعله يعتمد على عدد الحقبة ، ثم سرعات الجسيمات سوف تتغير بسلاسة.


وظائف الاختبار


نظرًا لأنه من الجيد جدًا إلقاء نظرة على عمل الطريقة ، فسنحصل على المزيد من وظائف الاختبار:


مدسوس بعيدا
parabol(x) = sum(u->u*u, x) # f(0,0) = 0, x_i ∈ [-10,10] shvefel(x) = sum(u-> -u*sin(sqrt(abs(u))), x) # f(420.9687,420.9687) = -819?, x_i ∈ [-500,500] rastrigin(x) = 10*length(x) + sum(u->u*u-10*cos(2*pi*u), x) # f(0,0) = 0, x_i ∈ [-5,5] ekly(x) = -20exp(-0.2sqrt(0.5(x[1]*x[1]+x[2]*x[2]))) - exp(0.5(cospi(2x[1])+cospi(2x[2]))) + 20 + ℯ # f(0,0) = 0, x_i ∈ [-5,5] rosenbrok(x) = 100(x[2]-x[1]*x[1])^2 + (x[1]-1)^2 # f(0,0) = 0, x_i ∈ [-5,5] bill(x) = (1.5-x[1]+x[1]*x[2])^2 + (2.25-x[1]+x[1]*x[2]*x[2])^2 + (2.625-x[1]+x[1]*x[2]^3)^2 # f(3,0.5) = 0, x_i ∈ [-5,5] boot(x) = (x[1]+2x[2]-7)^2 + (2x[1]+x[2]-5)^2 # f(1,3) = 0, x_i ∈ [-10,10] bukin6(x) = 100sqrt(abs(x[2]-0.01x[1]*x[1])) + 0.01abs(x[1]+10) # f(-10,1) = 0, x_i ∈ [-15,-5; -3,3] levy13(x) = sinpi(3x[1])^2 + (1+sinpi(3x[2])^2)*(x[1]-1)^2 + (1+sinpi(2x[2])^2)*(x[2]-1)^2 # f(1,1) = 0, x_i ∈ [-10,10] himmelblau(x) = (x[1]^2+x[2]-11)^2 + (x[1]+x[2]^2-7)^2 # f(3,2)... = 0, x_i ∈ [-5,5] camel3humped(x) = 2x[1]^2 - 1.05x[1]^4 + x[1]^6 /6 + x[1]*x[2] + x[2]^2 # f(0,0) = 0, x_i ∈ [-5,5] izom(x) = -cos(x[1])*cos(x[2])*exp(-( (x[1]-pi)^2 + (x[2]-pi)^2 )) # f(π,π) = -1, x_i ∈ [-100,100] holdertable(x) = -abs(sin(x[1])*cos(x[2])exp(abs( 1-sqrt(x[1]^2+x[2]^2)/pi ))) # f(±8.05502,±9.66459) = -19.2085, x_i ∈ [-10,10] shaffer4(x) = 0.5 + (cos(sin(abs(x[1]^2-x[2]^2)))^2-0.5) / (1+0.001(x[1]^2+x[2]^2))^2 # f(0,1.25313) = 0.292579, x_i ∈ [-100,100] 

وفي الواقع ، فإن مركز موارد المهاجرين نفسه:


 function mdpso(; nparts = 50, ndimes = 2, ages = 50, #   lover = [-10 -10], upper = [10 10], C1 = [1.9 1.9], #  - C2 = [1.8 1.8], Ac = [0.1 0.1], ) minind = 0 V = zeros(nparts,ndimes) #   n  n X = zeros(nparts,ndimes) funmin = -Inf Fmin = Inf Fbest = fill(Fmin, nparts) funx = zeros(nparts) xmem = zeros(nparts,ndimes) xbest = zeros(ndimes) #   #      for i in 1:nparts, j in 1:ndimes X[i,j] = randomer(lover[j], upper[j]) end for i in 1:ages for j in 1:nparts funx[j] = fun(X[j,:]) if funx[j] < Fbest[j] Fbest[j] = funx[j] xmem[j,:] = X[j,:] end end #      ploter(lover, upper, X, funx, i); funmin = minimum(funx) minind = argmin(funx) if funmin < Fmin Fmin = funmin xbest[:] = X[minind,:] end for j in 1:nparts, k in 1:ndimes R1 = rand() R2 = rand() V[j,k] = Ac[k]*V[j,k] + C1[k]*R1*(xmem[j,k] - X[j,k]) + C2[k]*R2*(xbest[k] - X[j,k]) X[j,k] += V[j,k] end println("Age № $i\n xbest:\n $(xbest[1]) $(xbest[2])") println("Fmin: $Fmin\n") end f = open("$fun.txt","w") #      write(f,"C1 = $C1, C2 = $C2, Ac = $Ac, lower = $lover, upper = $upper, ages = $ages, parts = $nparts") close(f) end 

مع ذلك ، مع الحساب ، يستغرق الحساب وقتًا أطول ، لكنه أجمل:


 using Plots pyplot() function ploter(l, u, xy, z, n_age ) contour(Xs, Ys, Zs, fill = true); #    (   ) #   ,       scatter!(xy[:,1], xy[:,2], legend = false, xaxis=( (l[1], u[1])), yaxis=( (l[2], u[2])) ); #savefig("$fun $n_age.png") #        end fun = bill low = [-4 -4] up = [4 4] Xs = range(low[1], stop = up[1], length = 80) Ys = range(low[2], stop = up[2], length = 80) Zs = [ fun([xy]) for y in Ys, x in Xs ] surface(Xs, Ys, Zs) xaxis!( (low[1], up[1]), low[1]:(up[1]-low[1])/5:up[1] ) yaxis!( (low[2], up[2]), low[2]:(up[2]-low[2])/5:up[2] ) 


 mdpso(C1 = [1.2 1.2], C2 = [1.1 1.1], Ac = [0.08 0.08], lower = [-4 -4], upper = [4 4], ages = 30) 


 fun = ekly mdpso(C1 = [1.7 1.7], C2 = [1.7 1.7], Ac = [0.07 0.07], lower = [-5 -5], upper = [5 5], ages = 15) 



 fun = himmelblau mdpso(C1 = [1.1 1.1], C2 = [1.0 1.0], Ac = [0.09 0.09], lower = [-5 -5], upper = [5 5], ages = 20, parts = 50) 



 fun = holdertable mdpso(C1 = [1.1 1.1], C2 = [1.0 1.0], Ac = [0.09 0.09], lower = [-10 -10], upper = [10 10], ages = 20, parts = 50) 



 fun = levy13 mdpso(C1 = [1.1 1.1], C2 = [1.0 1.0], Ac = [0.09 0.09], lower = [-10 -10], upper = [10 10], ages = 20, parts = 50) 



 fun = shaffer4 mdpso(C1 = [1.1 1.1], C2 = [1.0 1.0], Ac = [0.09 0.09], lower = [-100 -100], upper = [100 100], ages = 20, parts = 50) 



إن ما يخفف حقًا هو التركيز على المعلمات وعنصر العشوائية: إن لم يكن الذباب يطير بالقرب من الحد الأدنى العالمي ، فإن هذا الأمر برمته قد يقع في محيط محلي:


 fun = rastrigin mdpso(mdpso(C1 = [0.1 0.1], C2 = [1 1], Ac = [0.08 0.08], lover = low, upper = up, ages = 30)) 



نعم وقديم ليس جدا جيد روزنبروك لا يزال لا يسمح لفهم نفسه:


 fun = rosenbrok mdpso(C1 = [1.7 1.7], C2 = [1.5 1.5], Ac = [0.15 0.15], lover = low, upper = up, ages = 20, nparts = 50) ... Age № 20 xbest: 0.37796421341886866 0.12799160066705667 Fmin: 0.409026370833564 



ولكن كما قلت سابقًا ، يمكنك استخدام FDM للبحث عن تقريب جيد للحد الأدنى العالمي ، ثم تحديد ، على سبيل المثال ، Nelder-Mead:


طريقة البسيط
 vecl(x) = sqrt( sum(u -> u*u, x) ) function sortcoord(Mx) N = size(Mx,2) f = [fun(Mx[:,i]) for i in 1:N] #     Mx[:, sortperm(f)] end function normx(Mx) m = size(Mx,2) D = zeros(m-1,m) for i = 1:m, j = i+1:m D[i,j] = vecl(Mx[:,i] - Mx[:,j]) #     end D sqrt(maximum(D)) end function ofNelderMid(; ndimes = 2, ε = 1e-4, fit = [.1, .1], low = [-1 -1], up = [1 1]) k = 0 N = ndimes Xx = zeros(N, N+1) for i = 1:N+1 Xx[:,i] = fit end for i = 1:N Xx[i,i] += 0.5*vecl(fit) + ε end p = normx(Xx) while p > ε k += 1 Xx = sortcoord(Xx) Xo = [ sum(Xx[i,1:N])/N for i = 1:N ] #  - i-  Ro = 2Xo - Xx[:,N+1] FR = fun(Ro) if FR > fun(Xx[:,N+1]) for i = 2:N+1 Xx[:,i] = 0.5(Xx[:,1] + Xx[:,i]) end else if FR < fun(Xx[:,1]) Eo = Xo + 2(Xo - Xx[:,N+1]) if FR > fun(Eo) Xx[:,N+1] = Eo else Xx[:,N+1] = Ro end else if FR <= fun(Xx[:,N]) Xx[:,N+1] = Ro else Co = Xo + 0.5(Xo - Xx[:,N+1]) if FR > fun(Co) Xx[:,N+1] = Co else Xx[:,N+1] = Ro end end end end println(k, " ", p, " ", Xx[:,1]) p = normx(Xx) end #while fit = Xx[:,1] end 

 ofNelderMid(fit = [0.37796 0.127992]) ... 92 0.00022610400555036366 [1.0, 1.0] 93 0.00015987967588703512 [1.0, 1.0] 94 0.00011305200343052599 [1.0, 1.0] 2-element Array{Float64,1}: 0.9999999996645973 0.9999999995466575 

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


الأخلاقية


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

هذا كل شيء لهذا اليوم ، شكرا لك على اهتمامك وأتمنى لكم كل التوفيق!

Source: https://habr.com/ru/post/ar440234/


All Articles