Scythe Bugs and Missing Features.

Bugs:
- mean, sd, var, matrix methods etc are templated to be the same type as the
  matrix.  However these functions should return double, never int.

- Some matrix functions give compile errors because for int.  For instance
  sd() on a int matrix type fails because std:pow(int, int) is not valid.

- Something is messed up with rexp(). It's missing in the help file but:
  http://scythe.wustl.edu/api/classscythe_1_1rng_55184ef69537f8a47f82c1df4867d5dc.html#55184ef69537f8a47f82c1df4867d5dc
  Compiling the help using the latest doxygen seems to correct the problem, 
  although macros don't get documented really.

- Docs for scythe::min refer to max.

- What is linesearch? Docs are not very clear.

- Reading in from file is fragile.  outsheet in stata on OSX does not
  use the proper newline characters.  newlines across dos/unix/osx differ.

- col_interchange doesn't sound like it is what it is in help file.  It swaps
	rows rather than reordering.

- qnorm1 has return type "double double sd double" in docs.

Missing / Potential Features:

- rtnorm is parameterized using the variance whereas rnorm is parameterized
  with sd.  Possibility of errors.

- Fast versions of functions proportional to log density.  In most MCMC
  applications sum(log(scythe::dnorm)), for instance, takes up a significant
  amount of time.  However we usually just need the log of something
  proportional to density (for Slice or MH).  For instance, 
  dnorm(x, mean, sd, fastlog=true) = -(x-mean)^2/(2 * sd^2), which avoids
  and exp, sqrt, log, and a few other operations.  I get a 10
  fold increase in speed (of my entire MCMC algorithm) by replacing a
  sum(log(dnorm)) with 

	// faster sum log dnorm.
	template <typename T, matrix_order PO, matrix_style PS>
	double sum_fast_log_dnorm(const Matrix<T,PO,PS> &A, double mean, double sd)
	{
	  double result = 0;
	  for(int i = 0; i < A.size(); ++i) {
	    result += fast_log_dnorm(A[i], mean, sd);
	  }
	  return result;
	}

	double fast_log_dnorm(double x, double mean, double sd) {
	  return -pow((x - mean),2)/(2*pow(sd,2));
	}

  I haven't figured it out, but this should also be vectorizable using SSE
  and multiple cores.  Could lead to a lead to a 40 fold total increase in
  speed.

- qnorm1 is the only inverse cdf available.  It would be nice to have
  all the q functions available.  These should be in R, no?

- Missing some common random types, like multinomial.

- Upgrade MT to use SSE instructions. Honestly, however, random number
  generation speed does not seem to matter at all.  It's all the density
  calls, exp, log, pow, sqrt, which is where optimization should focus.

- How hard would it be to add SSE instructions to log and exp, at the very 
  least?  Might get 5 fold increase in speed, without MKL.
  See: http://gruntthepeon.free.fr/ssemath/

- The easiest path towards performance would be to fix things allow
  the Intel MKL to do optimizations.  I got 2-3x performance increase using
  Intel C++ compiler over g++4.01 without changing anything.
  Using getArray(), it would
  be easy to wrap MKL functions into Scythe with a day or two of work.  
  I suspect gains would be large, since it takes advantage of SSE, multi cores
  etc. Then add a -DHAS_MKL flag. See:
	http://www.intel.com/software/products/mkl/docs/WebHelp/mkl.htm