Particle swarm optimization: Difference between revisions

Content added Content deleted
m (→‎{{header|REXX}}: added/changed whitespace and comments, added a comment to the REXX section header.)
Line 416: Line 416:
Classic REXX doesn't have a   '''sine'''   function, so a RYO version is included here.
Classic REXX doesn't have a   '''sine'''   function, so a RYO version is included here.


The numeric precision is only limited to the number of decimal digits defined in the &nbsp; <big> '''pi''' </big> &nbsp; variable &nbsp; (in this case, '''100''').
The numeric precision is only limited to the number of decimal digits defined in the &nbsp; <big> '''pi''' </big> &nbsp; variable &nbsp; (in this case, &nbsp; '''100''').


This REXX version supports the specifying of '''X''', '''Y''', and '''D''', &nbsp; as well as the number of particles, and the number of decimal digits to be displayed.
This REXX version supports the specifying of &nbsp; '''X''', &nbsp; '''Y''', &nbsp; and &nbsp; '''D''', &nbsp; as well as the number of particles, &nbsp; and the number of decimal digits to be displayed.


The refinement loop is stopped when the function value stabilizes.
The refinement loop is stopped when the function value stabilizes.

<lang rexx>/*REXX pgm calc. Particle Swarm Optimization as it migrates through a solution*/
Note that REXX used decimal floating point, not binary.
numeric digits length(pi()) /*sDigs: is the # of displayed digits.*/
<lang rexx>/*REXX program calculates Particle Swarm Optimization as it migrates through a solution.*/
parse arg x y d #part sDigs . /*obtain optional arguments from the CL*/
if x=='' | x==',' then x= -0.5 /*is X not defined?*/
numeric digits length(pi()) /*sDigs: is the # of displayed digits.*/
if y=='' | y==',' then y= -1.5 /* " Y " " */
parse arg x y d #part sDigs . /*obtain optional arguments from the CL*/
if d=='' | d==',' then d= 1 /* " D " " */
if x=='' | x==',' then x= -0.5 /*Not specified? Then use the default.*/
if #part=='' | #part==',' then #part=1e12 /* " the # particles " " */
if y=='' | y==',' then y= -1.5 /* " " " " " " */
if sDigs=='' | sDigs==',' then sDigs= 25 /* " the # of digits " " */
if d=='' | d==',' then d= 1 /* " " " " " " */
minF=#part; old= /*number of particles is one billion. */
if #part=='' | #part==',' then #part=1e12 /* " " " " " " */
if sDigs=='' | sDigs==',' then sDigs= 25 /* " " " " " " */
say center('X',sDigs+3,'═') center('Y',sDigs+3,'═') center('D',sDigs+3,'═')
minF=#part; old= /*number of particles is one billion. */
say center('X', sDigs+3, '═') center('Y', sDigs+3, '═') center('D', sDigs+3, '═')
call refine x,y
call refine x,y
do until refine(minX,minY); d=d*.2
do until refine(minX,minY); d=d*.2 /*increase the difference.*/
end /*until ···*/ /* [↑] stop refining if no difference.*/
end /*until ···*/ /* [↑] stop refining if no difference.*/
say
say
indent=1 + (sDigs+3) * 2 /*compute the indentation for alignment*/
indent=1 + (sDigs+3) * 2 /*compute the indentation for alignment*/
say right('The global minimum for f(-.54719, -1.54719) ───► ', indent) fmt(f(-.54719, -1.54719))
say right('The global minimum for f(-.54719, -1.54719) ───► ', indent) fmt(f(-.54719, -1.54719))
say right('The published global minimum is:' , indent) fmt( -1.9133 )
say right('The published global minimum is:' , indent) fmt( -1.9133 )
exit /*stick a fork in it, we're all done. */
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
/*────────────────────────────────────────────────────────────────────────────*/
refine: parse arg xx,yy; h=d * .5
refine: parse arg xx,yy; h=d * .5
do x=xx-d to xx+d by h
do x=xx-d to xx+d by h
do y=yy-d to yy+d by h; f=f(x,y); if f>=minF then iterate
do y=yy-d to yy+d by h; f=f(x,y); if f>=minF then iterate
new=fmt(x) fmt(y) fmt(f); if new=old then return 1
new=fmt(x) fmt(y) fmt(f); if new=old then return 1
say new; minF=f; minX=x; minY=y; old=new
say new; minF=f; minX=x; minY=y; old=new
end /*y*/
end /*y*/
end /*x*/
end /*x*/
return 0
return 0
/*──────────────────────────────────────────────────────────────────────────────────one─liner subroutines───────────────────────────────*/
/*──────────────────────────────────────────────────────────────────────────────────one─liner subroutines───────────────────────────────*/
f: procedure: parse arg a,b; return sin(a+b) + (a-b)**2 - 1.5*a + 2.5*b + 1
f: procedure: parse arg a,b; return sin(a+b) + (a-b)**2 - 1.5*a + 2.5*b + 1
fmt: parse arg ?; ?=format(?,,sDigs); L=length(?); if pos(.,?)\==0 then ?=strip(strip(?,'T',0),'T',.); return left(?,L)
fmt: parse arg ?; ?=format(?,,sDigs); L=length(?); if pos(.,?)\==0 then ?=strip(strip(?,'T',0),'T',.); return left(?,L)
pi: pi=3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068; return pi
pi: pi=3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068; return pi
r2r: return arg(1) // (pi()*2) /*normalize radians ───► a unit circle.*/
r2r: return arg(1) // (pi()*2) /*normalize radians ───► a unit circle.*/
sin: procedure; parse arg x; x=r2r(x); numeric fuzz 5; z=x; _=x; q=x*x; do k=2 by 2 until p=z; p=z; _=-_*q/(k*(k+1)); z=z+_; end; return z</lang>
sin: procedure; parse arg x; x=r2r(x); numeric fuzz 5; z=x; _=x; q=x*x; do k=2 by 2 until p=z; p=z; _=-_*q/(k*(k+1)); z=z+_; end; return z</lang>
'''output''' &nbsp; when using the default inputs:
'''output''' &nbsp; when using the default inputs: