From 9c64554612b8a45dec406b1ba3f3d1e4705515f8 Mon Sep 17 00:00:00 2001
From: Mauricio Zambrano-Bigiarini <hzambran@users.noreply.github.com>
Date: Mon, 17 Sep 2012 15:19:43 +0000
Subject: [PATCH] swarm size and Vini.type are now SPSO2011 capable

---
 DESCRIPTION     |  4 ++--
 NEWS            | 13 +++++++++++++
 R/PSO_v2012.R   | 39 +++++++++++++++++++++++++--------------
 man/hydroPSO.Rd |  6 ++++--
 4 files changed, 44 insertions(+), 18 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index 6ed5bba..5786b1b 100755
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: hydroPSO
 Type: Package
 Title: Model-Independent Particle Swarm Optimisation for Environmental Models
-Version: 0.1-58
-Date: 2012-09-14
+Version: 0.1-58-1
+Date: 2012-09-17
 Author: Mauricio Zambrano-Bigiarini [aut, cre] and Rodrigo Rojas [ctb]
 Authors@R: c(person("Mauricio", "Zambrano-Bigiarini", email="mzb.devel@gmail.com", role=c("aut","cre")), person("Rodrigo", "Rojas", email="Rodrigo.RojasMujica@gmail.com", role=c("ctb")) )
 Maintainer: Mauricio Zambrano-Bigiarini <mzb.devel@gmail.com>
diff --git a/NEWS b/NEWS
index 3e7bdf5..6b2c909 100755
--- a/NEWS
+++ b/NEWS
@@ -1,5 +1,18 @@
 NEWS/ChangeLog for hydroPSO
 --------------------------
+0.1-58	(under-development)
+        o 'hydroPSO'      : -) towards SPSO 2011 capable...
+        
+                            -) 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'
+                            -) argument 'Vini.type' now allows the following 2 new values: 'random2011' and 'lhs2011', which set the initial velocity 
+                               values for the particles to random values according to the equation specified in the SPSO 2011, with a uniform distribution 
+                               or with a Latin-hypercube sampling, respectively. The old values 'random' and 'lhs' were replaced by the new 
+                               Default ? (so far, SPSO 2011)
+                            -) argument 'Vini.type'. The old values 'random' and 'lhs' were replaced by 'random2007' and 'lhs2007', in order to make 
+                               clear that the follow the equation described in SPSO 2007.
+        
 0.1-58	14-Sep-2012
         o 'hydroPSO'      : -) 'random.update' is now ONLY used when 'best.update="async". In hydroPSO 0.1-57 'random.update' was set to TRUE 
                                by default, independent of the 'best.update="sync" value.
