diff --git a/.gitignore b/.gitignore
index f1e0e11a..1c653b27 100644
--- a/.gitignore
+++ b/.gitignore
@@ -4,3 +4,8 @@
 *.Rproj
 
 inst/doc
+vignettes/*.R
+vignettes/*.html
+vignettes/*.png
+
+
diff --git a/NEWS.md b/NEWS.md
index 0ed909ab..f34188c9 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,3 +1,7 @@
+# matlib 0.5.2
+
+- added `swp()` function
+
 # matlib 0.5.1
 
 - added `len()` convenience function for Euclidean lengths
diff --git a/R/swp.R b/R/swp.R
new file mode 100644
index 00000000..755e9760
--- /dev/null
+++ b/R/swp.R
@@ -0,0 +1,56 @@
+
+
+#' The Matrix Sweep Operator
+#'
+#' The \code{swp} function "sweeps" a matrix on the rows and columns given in \code{index} to produce a new matrix
+#' with those rows and columns partialled out. This was defined as a fundamental statistical operation in
+#' multivariate methods by Beaton (1964) and expanded by Dempster (1969).
+#'
+#' If \code{M} is the partitioned matrix
+#' \deqn{\left[ \begin{array}{cc} \mathbf {R} &  \mathbf {S} \\ \mathbf {T} &  \mathbf {U} \end{array} \right]}
+#' where \eqn{R} is \eqn{q \times q} then \code{swp(M, 1:q)} gives
+#' \deqn{\left[ \begin{array}{cc} \mathbf {R}^{-1} &  \mathbf {R}^{-1}\mathbf {S} \\ -\mathbf {TR}^{-1} &  \mathbf {U}-\mathbf {TR}^{-1}\mathbf {S} \\ \end{array} \right]}
+#'
+#' @param M a numeric matrix
+#' @param index a numeric vector indicating the rows/columns to be swept.  The entries must be less than or equal
+#'     to the number or rows or columns in \code{M}.  If missing, the function sweeps on all rows/columns \code{1:min(dim(M))}.
+#'
+#' @return the matrix \code{M} with rows and columns in \code{indices}
+#' @references Beaton, A. E. (1964), \emph{The Use of Special Matrix Operations in Statistical Calculus}, Princeton, NJ: Educational Testing Service.
+#'
+#'      Dempster, A. P. (1969) \emph{Elements of Continuous Multivariate Analysis}. Addison-Wesley Publ. Co., Reading, Mass.
+#'
+#' @examples
+#' data(therapy)
+#' mod3 <- lm(therapy ~ perstest + IE + sex, data=therapy)
+#' X <- model.matrix(mod3)
+#' XY <- cbind(X, therapy=therapy$therapy)
+#' XY
+#' M <- crossprod(XY)
+#' swp(M, 1)
+#' swp(M, 1:2)
+
+swp <- function (M, index)
+{
+  p <- ncol(M)
+  if (missing(index)) index <- 1:min(dim(M))
+  u <- is.na(match(1:p, index))
+  a <- (1:p)[u]
+  out <- 0 * M
+  dimnames(out) <- dimnames(M)
+  if (length(a) == 0)
+    return(-solve(M))
+  else if (length(a) == p)
+    return(M)
+  else {
+    Saa <- M[a, a, drop = FALSE]
+    Sab <- M[a, index, drop = FALSE]
+    Sbb <- M[index, index, drop = FALSE]
+    B <- Sab %*% solve(Sbb)
+    out[a, a] <- Saa - B %*% t(Sab)
+    out[a, index] <- B
+    out[index, a] <- t(B)
+    out[index, index] <- -solve(Sbb)
+    return(out)
+  }
+}
diff --git a/man/swp.Rd b/man/swp.Rd
new file mode 100644
index 00000000..87b7f599
--- /dev/null
+++ b/man/swp.Rd
@@ -0,0 +1,44 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/swp.R
+\name{swp}
+\alias{swp}
+\title{The Matrix Sweep Operator}
+\usage{
+swp(M, index)
+}
+\arguments{
+\item{M}{a numeric matrix}
+
+\item{index}{a numeric vector indicating the rows/columns to be swept.  The entries must be less than or equal
+to the number or rows or columns in \code{M}.  If missing, the function sweeps on all rows/columns \code{1:min(dim(M))}.}
+}
+\value{
+the matrix \code{M} with rows and columns in \code{indices}
+}
+\description{
+The \code{swp} function "sweeps" a matrix on the rows and columns given in \code{index} to produce a new matrix
+with those rows and columns partialled out. This was defined as a fundamental statistical operation in
+multivariate methods by Beaton (1964) and expanded by Dempster (1969).
+}
+\details{
+If \code{M} is the partitioned matrix
+\deqn{\left[ \begin{array}{cc} \mathbf {R} &  \mathbf {S} \\ \mathbf {T} &  \mathbf {U} \end{array} \right]}
+where \eqn{R} is \eqn{q \times q} then \code{swp(M, 1:q)} gives
+\deqn{\left[ \begin{array}{cc} \mathbf {R}^{-1} &  \mathbf {R}^{-1}\mathbf {S} \\ -\mathbf {TR}^{-1} &  \mathbf {U}-\mathbf {TR}^{-1}\mathbf {S} \\ \end{array} \right]}
+}
+\examples{
+data(therapy)
+mod3 <- lm(therapy ~ perstest + IE + sex, data=therapy)
+X <- model.matrix(mod3)
+XY <- cbind(X, therapy=therapy$therapy)
+XY
+M <- crossprod(XY)
+swp(M, 1)
+swp(M, 1:2)
+}
+\references{
+Beaton, A. E. (1964), \emph{The Use of Special Matrix Operations in Statistical Calculus}, Princeton, NJ: Educational Testing Service.
+
+     Dempster, A. P. (1969) \emph{Elements of Continuous Multivariate Analysis}. Addison-Wesley Publ. Co., Reading, Mass.
+}
+