---
title: "Exceedingly Simple Sorting with Indices"
author: "Jan de Leeuw"
date: "Version 01, March 22, 2017"
output:
pdf_document:
keep_tex: yes
number_sections: yes
toc: yes
toc_depth: 3
html_document:
keep_md: yes
number_sections: yes
toc: yes
graphics: yes
fontsize: 12pt
abstract: We use the system `qsort` to write a routine that produces both the sort
an the order of a vector of doubles.
---
```{r function_code, echo = FALSE}
source("mySort.R")
```
```{r packages, echo = FALSE}
options (digits = 10)
library (microbenchmark)
```
Note: This is a working paper which will be expanded/updated frequently. All suggestions for improvement are welcome. The directory [gifi.stat.ucla.edu/mysort](http://gifi.stat.ucla.edu/mysort) has a pdf version, the complete Rmd file with the code chunks,
and the R source code.
#Introduction
This may not be of much use to anyone else, but I'll make it available anyway. For a larger project I needed something simple in C to sort a vector and return both
the sorted vector and the corresponding index vector. Thus it combines `sort` and `order` in R. It creates an array of structures
in C, where each structure has a value member and an index member. It then sorts the structures using the system `qsort`, by comparison
of the values. Finally it splits the sorted array of structures into values and indices and returns those to R in a two-element list.
#Example
```{r example}
set.seed (12345)
x <- rnorm (10)
mySort (x)
```
as compared to
```{r standard}
list (values = sort (x), indices = order (x))
```
A quick comparison shows very little timing difference.
```{r compare}
f <- function () {
x <- rnorm(10000)
h <- mySort (x)
}
g <- function () {
x <- rnorm (10000)
h <- list (values = sort (x), indices = order (x))
}
microbenchmark (f, g, times = 1000)
```
#C code
```{c c_code, engine='c', results='hide'}
#include
#include
int myComp(const void *, const void *);
void mySort(double *, int *, int *);
struct couple {
double value;
int index;
};
int myComp(const void *px, const void *py) {
double x = ((struct couple *)px)->value;
double y = ((struct couple *)py)->value;
return (int)copysign(1.0, x - y);
}
void mySort(double *x, int *k, int *n) {
int nn = *n;
struct couple *xi =
(struct couple *)calloc((size_t)nn, (size_t)sizeof(struct couple));
for (int i = 0; i < nn; i++) {
xi[i].value = x[i];
xi[i].index = i + 1;
}
(void)qsort(xi, (size_t)nn, (size_t)sizeof(struct couple), myComp);
for (int i = 0; i < nn; i++) {
x[i] = xi[i].value;
k[i] = xi[i].index;
}
free(xi);
}
```
#R code
```{r r_code}
dyn.load("mySort.so")
mySort <- function (x) {
h <-
.C(
"mySort",
values = as.double (x),
indices = as.integer(1:length(x)),
as.integer(length(x))
)
return (list (values = h$values, indices = h$indices))
}
```