diff --git a/R/PSO_v2012.R b/R/PSO_v2012.R
index e9db594..55cfe18 100755
--- a/R/PSO_v2012.R
+++ b/R/PSO_v2012.R
@@ -121,7 +121,7 @@ compute.CF <- function(c1, c2) {
 # 'topology'  : character, with the topology to be used in PSO. Valid values
 #               are in c('gbest', 'lbest')
 # 'method'    : character, with the method to be used as PSO algorithm. Valid values
-#               are in c('pso', 'ipso', 'fips', 'wfips')
+#               are in c('spso2007', 'spso2011', 'ipso', 'fips', 'wfips')
 
 # Result      : vector of 'n' velocities, one for each parameter, corresponding to the current particle
 ################################################################################
@@ -145,7 +145,7 @@ compute.veloc <- function(x, v, w, c1, c2, CF, Pbest, part.index, gbest,
   r1 <- runif(n, min=0, max=1)
   r2 <- runif(n, min=0, max=1)
   
-  if ( method=="pso" ) {
+  if ( method=="spso2007" ) {
   
     ifelse(part.index != localBest.pos,  
            vn <- CF * ( w*v + r1*c1*(pbest-x) + r2*c2*(localBest-x) ),
@@ -807,6 +807,7 @@ InitializateX <- function(npart, param.IDs, x.MinMax, x.ini.type) {
 # Author : Mauricio Zambrano-Bigiarini
 # Started: 24-Dec-2010
 # Updates: 24-Nov-2011
+#          17-Sep-2012
 ################################################################################
 # Purpose: Function for the initialization of the position and the velocities 
 #          of all the particles in the swarm
@@ -822,7 +823,7 @@ InitializateX <- function(npart, param.IDs, x.MinMax, x.ini.type) {
 #                 Second column has the maximum possible value for each parameter
 # 'v.ini'      : character, indicating how to carry out the initialization 
 #                of the velocitites of all the particles in the swarm
-#                valid values are in c(0, 'random', 'lhs') 
+#                valid values are in c('zero', 'random2007', 'lhs2007', 'random2011', 'lhs2011') 
 InitializateV <- function(npart, param.IDs, x.MinMax, v.ini.type, Xini) {
 
   # Number of parameters
@@ -833,13 +834,17 @@ InitializateV <- function(npart, param.IDs, x.MinMax, v.ini.type, Xini) {
   # Rows = 'npart'; 
   # Columns = 'n' (Dimension of the Solution Space)
   # Random bounded values are assigned to each dimension
-  if ( v.ini.type=="random" ) {
+  if ( v.ini.type=="random2007" ) {
       V <- ( Random.Bounded.Matrix(npart, x.MinMax) - Xini)/2
-  } else if ( v.ini.type=="lhs" ) {
+  } else if ( v.ini.type=="lhs2007" ) {
       V <- ( rLHS(npart, x.MinMax) - Xini)/2
     } else if ( v.ini.type=="zero" ) {
         V <- matrix(0, ncol=n, nrow=npart, byrow=TRUE)    
-      }
+      } else if ( v.ini.type=="random2011" ) {
+          V <- Random.Bounded.Matrix(npart, (x.MinMax - cbind(Xini, Xini) ) 
+        } else if ( v.ini.type=="lhs2011" ) {
+            V <- rLHS(npart, x.MinMax - cbind(Xini, Xini) )
+          } 
   colnames(V) <- param.IDs 
   rownames(V) <- paste("Part", 1:npart, sep="")
 
@@ -1378,7 +1383,7 @@ hydroPSO <- function(
                     par, 
                     fn="hydromod",  
                     ...,
-                    method=c("pso", "ipso", "fips", "wfips"),
+                    method=c("spso2011", "spso2007", "ipso", "fips", "wfips"),
                     lower=-Inf,
                     upper=Inf,                
                     control=list(),
@@ -1432,15 +1437,18 @@ hydroPSO <- function(
 	    c1= 0.5+log(2), 
 	    c2= 0.5+log(2), 
 	    use.CF= FALSE, 
-	    lambda= 1, 
+	    lambda= 1,
+
 	    abstol= NULL,    
 	    reltol=sqrt(.Machine$double.eps),             
 	    Xini.type=c("lhs", "random"),  
-	    Vini.type=c("lhs", "random", "zero"), 
+	    Vini.type=c("lhs2007", "random2007", "zero", "lhs2011", "random2011"), 
 	    best.update=c("sync", "async"),
 	    random.update=TRUE,
 	    boundary.wall=c("reflecting", "damping", "absorbing", "invisible"),
-	    topology=c("random", "gbest", "lbest", "vonNeumann"), K=3, iter.ini=0, ngbest=4, # only used when 'method=ipso'   
+	    topology=c("random", "gbest", "lbest", "vonNeumann"), K=3, 
+	    iter.ini=0, # only used when 'topology=lbest'   
+	    ngbest=4,   # only used when 'method=ipso'   
 
 	    use.IW = TRUE, IW.type=c("linear", "non-linear", "runif", "aiwf", "GLratio"), IW.w=1/(2*log(2)), IW.exp= 1, 
 	    use.TVc1= FALSE, TVc1.type=c("non-linear", "linear", "GLratio"), TVc1.rng= c(1.28, 1.05), TVc1.exp= 1.5, 
@@ -1457,7 +1465,7 @@ hydroPSO <- function(
 	    REPORT=100 
 	       )
 
-    MinMax        <- match.arg(control[["MinMax"]], con[["MinMax"]])
+    MinMax        <- match.arg(control[["MinMax"]], con[["MinMax"]])    
     Xini.type     <- match.arg(control[["Xini.type"]], con[["Xini.type"]]) 
     Vini.type     <- match.arg(control[["Vini.type"]], con[["Vini.type"]]) 
     best.update   <- match.arg(control[["best.update"]], con[["best.update"]]) 
@@ -1476,9 +1484,12 @@ hydroPSO <- function(
     drty.in           <- con[["drty.in"]]
     drty.out          <- con[["drty.out"]]
     param.ranges      <- con[["param.ranges"]]         
-    digits            <- con[["digits"]]
-
-    npart             <- ifelse(is.na(con[["npart"]]),ceiling(10+2*sqrt(n)),con[["npart"]])
+    digits            <- con[["digits"]]                
+    npart             <- ifelse(is.na(con[["npart"]]), 
+                                ifelse(method %in% c("spso2007", "spso2011"), 
+                                       ifelse(method=="spso2007", ceiling(10+2*sqrt(n)), 40),
+                                       40), 
+                                con[["npart"]] )                                 
     maxit             <- con[["maxit"]] 
     maxfn             <- con[["maxfn"]] 
     c1                <- con[["c1"]] 
diff --git a/man/hydroPSO.Rd b/man/hydroPSO.Rd
index 8b2da35..ce0d998 100755
--- a/man/hydroPSO.Rd
+++ b/man/hydroPSO.Rd
@@ -14,7 +14,7 @@ Enhanced Particle Swarm Optimisation algorithm
 Particle Swarm Optimisation algorithm to calibrate environmental models. It includes a series of controlling options and PSO variants to improve the performance of the algorithm and customize it to different calibration problems
 }
 \usage{
-hydroPSO(par, fn= "hydromod", ..., method=c("pso", "ipso", "fips", "wfips"),
+hydroPSO(par, fn= "hydromod", ..., method=c("spso2011", "spso2007", "ipso", "fips", "wfips"),
          lower=-Inf, upper=Inf, control=list(), 
          model.FUN=NULL, model.FUN.args=list() )
 }
@@ -34,13 +34,15 @@ OPTIONAL. Only used when \code{fn!='hydromod'}. \cr
 further arguments to be passed to \code{fn}.
 }
   \item{method}{
-character, variant of the PSO algorithm to be used. Valid values are in \code{c('pso', 'ipso', 'fips', 'wfips')}: \cr
+character, variant of the PSO algorithm to be used. Valid values are in \code{c('spso2011', 'spso2007', 'ipso', 'fips', 'wfips')}: \cr
 
 \kbd{pso}: at each iteration particles are attracted to its own best-known personal and to the best-known global position. Each particle is connected to a neighbourhood of particles depending on the \code{topology} value \cr
 \kbd{ipso}: at each iteration particles in the swarm are rearranged in descending order according to their goodness-of-fit and the best \code{ngbest} particles are used to modify particles' position and velocity (see Zhao, 2006). Each particle is connected to a neighbourhood of particles depending on the \code{topology} value \cr
 \kbd{fips}: at each iteration ALL particles contribute to modify the particles' position and velocity (see Mendes et al., 2004). Each particle is connected to a neighbourhood of particles depending on the \code{topology} value \cr
 \kbd{wfips}: same implementation as \kbd{fips} method, but the contribution of each particle is weighted according to their goodness-of-fit value (see Mendes et al., 2004) \cr
 By default \code{method=pso}
+\kbd{spso2007}: As in \code{method='pso'}, but with values of the PSO engine set to the values defined in the Standard PSO 2007 (SPSO 2007, see Clerc 2011) \cr
+\kbd{spso2011}: As in \code{method='pso'}, but with values of the PSO engine set to the values defined in the Standard PSO 2011 (SPSO 2011, see Clerc 2011)
 }
   \item{lower}{
 numeric, lower boundary for each parameter \cr
-- 
GitLab