Commit 2a877243 authored by Liam M. Healy's avatar Liam M. Healy
Browse files

Merge branch 'antik' 

Conflicts:
	gsll.asd
	ordinary-differential-equations/ode-system.lisp
	random/dirichlet.lisp
parents dfdabff8 13e65854

basis-splines.lisp

+7 −7
Original line number Diff line number Diff line 
;; Basis splines.

;; Liam Healy 2008-02-18 14:43:20EST basis-splines.lisp

;; Time-stamp: <2010-06-30 19:57:28EDT basis-splines.lisp>

;; Time-stamp: <2011-05-26 12:37:29EDT basis-splines.lisp>

;;

;; Copyright 2008, 2009 Liam M. Healy

;; Copyright 2008, 2009, 2011 Liam M. Healy

;; Distributed under the terms of the GNU General Public License

;;

;; This program is free software: you can redistribute it and/or modify
@@ -125,18 +125,18 @@ 
      (let* ((xi (coerce (* i (/ 15 (1- ndata))) 'double-float))

	     (yi (+ (* (cos xi) (exp (* -0.1d0 xi)))

		    (sample rng :gaussian :sigma sigma))))

	(setf (grid:gref x i) xi

	      (grid:gref y i) yi

	      (grid:gref w i) (/ (expt sigma 2)))))

	(setf (grid:aref x i) xi

	      (grid:aref y i) yi

	      (grid:aref w i) (/ (expt sigma 2)))))

    ;; Uniform breakpoints [0, 15]

    (uniform-knots 0.0d0 15.0d0 bw)

    ;; Fit matrix

    (dotimes (i ndata)

      ;; Compute B_j

      (evaluate bw (grid:gref x i) :b B)

      (evaluate bw (grid:aref x i) :b B)

      ;; Fill in row i of X

      (dotimes (j ncoeffs)

	(setf (grid:gref Xmatrix i j) (grid:gref B j))))

	(setf (grid:aref Xmatrix i j) (grid:aref B j))))

    ;; Do the fit

    (linear-mfit Xmatrix y c w nil cov mw)

    ;; Return the smoothed curve

calculus/monte-carlo.lisp

+5 −5
Original line number Diff line number Diff line 
;; Monte Carlo Integration

;; Liam Healy Sat Feb  3 2007 - 17:42

;; Time-stamp: <2010-08-07 21:41:14EDT monte-carlo.lisp>

;; Time-stamp: <2011-01-10 17:59:20EST monte-carlo.lisp>

;;

;; Copyright 2007, 2008, 2009 Liam M. Healy

;; Copyright 2007, 2008, 2009, 2011 Liam M. Healy

;; Distributed under the terms of the GNU General Public License

;;

;; This program is free software: you can redistribute it and/or modify
@@ -51,7 +51,7 @@ 
	      (scalars t))

  "gsl_monte_plain_integrate"

  ((callback :pointer)

   ((foreign-pointer lower-limits) :pointer) ((foreign-pointer upper-limits) :pointer)

   ((grid:foreign-pointer lower-limits) :pointer) ((grid:foreign-pointer upper-limits) :pointer)

   ((dim0 lower-limits) sizet) (number-of-samples sizet)

   ((mpointer generator) :pointer)

   ((mpointer state) :pointer)
@@ -112,7 +112,7 @@ 
	      (scalars t))

  "gsl_monte_miser_integrate"

  ((callback :pointer)

   ((foreign-pointer lower-limits) :pointer) ((foreign-pointer upper-limits) :pointer)

   ((grid:foreign-pointer lower-limits) :pointer) ((grid:foreign-pointer upper-limits) :pointer)

   ((dim0 lower-limits) sizet) (number-of-samples sizet)

   ((mpointer generator) :pointer)

   ((mpointer state) :pointer)
@@ -173,7 +173,7 @@ 
	      (scalars t))

  "gsl_monte_vegas_integrate"

  ((callback :pointer)

   ((foreign-pointer lower-limits) :pointer) ((foreign-pointer upper-limits) :pointer)

   ((grid:foreign-pointer lower-limits) :pointer) ((grid:foreign-pointer upper-limits) :pointer)

   ((dim0 lower-limits) sizet) (number-of-samples sizet)

   ((mpointer generator) :pointer)

   ((mpointer state) :pointer)

calculus/numerical-integration.lisp

+3 −3
Original line number Diff line number Diff line 
;; Numerical integration

;; Liam Healy, Wed Jul  5 2006 - 23:14

;; Time-stamp: <2010-08-07 21:46:21EDT numerical-integration.lisp>

;; Time-stamp: <2011-01-10 17:59:21EST numerical-integration.lisp>

;;

;; Copyright 2006, 2007, 2008, 2009, 2010 Liam M. Healy

;; Copyright 2006, 2007, 2008, 2009, 2010, 2011 Liam M. Healy

;; Distributed under the terms of the GNU General Public License

;;

;; This program is free software: you can redistribute it and/or modify
@@ -157,7 +157,7 @@ 
	      (limit 1000) (workspace (make-integration-workspace limit)))

  "gsl_integration_qagp"

  ((callback :pointer)

   ((foreign-pointer points) :pointer) ((dim0 points) sizet)

   ((grid:foreign-pointer points) :pointer) ((dim0 points) sizet)

   (absolute-error :double) (relative-error :double) (limit sizet)

   ((mpointer workspace) :pointer)

   (result (:pointer :double)) (abserr (:pointer :double)))

