diff --git a/NEWS b/NEWS
index 69ca1be7ba8681763d2ecc077cfc018de630856a..158d6a0d9c0b5880e67431195df92aef9bec4184 100755
--- a/NEWS
+++ b/NEWS
@@ -1,8 +1,14 @@
NEWS/ChangeLog for hydroPSO
--------------------------
0.1-59 (under-development & not working yet)
- o 'hydroPSO' : -) towards SPSO 2011 capable...
+ o 'hydroPSO' : -) towards SPSO 2011 capable, which includes:
+ * new default value for swarm size (npart=40),
+ * new equation for updating the velocity ( V[t+1] = w*V + Gr -x + r*alea_sphere(Gr, ||Gr-x||) ),
+ * new initialization of the velocity ( Vini = U(lower-Xini, upper-Xini) )
+ * optional normalisation of parameter values
+ -) new 'normalise' parameter for the 'control' variable, in order to improve the performance when the search space is not
+ a hypercube
-) argument 'method' now allows the following 2 new values: 'spso2007' and 'spso2011', which set the values of the PSO
engine to the ones of the SPSO 2007 and SPSO 2011, respectively.
-) argument 'method'. The old (default) value 'pso' was replaced by the new (default) value 'spso2011'
@@ -16,8 +22,12 @@ NEWS/ChangeLog for hydroPSO
method == 'spso2007' => npart=10+2*[sqrt(n)]
method != 'spso2007' => npart=40
-) in the documentation, default values are now mentioned for each argument of 'control'.
+ -) the 'PSO_logfile.txt' now contains information about: 'best.update', 'random.update', 'Xini.type', 'Vini.type'.
+ In addition, information about 'IW.type' and 'IW.exp' is now written only when 'length(IW.w) > 1'
- o Results obtained when running 'hydroPSO' >= 0.1-59 with default settings WILL BE DIFFERENT from those obtained with 'hydroPSO' <= 0.1-58, due to the following changes in default values:
+ o Running 'hydroPSO' >= 0.1-59 with default settings will produce DIFFERENT RESULTS from those obtained with 'hydroPSO' <= 0.1-58, due to the
+ following changes in default values:
+ ver <= 0.1-58 -> ver >= 0.1-59
-) npart : 10+2*[sqrt(n)] -> 40
-) Vini.type : 'lhs2007' -> 'random2011'
-) boundary.wall: 'reflecting' -> 'absorbing'
diff --git a/R/PSO_v2012.R b/R/PSO_v2012.R
index 938409de94e4a6f1a063e752ee74a7caeb1b9674..757ac78a2f92a9273db479041c6da96df2a1e1ee 100755
--- a/R/PSO_v2012.R
+++ b/R/PSO_v2012.R
@@ -537,7 +537,7 @@ position.update.and.boundary.treatment <- function(x, v, x.MinMax, boundary.wall
if ( length(byd.min.pos) > 0) {
if ( boundary.wall == "absorbing") {
x.new[byd.min.pos] <- x.min[byd.min.pos]
- v.new[byd.min.pos] <- 0*v[byd.min.pos]
+ v.new[byd.min.pos] <- -0.0*v[byd.min.pos]
} else if ( boundary.wall == "reflecting") {
x.new[byd.min.pos] <- 2*x.min[byd.min.pos] - x.new[byd.min.pos]
v.new[byd.min.pos] <- v[byd.min.pos]
@@ -555,7 +555,7 @@ position.update.and.boundary.treatment <- function(x, v, x.MinMax, boundary.wall
if ( length(byd.max.pos) > 0 ) {
if ( boundary.wall == "absorbing") {
x.new[byd.max.pos] <- x.max[byd.max.pos]
- v.new[byd.max.pos] <- 0*v[byd.max.pos]
+ v.new[byd.max.pos] <- -0.0*v[byd.max.pos]
} else if ( boundary.wall == "reflecting") {
x.new[byd.max.pos] <- 2*x.max[byd.max.pos] - x.new[byd.max.pos]
v.new[byd.max.pos] <- v[byd.max.pos]
@@ -1559,6 +1559,7 @@ hydroPSO <- function(
maxfn=Inf,
c1= 0.5+log(2),
c2= 0.5+log(2),
+ use.IW= TRUE, IW.w=1/(2*log(2)), IW.type=c("linear", "non-linear", "runif", "aiwf", "GLratio"), IW.exp= 1,
use.CF= FALSE,
lambda= 1,
@@ -1572,8 +1573,9 @@ hydroPSO <- function(
topology=c("random", "gbest", "lbest", "vonNeumann"), K=3,
iter.ini=0, # only used when 'topology=lbest'
ngbest=4, # only used when 'method=ipso'
+
+ normalise=FALSE,
- use.IW = TRUE, IW.w=1/(2*log(2)), IW.type=c("linear", "non-linear", "runif", "aiwf", "GLratio"), IW.exp= 1,
use.TVc1= FALSE, TVc1.rng= c(1.28, 1.05), TVc1.type=c("linear", "non-linear", "GLratio"), TVc1.exp= 1.5,
use.TVc2= FALSE, TVc2.rng= c(1.05, 1.28), TVc2.type=c("linear", "non-linear"), TVc2.exp= 1.5,
use.TVlambda=FALSE, TVlambda.rng= c(1, 0.25), TVlambda.type=c("linear", "non-linear"), TVlambda.exp= 1,
@@ -1623,6 +1625,9 @@ hydroPSO <- function(
maxfn <- con[["maxfn"]]
c1 <- con[["c1"]]
c2 <- con[["c2"]]
+ use.IW <- as.logical(con[["use.IW"]])
+ IW.w <- con[["IW.w"]]
+ IW.exp <- con[["IW.exp"]]
use.CF <- con[["use.CF"]]
lambda <- con[["lambda"]]
abstol <- con[["abstol"]]
@@ -1630,10 +1635,8 @@ hydroPSO <- function(
random.update <- as.logical(con[["random.update"]])
K <- con[["K"]]
iter.ini <- con[["iter.ini"]]
- ngbest <- con[["ngbest"]]
- use.IW <- as.logical(con[["use.IW"]])
- IW.w <- con[["IW.w"]]
- IW.exp <- con[["IW.exp"]]
+ ngbest <- con[["ngbest"]]
+ normalise <- as.logical(con[["normalise"]])
use.TVc1 <- as.logical(con[["use.TVc1"]])
TVc1.rng <- con[["TVc1.rng"]]
TVc1.exp <- con[["TVc1.exp"]]
@@ -1744,6 +1747,14 @@ hydroPSO <- function(
} # ELSE end
n <- nrow(X.Boundaries)
+
+ if (normalise) {
+ X.Boundaries[,1] <- rep(0, n)
+ X.Boundaries[,2] <- rep(1, n)
+
+ lower.mat <- matrix( rep(X.Boundaries[,1], npart), nrow=npart, byrow=TRUE)
+ upper.mat <- matrix( rep(X.Boundaries[,2], npart), nrow=npart, byrow=TRUE)
+ } # IF end
if (is.null(rownames(X.Boundaries))) {
param.IDs <- paste("Param", 1:n, sep="")
@@ -1890,7 +1901,7 @@ hydroPSO <- function(
} else if ( (class(par)=="matrix") | (class(par)=="data.frame") ) {
tmp <- ncol(par)
if ( tmp != n )
- stop( paste("Invalid argument: 'ncol(par) != n' (",tmp, "!=", n, ")", sep="") )
+ stop( "Invalid argument: 'ncol(par) != n' (",tmp, "!=", n, ")" )
tmp <- nrow(par)
X[1:tmp,] <- par
} # ELSE end
@@ -2034,6 +2045,18 @@ hydroPSO <- function(
writeLines(c("c2 :", c2), PSOparam.TextFile, sep=" ")
writeLines("", PSOparam.TextFile)
} # ELSE end
+ writeLines(c("use.IW :", use.IW), PSOparam.TextFile, sep=" ")
+ writeLines("", PSOparam.TextFile)
+ if (use.IW) {
+ writeLines(c("IW.w :", IW.w), PSOparam.TextFile, sep=" ")
+ writeLines("", PSOparam.TextFile)
+ if ( length(IW.w) > 1 ) {
+ writeLines(c("IW.type :", IW.type), PSOparam.TextFile, sep=" ")
+ writeLines("", PSOparam.TextFile)
+ writeLines(c("IW.exp :", IW.exp), PSOparam.TextFile, sep=" ")
+ writeLines("", PSOparam.TextFile)
+ } # IF end
+ } # IF end
if (use.TVlambda) {
writeLines(c("use.TVlambda :", use.TVlambda), PSOparam.TextFile, sep=" ")
writeLines("", PSOparam.TextFile)
@@ -2047,16 +2070,6 @@ hydroPSO <- function(
writeLines(c("lambda :", lambda), PSOparam.TextFile, sep=" ")
writeLines("", PSOparam.TextFile)
} # ELSE end
- if (use.IW) {
- writeLines(c("use.IW :", use.IW), PSOparam.TextFile, sep=" ")
- writeLines("", PSOparam.TextFile)
- writeLines(c("IW.type :", IW.type), PSOparam.TextFile, sep=" ")
- writeLines("", PSOparam.TextFile)
- writeLines(c("IW.w :", IW.w), PSOparam.TextFile, sep=" ")
- writeLines("", PSOparam.TextFile)
- writeLines(c("IW.exp :", IW.exp), PSOparam.TextFile, sep=" ")
- writeLines("", PSOparam.TextFile)
- } # IF end
writeLines(c("maxfn :", maxfn), PSOparam.TextFile, sep=" ")
writeLines("", PSOparam.TextFile)
writeLines(c("abstol :", abstol), PSOparam.TextFile, sep=" ")
@@ -2284,9 +2297,12 @@ hydroPSO <- function(
} # IF end
##########################################################################
+
+ if (normalise) X <- X * (upper.mat - lower.mat) + lower.mat
+
# 3.a) Evaluate the particles fitness
if ( fn.name != "hydromod" ) {
-
+
# Evaluating an R Function
GoF <- apply(X, fn, MARGIN=1, ...)
@@ -2322,7 +2338,6 @@ hydroPSO <- function(
} # ELSE end
-
if ( best.update == "sync" ) {
tmp <- sync.update.pgbests(x=X,
xt.fitness= Xt.fitness[iter, ],
@@ -2574,7 +2589,7 @@ hydroPSO <- function(
gbest.fit.bak <- ngbest.fit
} # IF end
- if (verbose) message("[Re-grouping particles in the swarm (iter: ", iter, ") ...]")
+ if (verbose) message("[ Re-grouping particles in the swarm (iter: ", iter, ") ... ]")
tmp <- RegroupingSwarm(x=X,
gbest= X.best.part[gbest.pos, ],
@@ -2665,13 +2680,17 @@ hydroPSO <- function(
# File 'BestParamPerIter.txt' #
GoF <- gbest.fit
if(is.finite(GoF)) {
+
+ ifelse(normalise, temp <- X.best.part[gbest.pos, ] * (upper - lower) + lower,
+ temp <- X.best.part[gbest.pos, ] )
+
writeLines( as.character( c(iter,
formatC(GoF, format="E", digits=digits, flag=" "),
- formatC(X.best.part[gbest.pos, ], format="E", digits=digits, flag=" ")
+ formatC(temp, format="E", digits=digits, flag=" ")
) ), BestParamPerIter.TextFile, sep=" ")
} else writeLines( as.character( c(iter,
"NA",
- formatC(X.best.part[gbest.pos, ], format="E", digits=digits, flag=" ")
+ formatC(temp, format="E", digits=digits, flag=" ")
) ), BestParamPerIter.TextFile, sep=" ")
writeLines("", BestParamPerIter.TextFile)
flush(BestParamPerIter.TextFile)
@@ -2700,6 +2719,8 @@ hydroPSO <- function(
} # WHILE end
########################## END Main Iteration Loop #########################
+
+ if (normalise) X.best.part <- X.best.part * (upper.mat - lower.mat) + lower.mat
if (write2disk) {
close(OFout.Text.file)
diff --git a/man/hydroPSO.Rd b/man/hydroPSO.Rd
index e05f60e031f13656f680c7f1673fc575f837d3b6..9df5adcb2417fbb0ecec9048a16586ef4447e749 100755
--- a/man/hydroPSO.Rd
+++ b/man/hydroPSO.Rd
@@ -107,6 +107,15 @@ By default \code{c1= 0.5 + log(2)}
\item{c2}{
numeric, social acceleration coefficient. Encourages the exploration of the current global best and reflects how much the particle is influenced by the best-known optimum of the swarm \cr
By default \code{c2= 0.5 + log(2)}
+}
+ \item{use.IW}{
+logical, indicates if an inertia weight (\env{w}) will be used to avoid swarm explosion, i.e. particles flying around their best position without converging into it (see Shi and Eberhart, 1998) \cr
+By default \code{use.IW=TRUE}
+}
+ \item{IW.w}{
+OPTIONAL. Used only when \code{use.IW= TRUE \& IW.type!='GLratio'} \cr
+numeric, value of the inertia weight(s) (\env{w} or \env{[w.ini, w.fin]}). It can be a single number which is used for all iterations, or it can be a vector of length 2 with the initial and final values (in that order) that \env{w} will take along the iterations \cr
+By default \code{IW.w=1/(2*log(2))}
}
\item{use.CF}{
logical, indicates if the Clerc's Constriction Factor (see Clerc, 1999; Eberhart and Shi, 2000; Clerc and Kennedy, 2002) is used to avoid swarm explosion \cr
@@ -180,12 +189,13 @@ OPTIONAL. Only used when \code{method=='ipso'} \cr
numeric, number of particles considered in the computation of the global best \cr
By default \code{ngbest=4} (see Zhao, 2006)
}
- \item{use.IW}{
-logical, indicates if an inertia weight (\env{w}) will be used to avoid swarm explosion, i.e. particles flying around their best position without converging into it (see Shi and Eberhart, 1998) \cr
-By default \code{use.IW=TRUE}
+ \item{randomise}{
+logical, indicates whether the parameter values have to be normalised to the [0,1] interval during the optimisation or not\cr
+This option appears in the C and Matlab version of SPSO-2011 (See \url{http://www.particleswarm.info/standard_pso_2011_c.zip}) and there it is recommended to use this option when the search space is not an hypercube. If the search space is an hypercube, it is better not normalise (there is a small difference between the position without any normalisation and the de-normalised one)\cr
+By default \code{randomise=FALSE}
}
\item{IW.type}{
-OPTIONAL. Used only when \code{use.IW= TRUE} \cr
+OPTIONAL. Used only when \code{use.IW= TRUE} \& \code{length(IW.w)>1}\cr
character, defines how the inertia weight \env{w} will vary along iterations. Valid values are: \cr
-)\kbd{linear}: \env{w} varies linearly between the initial and final values specified in \code{IW.w} (see Shi and Eberhart, 1998; Zheng et al., 2003). This is the DEFAULT option \cr
-)\kbd{non-linear}: \env{w} varies non-linearly between the initial and final values specified in \code{IW.w} with exponential factor \kbd{IW.exp} (see Chatterjee and Siarry, 2006) \cr
@@ -193,11 +203,6 @@ character, defines how the inertia weight \env{w} will vary along iterations. Va
-)\kbd{aiwf}: adaptive inertia weight factor, where the inertia weight is varied adaptively depending on the goodness-of-fit values of the particles (see Liu et al., 2005) \cr
-)\kbd{GLratio}: \env{w} varies according to the ratio between the global best and the average of the particle's local best (see Arumugam and Rao, 2008) \cr
By default \code{IW.type='linear'}
-}
- \item{IW.w}{
-OPTIONAL. Used only when \code{use.IW= TRUE \& IW.type!='GLratio'} \cr
-numeric, value of the inertia weight(s) (\env{w} or \env{[w.ini, w.fin]}). It can be a single number which is used for all iterations, or it can be a vector of length 2 with the initial and final values (in that order) that \env{w} will take along the iterations \cr
-By default \code{IW.w=1/(2*log(2))}
}
\item{IW.exp}{
OPTIONAL. Used only when \code{use.IW= TRUE} AND \code{IW.type= 'non-linear'} \cr