diff --git a/DESCRIPTION b/DESCRIPTION
index c6686ca..a5202b9 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -9,6 +9,7 @@ Imports:
   phytools,
   subplex,
   DEoptim,
+  pso,
   Matrix,
   expm,
   SparseM,
@@ -60,6 +61,6 @@ Description:
  <DOI:10.1093/sysbio/syaa048>.
  Also contains functions to simulate the diversity-dependent
  process.
-RoxygenNote: 7.3.2
+RoxygenNote: 7.3.3
 Config/testthat/edition: 3
 URL: https://rsetienne.github.io/DDD/
diff --git a/R/dd_utils.R b/R/dd_utils.R
index 3577199..6343539 100644
--- a/R/dd_utils.R
+++ b/R/dd_utils.R
@@ -694,15 +694,18 @@ simplex = function(fun, trparsopt, optimpars, ...)
   invisible(out)
 }
 
-#' Carries out optimization (finding a minimum)
-#' 
-#' A wrapper to use several optimization routines, currently only 'simplex' (a
-#' method adopted from Matlab, or 'subplex', from the R package subplex). The
-#' function is called from several packages by the same author.
-#' 
-#' 
-#' @param optimmethod The method to use for optimization, either 'simplex' or
-#' 'subplex'
+#' Carries out optimization (finding a maximum)
+#' 
+#' A wrapper to use several optimization routines 
+#' 
+#' @param optimmethod The method to use for optimization. There is a choice of:
+#' 'simplex', which is a simplex algorithm that shows (when verbose > 0) the
+#' optimization trajectory.
+#' 'subplex', from the package with the same name
+#' 'DEoptim', from the package with the same name
+#' 'pso', from the package with the same name
+#' 'optim::name_of_algorithm', the algorithm from the optim package (e.g.
+#' "Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN", or "Brent")
 #' @param optimpars Parameters of the optimization: 1) relative tolerance in
 #' function arguments, 2) relative tolerance in function value, 3) absolute 
 #' tolerance in function arguments as well as the function value, 4) 
@@ -809,6 +812,23 @@ optimizer <- function(
                                               fun = fun,
                                               ...))
       outnew <- list(par = outnew$par, fvalues = -outnew$value, conv = outnew$convergence)
+    } else if(optimmethod == 'pso') {
+      minfun4 <- function(trparsopt, fun, ...)
+      {           
+        return(-fun(trparsopt = trparsopt, ...))
+      }
+      outnew <- suppressWarnings(pso::psoptim(par = trparsopt,
+                                              fn = minfun4,
+                                              lower = 0,
+                                              upper = 1,
+                                              control = list(reltol = optimpars[2],
+                                                             abstol = optimpars[3],
+                                                             maxit = optimpars[4]),
+                                              fun = fun,
+                                              ...))
+      outnew <- list(par = outnew$par, fvalues = -outnew$value, conv = outnew$convergence)
+    } else {
+      stop('The specified optimization method is not supported')
     }
     if(cy > 1 & (any(is.na(outnew$par)) | any(is.nan(outnew$par)) | is.na(outnew$fvalues) | is.nan(outnew$fvalues) | outnew$conv != 0))
     {
diff --git a/man/optimizer.Rd b/man/optimizer.Rd
index 0ef3346..eb2be3b 100644
--- a/man/optimizer.Rd
+++ b/man/optimizer.Rd
@@ -2,7 +2,7 @@
 % Please edit documentation in R/dd_utils.R
 \name{optimizer}
 \alias{optimizer}
-\title{Carries out optimization (finding a minimum)}
+\title{Carries out optimization (finding a maximum)}
 \usage{
 optimizer(
   optimmethod = "simplex",
@@ -15,8 +15,14 @@ optimizer(
 )
 }
 \arguments{
-\item{optimmethod}{The method to use for optimization, either 'simplex' or
-'subplex'}
+\item{optimmethod}{The method to use for optimization. There is a choice of:
+'simplex', which is a simplex algorithm that shows (when verbose > 0) the
+optimization trajectory.
+'subplex', from the package with the same name
+'DEoptim', from the package with the same name
+'pso', from the package with the same name
+'optim::name_of_algorithm', the algorithm from the optim package (e.g.
+"Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN", or "Brent")}
 
 \item{optimpars}{Parameters of the optimization: 1) relative tolerance in
 function arguments, 2) relative tolerance in function value, 3) absolute 
@@ -45,9 +51,7 @@ of the optimization routine}
 optimization converged (\code{conv})}.
 }
 \description{
-A wrapper to use several optimization routines, currently only 'simplex' (a
-method adopted from Matlab, or 'subplex', from the R package subplex). The
-function is called from several packages by the same author.
+A wrapper to use several optimization routines
 }
 \examples{