data/both.lisp

+17 −17
Original line number Diff line number Diff line 
;; Functions for both vectors and matrices.

;; Liam Healy 2008-04-26 20:48:44EDT both.lisp

;; Time-stamp: <2010-11-25 09:31:17EST both.lisp>

;; Time-stamp: <2011-05-26 12:37:36EDT both.lisp>

;;

;; Copyright 2008, 2009, 2010 Liam M. Healy

;; Copyright 2008, 2009, 2010, 2011 Liam M. Healy

;; Distributed under the terms of the GNU General Public License

;;

;; This program is free software: you can redistribute it and/or modify
@@ -35,13 +35,13 @@ 
  :export nil

  :documentation "Allocate memory for the GSL struct given a block pointer.")



(defmfun alloc-from-block ((object matrix) blockptr)

(defmfun alloc-from-block ((object grid:matrix) blockptr)

  ("gsl_" :category :type "_alloc_from_block")

  ((blockptr :pointer)

   (0 sizet)				; offset

   ((first (dimensions object)) sizet)	; number of rows

   ((second (dimensions object)) sizet)	; number of columns

   ((second (dimensions object)) sizet))	; "tda" = number of columns for now

   ((first (grid:dimensions object)) sizet)	; number of rows

   ((second (grid:dimensions object)) sizet)	; number of columns

   ((second (grid:dimensions object)) sizet))	; "tda" = number of columns for now

  :definition :methods

  :c-return :pointer

  :export nil)
@@ -81,7 +81,7 @@ 
;;;;****************************************************************************

;;;; Array elements; used in callbacks scalarsp=T only

;;;;****************************************************************************

;;; Normal foreign array access is with grid:gref, but in order to

;;; Normal foreign array access is with grid:aref, but in order to

;;; avoid the overhead of instantiating a foreign-array object to

;;; access components, we use these macros which expand to gsl_*_get

;;; and gsl_*_set.
@@ -113,7 +113,7 @@ 
    `(cffi:foreign-funcall

      ,(actual-gsl-function-name

	`("gsl_" :category :type ,(if value "_set" "_get"))

	(if matrixp 'matrix 'vector)

	(if matrixp 'grid:matrix 'vector)

	element-type)

      :pointer ,mpointer

      sizet ,(first indices)
@@ -192,7 +192,7 @@ 
  :documentation			; FDL

  "Add the scalar complex x to all the elements of array a.")



