Tweet
' Simulated annealing (SA) demonstration: Traveling Salesman Problem
'
' See this link for an introduction:
' http://en.wikipedia.org/wiki/Simulated_annealing
'
' The traveling salesman problem is solved very effectively using
' this method. The following code is a Powerbasic translation of
' the Fortran code found in the following link, which explains in
' details how the method works:
' http://www.nrbook.com/a/bookfpdf/f10-9.pdf
'
' The present code is set up with 200 cities. You may well solve
' a problem with 500 cities in less than a minute. Amazing!
'
' Best regards,
'
' Erik Christensen ------------ April 24, 2008
'
' See this link for an introduction:
' http://en.wikipedia.org/wiki/Simulated_annealing
'
' Simulated annealing (SA) is a generic probabilistic meta-algorithm
' for the global optimization problem, namely locating a good
' approximation to the global optimum of a given function in a large
' search space. It is often used when the search space is discrete
' (e.g., all tours that visit a given set of cities). For certain
' problems, simulated annealing may be more effective than
' exhaustive enumeration - provided that the goal is merely to find
' an acceptably good solution in a fixed amount of time, rather than
' the best possible solution.
' The name and inspiration come from annealing in metallurgy, a
' technique involving heating and controlled cooling of a material
' to increase the size of its crystals and reduce their defects. The
' heat causes the atoms to become unstuck from their initial
' positions (a local minimum of the internal energy) and wander
' randomly through states of higher energy; the slow cooling gives
' them more chances of finding configurations with lower internal
' energy than the initial one.
' By analogy with this physical process, each step of the SA
' algorithm replaces the current solution by a random "nearby"
' solution, chosen with a probability that depends on the difference
' between the corresponding function values and on a global
' parameter T (called the temperature), that is gradually decreased
' during the process. The dependency is such that the current
' solution changes almost randomly when T is large, but increasingly
' "downhill" as T goes to zero. The allowance for "uphill" moves
' saves the method from becoming stuck at local minima-which are the
' bane of greedier methods.
' The method was independently described by S. Kirkpatrick, C. D.
' Gelatt and M. P. Vecchi in 1983, and by V. Cerný in 1985. The
' method is an adaptation of the Metropolis-Hastings algorithm, a
' Monte Carlo method to generate sample states of a thermodynamic
' system, invented by N. Metropolis et al in 1953.
'
' for the global optimization problem, namely locating a good
' approximation to the global optimum of a given function in a large
' search space. It is often used when the search space is discrete
' (e.g., all tours that visit a given set of cities). For certain
' problems, simulated annealing may be more effective than
' exhaustive enumeration - provided that the goal is merely to find
' an acceptably good solution in a fixed amount of time, rather than
' the best possible solution.
' The name and inspiration come from annealing in metallurgy, a
' technique involving heating and controlled cooling of a material
' to increase the size of its crystals and reduce their defects. The
' heat causes the atoms to become unstuck from their initial
' positions (a local minimum of the internal energy) and wander
' randomly through states of higher energy; the slow cooling gives
' them more chances of finding configurations with lower internal
' energy than the initial one.
' By analogy with this physical process, each step of the SA
' algorithm replaces the current solution by a random "nearby"
' solution, chosen with a probability that depends on the difference
' between the corresponding function values and on a global
' parameter T (called the temperature), that is gradually decreased
' during the process. The dependency is such that the current
' solution changes almost randomly when T is large, but increasingly
' "downhill" as T goes to zero. The allowance for "uphill" moves
' saves the method from becoming stuck at local minima-which are the
' bane of greedier methods.
' The method was independently described by S. Kirkpatrick, C. D.
' Gelatt and M. P. Vecchi in 1983, and by V. Cerný in 1985. The
' method is an adaptation of the Metropolis-Hastings algorithm, a
' Monte Carlo method to generate sample states of a thermodynamic
' system, invented by N. Metropolis et al in 1953.
'
' this method. The following code is a Powerbasic translation of
' the Fortran code found in the following link, which explains in
' details how the method works:
' http://www.nrbook.com/a/bookfpdf/f10-9.pdf
'
' The present code is set up with 200 cities. You may well solve
' a problem with 500 cities in less than a minute. Amazing!
'
' Best regards,
'
' Erik Christensen ------------ April 24, 2008
Code:
' Simulated annealing (SA) demonstration: Traveling Salesman Problem
'
' See this link for an introduction:
' http://en.wikipedia.org/wiki/Simulated_annealing
' quote:
' "Simulated annealing (SA) is a generic probabilistic meta-algorithm
' for the global optimization problem, namely locating a good
' approximation to the global optimum of a given function in a large
' search space. It is often used when the search space is discrete
' (e.g., all tours that visit a given set of cities). For certain
' problems, simulated annealing may be more effective than
' exhaustive enumeration - provided that the goal is merely to find
' an acceptably good solution in a fixed amount of time, rather than
' the best possible solution.
' The name and inspiration come from annealing in metallurgy, a
' technique involving heating and controlled cooling of a material
' to increase the size of its crystals and reduce their defects. The
' heat causes the atoms to become unstuck from their initial
' positions (a local minimum of the internal energy) and wander
' randomly through states of higher energy; the slow cooling gives
' them more chances of finding configurations with lower internal
' energy than the initial one.
' By analogy with this physical process, each step of the SA
' algorithm replaces the current solution by a random "nearby"
' solution, chosen with a probability that depends on the difference
' between the corresponding function values and on a global
' parameter T (called the temperature), that is gradually decreased
' during the process. The dependency is such that the current
' solution changes almost randomly when T is large, but increasingly
' "downhill" as T goes to zero. The allowance for "uphill" moves
' saves the method from becoming stuck at local minima-which are the
' bane of greedier methods.
' The method was independently described by S. Kirkpatrick, C. D.
' Gelatt and M. P. Vecchi in 1983, and by V. Cerný in 1985. The
' method is an adaptation of the Metropolis-Hastings algorithm, a
' Monte Carlo method to generate sample states of a thermodynamic
' system, invented by N. Metropolis et al in 1953."
'
' The traveling salesman problem is solved very effectively using
' this method. The following code is a Powerbasic translation of
' the Fortran code found in the following link, which explains in
' details how the method works:
' http://www.nrbook.com/a/bookfpdf/f10-9.pdf
'
' The present code is set up with 200 cities. You may well solve
' a problem with 500 cities in less than a minute. Amazing!
'
' Best regards,
'
' Erik Christensen ------------ April 24, 2008
#COMPILE EXE
#DIM ALL
#INCLUDE "WIN32API.INC"
#INCLUDE "COMMCTRL.INC"
#INCLUDE "PBForms.INC"
'
%IDC_BUTTON_EXIT = 1002
%IDC_BUTTON_RUN = 1001
%IDC_GRAPHIC1 = 1003
%IDC_TEXTBOX1 = 1004
%IDD_DIALOG1 = 101
'
GLOBAL s1 AS STRING
'
FUNCTION ALEN (BYVAL X1 AS SINGLE, BYVAL X2 AS SINGLE, BYVAL Y1 AS SINGLE, BYVAL Y2 AS SINGLE) AS SINGLE
FUNCTION = SQR((X2 - X1) ^ 2 + (Y2 - Y1) ^ 2)
END FUNCTION
'
SUB ANNEAL (BYREF X() AS SINGLE, BYREF Y() AS SINGLE, BYREF IORDER() AS LONG, BYVAL NCITY AS LONG, BYVAL hDlg AS LONG)
DIM N(6) AS LONG
LOCAL NOVER AS LONG, NLIMIT AS LONG, TFACTR AS SINGLE, PATH AS SINGLE, T AS SINGLE
LOCAL I AS LONG, J AS LONG, I1 AS LONG, I2 AS LONG, K AS LONG, NN AS LONG, II AS LONG
LOCAL NSUCC AS LONG, IDEC AS LONG, DE AS SINGLE, ANS AS LONG, lResult AS LONG, lResult2 AS LONG
NOVER = 100 * NCITY
NLIMIT = 10 * NCITY
TFACTR = .9
PATH = 0!
T = .5
FOR I = 1 TO NCITY - 1
I1 = IORDER(I)
I2 = IORDER(I + 1)
PATH = PATH + ALEN(X(I1), X(I2), Y(I1), Y(I2))
NEXT I
I1 = IORDER(NCITY)
I2 = IORDER(1)
PATH = PATH + ALEN(X(I1), X(I2), Y(I1), Y(I2))
FOR J = 1 TO 100
NSUCC = 0
FOR K = 1 TO NOVER
DO
N(1) = 1 + INT(NCITY * RND)
N(2) = 1 + INT((NCITY - 1) * RND)
IF N(2) >= N(1) THEN N(2) = N(2) + 1
IDEC = RND(0, 1)
NN = 1 + (N(1) - N(2) + NCITY - 1) MOD NCITY
LOOP WHILE NN < 3
IF IDEC = 0 THEN
N(3) = N(2) + INT(ABS(NN - 2) * RND) + 1
N(3) = 1 + N(3) - 1 - NCITY * INT((N(3) - 1) / NCITY)
CALL TRNCST(X(), Y(), IORDER(), NCITY, N(), DE)
CALL METROP(DE, T, ANS)
IF ANS THEN
NSUCC = NSUCC + 1
PATH = PATH + DE
CALL TRNSPT(IORDER(), NCITY, N())
END IF
ELSE
CALL REVCST(X(), Y(), IORDER(), NCITY, N(), DE)
CALL METROP(DE, T, ANS)
IF ANS THEN
NSUCC = NSUCC + 1
PATH = PATH + DE
CALL REVERS(IORDER(), NCITY, N())
END IF
END IF
IF NSUCC >= NLIMIT THEN EXIT FOR
NEXT K
GRAPHIC CLEAR %WHITE
GRAPHIC SET POS (X(IORDER(1)), Y(IORDER(1)))
FOR I = 1 TO NCITY
IF I < NCITY THEN I2 = IORDER(I+1) ELSE I2 = IORDER(1)
GRAPHIC LINE - (X(I2), Y(I2)), %RED
GRAPHIC ELLIPSE (X(I2-.1), Y(I2)-.06) - (X(I2)+.1, Y(I2)+.06), %RED, %RED
NEXT I
s1 = s1 + "T="+STR$(ROUND(T,3))+$TAB+" Length ="+STR$(ROUND(PATH,3))+$TAB+" Moves:"+STR$(NSUCC)+$CRLF
CONTROL SET TEXT hDlg, %IDC_TEXTBOX1, s1
CONTROL SEND hDlg, %IDC_TEXTBOX1, %EM_GETLINECOUNT, 0, 0 TO lResult
CONTROL SEND hDlg, %IDC_TEXTBOX1, %EM_GETFIRSTVISIBLELINE, 0, 0 TO lResult2
' scroll down to last calculation
CONTROL SEND hDlg, %IDC_TEXTBOX1, %EM_LINESCROLL, 0, lResult - lResult2 - 9
T = T * TFACTR
IF NSUCC = 0 THEN EXIT FOR
NEXT J
END SUB
'
SUB METROP (BYREF DE AS SINGLE, BYREF T AS SINGLE, BYREF ANS AS LONG) ' Metropolis algorithm
ANS = (DE < 0!) OR (RND < EXP(-DE / T))
END SUB
'
SUB REVCST (BYREF X() AS SINGLE, BYREF Y() AS SINGLE, BYREF IORDER() AS LONG, BYVAL NCITY AS LONG, BYREF N() AS LONG, BYREF DE AS SINGLE)
DIM XX(4) AS SINGLE, YY(4) AS SINGLE, J AS LONG, II AS LONG
N(3) = 1 + N(1) + NCITY - 2 - NCITY * INT((N(1) + NCITY - 2) / NCITY)
N(4) = 1 + N(2) - NCITY * INT(N(2) / NCITY)
FOR J = 1 TO 4
II = IORDER(N(J))
XX(J) = X(II)
YY(J) = Y(II)
NEXT J
DE = -ALEN(XX(1), XX(3), YY(1), YY(3))
DE = DE - ALEN(XX(2), XX(4), YY(2), YY(4))
DE = DE + ALEN(XX(1), XX(4), YY(1), YY(4))
DE = DE + ALEN(XX(2), XX(3), YY(2), YY(3))
END SUB
'
SUB REVERS (BYREF IORDER() AS LONG, BYVAL NCITY AS LONG, BYREF N() AS LONG)
LOCAL NN AS LONG, J AS LONG, K AS LONG, L AS LONG
NN = (1 + (N(2) - N(1) + NCITY) MOD NCITY) / 2
FOR J = 1 TO NN
K = 1 + N(1) + J - 2 - NCITY * INT((N(1) + J - 2) / NCITY)
L = 1 + (N(2) - J + NCITY) MOD NCITY
SWAP IORDER(K), IORDER(L)
NEXT J
END SUB
'
SUB TRNCST (BYREF X() AS SINGLE, BYREF Y() AS SINGLE, BYREF IORDER() AS LONG, BYVAL NCITY AS LONG, BYREF N() AS LONG, BYREF DE AS SINGLE)
DIM XX(6) AS SINGLE, YY(6) AS SINGLE, J AS LONG, II AS LONG
N(4) = 1 + N(3) - NCITY * INT(N(3) / NCITY)
N(5) = 1 + N(1) + NCITY - 2 - NCITY * INT((N(1) + NCITY - 2) / NCITY)
N(6) = 1 + N(2) - NCITY * INT(N(2) / NCITY)
FOR J = 1 TO 6
II = IORDER(N(J))
XX(J) = X(II)
YY(J) = Y(II)
NEXT J
DE = -ALEN(XX(2), XX(6), YY(2), YY(6))
DE = DE - ALEN(XX(1), XX(5), YY(1), YY(5))
DE = DE - ALEN(XX(3), XX(4), YY(3), YY(4))
DE = DE + ALEN(XX(1), XX(3), YY(1), YY(3))
DE = DE + ALEN(XX(2), XX(4), YY(2), YY(4))
DE = DE + ALEN(XX(5), XX(6), YY(5), YY(6))
END SUB
'
SUB TRNSPT (BYREF IORDER() AS LONG, BYVAL NCITY AS LONG, BYREF N() AS LONG)
DIM JORDER(NCITY) AS LONG, M1 AS LONG, M2 AS LONG, M3 AS LONG, NN AS LONG
LOCAL J AS LONG, JJ AS LONG
M1 = 1 + N(2) - N(1) + NCITY - NCITY * INT((N(2) - N(1) + NCITY) / NCITY)
M2 = 1 + N(5) - N(4) + NCITY - NCITY * INT((N(5) - N(4) + NCITY) / NCITY)
M3 = 1 + N(3) - N(6) + NCITY - NCITY * INT((N(3) - N(6) + NCITY) / NCITY)
NN = 1
FOR J = 1 TO M1
JJ = 1 + J + N(1) - 2 - NCITY * INT((J + N(1) - 2) / NCITY)
JORDER(NN) = IORDER(JJ)
INCR NN
NEXT J
IF M2 > 0 THEN
FOR J = 1 TO M2
JJ = 1 + J + N(4) - 2 - NCITY * INT((J + N(4) - 2) / NCITY)
JORDER(NN) = IORDER(JJ)
INCR NN
NEXT J
END IF
IF M3 > 0 THEN
FOR J = 1 TO M3
JJ = 1 + J + N(6) - 2 - NCITY * INT((J + N(6) - 2) / NCITY)
JORDER(NN) = IORDER(JJ)
INCR NN
NEXT J
END IF
FOR J = 1 TO NCITY
IORDER(J) = JORDER(J)
NEXT J
END SUB
'
CALLBACK FUNCTION ShowDIALOG1Proc()
LOCAL st AS STRING, lResult AS LONG, lResult2 AS LONG
STATIC NCITY AS LONG, I AS LONG, II AS LONG, I2 AS LONG, hFont AS LONG
SELECT CASE AS LONG CBMSG
CASE %WM_INITDIALOG
' Initialization handler
hFont = CreateFont(16,0,0,0,400,0,0,0,0,3,2,1,82,"Ariel")
CONTROL SEND CBHNDL, %IDC_TEXTBOX1, %WM_SETFONT, hFont, %TRUE
NCITY = 200 ' number of cities
DIM X(NCITY) AS STATIC SINGLE, Y(NCITY) AS STATIC SINGLE, IORDER(NCITY) AS STATIC LONG
DIALOG SET TEXT CBHNDL, "Traveling Salesman Problem solved by Simulated Annealing:"+STR$(NCITY)+" Cities"
GRAPHIC ATTACH CBHNDL, %IDC_GRAPHIC1
GRAPHIC CLEAR %WHITE
CASE %WM_NCACTIVATE
STATIC hWndSaveFocus AS DWORD
IF ISFALSE CBWPARAM THEN
' Save control focus
hWndSaveFocus = GetFocus()
ELSEIF hWndSaveFocus THEN
' Restore control focus
SetFocus(hWndSaveFocus)
hWndSaveFocus = 0
END IF
CASE %WM_COMMAND
' Process control notifications
SELECT CASE AS LONG CBCTL
CASE %IDC_BUTTON_RUN
IF CBCTLMSG = %BN_CLICKED OR CBCTLMSG = 1 THEN
s1 = ""
RANDOMIZE TIMER
' Create points of sale: x and y coordinate of city
FOR I = 1 TO NCITY
X(I) = 0.1+RND*9.8
Y(I) = 0.1+RND*9.8
IORDER(I) = I
NEXT I
GRAPHIC SCALE (0, 10) - (10, 0)
'
CALL ANNEAL(X(), Y(), IORDER(), NCITY, CBHNDL)
'
s1 = s1+$CRLF+"********** Finished **********"+$CRLF+ "Final path:"+$CRLF+"City"+$TAB+" x"+$TAB+" y"+$CRLF
CONTROL SET TEXT CBHNDL, %IDC_TEXTBOX1, s1
CONTROL SEND CBHNDL, %IDC_TEXTBOX1, %EM_GETLINECOUNT, 0, 0 TO lResult
CONTROL SEND CBHNDL, %IDC_TEXTBOX1, %EM_GETFIRSTVISIBLELINE, 0, 0 TO lResult2
FOR I = 1 TO NCITY
II = IORDER(I)
s1 = s1 + STR$(II) + $TAB + STR$(ROUND(X(II),3))+ $TAB + STR$(ROUND(Y(II),3))+$CRLF
NEXT I
CONTROL SET TEXT CBHNDL, %IDC_TEXTBOX1, s1
' scroll down to "Finished" and "Final Path"
CONTROL SEND CBHNDL, %IDC_TEXTBOX1, %EM_LINESCROLL, 0, lResult - lResult2 - 22
END IF
CASE %IDC_BUTTON_EXIT
IF CBCTLMSG = %BN_CLICKED OR CBCTLMSG = 1 THEN DIALOG END CBHNDL
CASE %WM_DESTROY
DeleteObject hFont
PostQuitMessage 0
FUNCTION = 0 : EXIT FUNCTION
END SELECT
END SELECT
END FUNCTION
'
FUNCTION ShowDIALOG1(BYVAL hParent AS DWORD) AS LONG
LOCAL lRslt AS LONG
LOCAL hDlg AS DWORD
DIALOG NEW PIXELS, hParent, "Traveling Salesman Problem solved by " + _
"Simulated Annealing", , , 690, 421, %WS_POPUP OR %WS_BORDER _
OR %WS_DLGFRAME OR %WS_SYSMENU OR %WS_MINIMIZEBOX OR _
%WS_CLIPSIBLINGS OR %WS_VISIBLE OR %DS_MODALFRAME OR %DS_3DLOOK OR _
%DS_NOFAILCREATE OR %DS_SETFONT, %WS_EX_CONTROLPARENT OR %WS_EX_LEFT _
OR %WS_EX_LTRREADING OR %WS_EX_RIGHTSCROLLBAR, TO hDlg
CONTROL ADD BUTTON, hDlg, %IDC_BUTTON_RUN, "&Run ", 276, 390, 162, 26, _
%WS_CHILD OR %WS_VISIBLE OR %WS_TABSTOP OR %BS_TEXT OR _
%BS_PUSHBUTTON OR %BS_CENTER OR %BS_VCENTER, %WS_EX_STATICEDGE OR _
%WS_EX_LEFT OR %WS_EX_LTRREADING
CONTROL ADD BUTTON, hDlg, %IDC_BUTTON_EXIT, "E&xit", 612, 390, 72, 26, _
%WS_CHILD OR %WS_VISIBLE OR %WS_TABSTOP OR %BS_TEXT OR _
%BS_PUSHBUTTON OR %BS_CENTER OR %BS_VCENTER, %WS_EX_STATICEDGE OR _
%WS_EX_LEFT OR %WS_EX_LTRREADING
CONTROL ADD GRAPHIC, hDlg, %IDC_GRAPHIC1, "", 6, 6, 384, 378, %WS_CHILD _
OR %WS_VISIBLE, %WS_EX_DLGMODALFRAME
CONTROL ADD TEXTBOX, hDlg, %IDC_TEXTBOX1, "", 396, 6, 288, 378, %WS_CHILD _
OR %WS_VISIBLE OR %WS_TABSTOP OR %WS_HSCROLL OR %WS_VSCROLL OR _
%ES_LEFT OR %ES_MULTILINE, %WS_EX_CLIENTEDGE OR %WS_EX_LEFT OR _
%WS_EX_LTRREADING OR %WS_EX_RIGHTSCROLLBAR
DIALOG SHOW MODAL hDlg, CALL ShowDIALOG1Proc TO lRslt
FUNCTION = lRslt
END FUNCTION
'
FUNCTION PBMAIN()
PBFormsInitComCtls (%ICC_WIN95_CLASSES OR %ICC_DATE_CLASSES OR _
%ICC_INTERNET_CLASSES)
ShowDIALOG1 %HWND_DESKTOP
END FUNCTION
Comment