(defmethod elt+ ((x float) (a foreign-array))

(defmethod elt+ ((x float) (a grid:foreign-array))

  (elt+ a x))

  

(defmfun elt- ((a both) (b both))
@@ -207,7 +207,7 @@ 
  "Subtract the elements of b from the elements of a.

   The two must have the same dimensions.")



(defmethod elt- ((a foreign-array) (x float))

(defmethod elt- ((a grid:foreign-array) (x float))

  (elt+ a (- x)))



(defmfun elt* ((a vector) (b vector))
@@ -222,7 +222,7 @@ 
  "Multiply the elements of a by the elements of b.

   The two must have the same dimensions.")



(defmfun elt* ((a matrix) (b matrix))

(defmfun elt* ((a grid:matrix) (b grid:matrix))

  ("gsl_" :category :type "_mul_elements")

  (((mpointer a) :pointer) ((mpointer b) :pointer))

  :definition :methods
@@ -243,7 +243,7 @@ 
  "Divide the elements of a by the elements of b.

   The two must have the same dimensions.")



(defmfun elt/ ((a matrix) (b matrix))

(defmfun elt/ ((a grid:matrix) (b grid:matrix))

  ("gsl_" :category :type "_div_elements")

  (((mpointer a) :pointer) ((mpointer b) :pointer))

  :definition :methods
@@ -252,7 +252,7 @@ 
  :outputs (a)

  :return (a))



(defmethod elt/ ((a foreign-array) (x number))

(defmethod elt/ ((a grid:foreign-array) (x number))

  (elt* a (/ x)))



(defmfun elt* ((a both) (x float))
@@ -279,7 +279,7 @@ 
  :documentation			; FDL

  "Multiply the elements of a by the scalar complex factor x.")



(defmethod elt* ((x float) (a foreign-array))

(defmethod elt* ((x float) (a grid:foreign-array))

  (elt* a x))



;;;;****************************************************************************
@@ -329,7 +329,7 @@ 
  "The index of the minimum value in a.  When there are several

  equal minimum elements, then the lowest index is returned.")



(defmfun min-index ((a matrix))

(defmfun min-index ((a grid:matrix))

  ("gsl_" :category :type "_min_index")

  (((mpointer a) :pointer) (imin (:pointer sizet)) (jmin (:pointer sizet)))

  :definition :methods
@@ -348,7 +348,7 @@ 
  "The index of the maximum value in a.  When there are several

  equal maximum elements, then the lowest index is returned.")



(defmfun max-index ((a matrix))

(defmfun max-index ((a grid:matrix))

  ("gsl_" :category :type "_max_index")

  (((mpointer a) :pointer) (imin (:pointer sizet)) (jmin (:pointer sizet)))

  :definition :methods
@@ -369,7 +369,7 @@ 
  returned.  Returned indices are minimum, maximum; for matrices

  imin, jmin, imax, jmax.")



(defmfun minmax-index ((a matrix))

(defmfun minmax-index ((a grid:matrix))

  ("gsl_" :category :type "_minmax_index")

  (((mpointer a) :pointer)

   (imin (:pointer sizet)) (jmin (:pointer sizet))

data/combination.lisp

+4 −4
Original line number Diff line number Diff line 
;; Combinations

;; Liam Healy, Sun Mar 26 2006 - 11:51

;; Time-stamp: <2010-07-16 17:14:03EDT combination.lisp>

;; Time-stamp: <2011-01-10 18:16:28EST combination.lisp>

;;

;; Copyright 2006, 2007, 2008, 2009, 2010 Liam M. Healy

;; Copyright 2006, 2007, 2008, 2009, 2010, 2011 Liam M. Healy

;; Distributed under the terms of the GNU General Public License

;;

;; This program is free software: you can redistribute it and/or modify
@@ -37,7 +37,7 @@ 
    (setf (grid:metadata-slot object 'mpointer)

	  mptr

	  (cffi:foreign-slot-value mptr 'gsl-combination-c 'data)

	  (foreign-pointer object)

	  (grid:foreign-pointer object)

	  (cffi:foreign-slot-value mptr 'gsl-combination-c 'range)

	  range

	  (cffi:foreign-slot-value mptr 'gsl-combination-c 'size)
@@ -56,7 +56,7 @@ 
	     (make-instance

	      'combination

	      :element-type '(unsigned-byte #+int64 64 #+int32 32)

	      :range (combination-range n) :dimensions (dimensions k))

	      :range (combination-range n) :dimensions (grid:dimensions k))

	     (make-instance

	      'combination

	      :element-type '(unsigned-byte #+int64 64 #+int32 32)