kclglue.lsp

   1 ;;;; kclglue -- Interface to C library
   2 ;;;;
   3 ;;;; Copyright (c) 1991, by Luke Tierney. Permission is granted for
   4 ;;;; unrestricted use.
   5
   6 (in-package 'lisp-stat-basics)
   7
   8 (eval-when (compile load eval)
   9            (set-macro-character
  10             #\%
  11             #'(lambda (stream char) (values (read-line stream)))))
  12
  13 (Clines
  14 %#define IN_KCL_GLUE
  15 %#include "lib/linalg.h"
  16 %extern double rcondest_front();
  17 %extern char *calloc();
  18 %char buf[1000];
  19 )
  20
  21 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  22 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  23 ;;;;
  24 ;;;;                           Basic Utilities
  25 ;;;;
  26 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  27 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  28
  29 ;;;;
  30 ;;;;                Storage Allocation Functions
  31 ;;;;
  32
  33 (defun null-ptr-p (p) (= p 0))
  34 (defun ptr-eq (p q) (= p q))
  35
  36 (Clines
  37 %int la_base_allocate(n, m)
  38 %       unsigned n, m;
  39 %{
  40 %  char *p = calloc(n, m);
  41 %  if (p == nil) xlfail("allocation failed");
  42 %  return((int) p);
  43 %}
  44 )
  45
  46 (defentry la-base-allocate (int int) (int "la_base_allocate"))
  47
  48 (Clines
  49 %int la_base_free_alloc(p)
  50 %       int p;
  51 %{
  52 %  if (p) free((char *) p);
  53 %  return(0);
  54 %}
  55 )
  56
  57 (defentry la-base-free (int) (int "la_base_free_alloc"))
  58
  59 (Clines
  60 %static int mode_size(mode)
  61 %       int mode;
  62 %{
  63 %  switch (mode) {
  64 %  case IN: return(sizeof(int));
  65 %  case RE: return(sizeof(double));
  66 %  case CX: return(sizeof(Complex));
  67 %  }
  68 %  return(0);
  69 %}
  70 )
  71
  72 (defentry la-mode-size (int) (int "mode_size"))
  73
  74 (defCfun "int la_allocate(n, m) int n, m;" 0
  75 %{
  76 %  int p;
  77    ((la-allocate (int "n") (int "m")) (int "p"))
  78 %  return(p);
  79 %}
  80 )
  81
  82 (defCfun "la_free_alloc(p) int p;" 0
  83 %{
  84    (la-free (int "p"))
  85 %}
  86 )
  87
  88 (defentry al (int int) (int "la_allocate"))
  89 (defentry fr (int) (int "la_free_alloc"))
  90
  91 ;;;;
  92 ;;;;                Storage Access Functions
  93 ;;;;
  94
  95 (Clines
  96 %static int get_integer(p, i)
  97 %       int p, i;
  98 %{
  99 %  return(*(((int *) p) + i));
 100 %}
 101 )
 102
 103 (defentry la-get-integer (int int) (int "get_integer"))
 104
 105 (Clines
 106 %static double get_double(p, i)
 107 %       int p, i;
 108 %{
 109 %  return(*(((double *) p) + i));
 110 %}
 111 )
 112
 113 (defentry la-get-double (int int) (double "get_double"))
 114
 115 (Clines
 116 %static double get_complex_real(p, i)
 117 %       int p, i;
 118 %{
 119 %  Complex *c = ((Complex *) p) + i;
 120 %  return(c->real);
 121 %}
 122 )
 123
 124 (defentry la-get-complex-real (int int) (double "get_complex_real"))
 125
 126 (Clines
 127 %static double get_complex_imag(p, i)
 128 %       int p, i;
 129 %{
 130 %  Complex *c = ((Complex *) p) + i;
 131 %  return(c->imag);
 132 %}
 133 )
 134
 135 (defentry la-get-complex-imag (int int) (double "get_complex_imag"))
 136
 137 (defun la-get-complex (p i)
 138   (complex (la-get-complex-real p i) (la-get-complex-imag p i)))
 139
 140 (defun la-get-pointer (p i) (la-get-integer p i))
 141
 142 ;;;;
 143 ;;;;                Storage Mutation Functions
 144 ;;;;
 145
 146 (Clines
 147 %static int put_integer(p, i, x)
 148 %       int p, i, x;
 149 %{
 150 %  *(((int *) p) + i) = x;
 151 %  return(0);
 152 %}
 153 )
 154
 155 (defentry la-put-integer (int int int) (int "put_integer"))
 156
 157 (Clines
 158 %static int put_double(p, i, x)
 159 %       int p, i;
 160 %       double x;
 161 %{
 162 %  *(((double *) p) + i) = x;
 163 %  return(0);
 164 %}
 165 )
 166
 167 (defentry la-put-double (int int double) (int "put_double"))
 168
 169 (Clines
 170 %static int put_complex(p, i, x, y)
 171 %       int p, i;
 172 %       double x, y;
 173 %{
 174 %  Complex *c = ((Complex *) p) + i;
 175 %  c->real = x;
 176 %  c->imag = y;
 177 %  return(0);
 178 %}
 179 )
 180
 181 (defentry la-put-complex (int int double double) (int "put_complex"))
 182
 183 (defun la-put-pointer (p i x) (la-put-integer p i x))
 184
 185 ;;;;
 186 ;;;; XLISP internal error message emulation
 187 ;;;;
 188
 189 (defvar *buf* (make-string 1000))
 190
 191 (defun set-buf-char (i c) (setf (elt *buf* i) (code-char c)))
 192
 193 (defun get-buf (&optional (n (position (code-char 0) *buf*)))
 194   (subseq *buf* 0 n))
 195
 196 (Clines
 197 %static int bufpos = 0;
 198 %
 199 %static resetbuf() { bufpos = 0; }
 200 %
 201 )
 202
 203 (defCfun "static prbuf(s) char *s;" 0
 204 %{
 205 %  object ch;
 206 %  int i, n;
 207 %
 208 %  n = strlen(s);
 209 %  for (i = 0; i <n; i++, bufpos++) {
 210      (set-buf-char (int "bufpos") (int "(int) s[i]"))
 211 %  }
 212    (set-buf-char (int "bufpos") (int "(int) 0"))
 213 %}
 214 )
 215
 216 (defCfun "xlfail(s) char *s;" 0
 217 %{
 218 %  object buf;
 219 %
 220 %  resetbuf();
 221 %  prbuf(s);
 222    ((get-buf (int "bufpos")) "buf")
 223    (error "buf")
 224 %}
 225 )
 226
 227 (defCfun "stdputstr(s) char *s;" 0
 228 %{
 229 %  object buf;
 230 %
 231 %  resetbuf();
 232 %  prbuf(s);
 233    ((get-buf (int "bufpos")) "buf")
 234    (princ "buf")
 235 %}
 236 )
 237
 238 (Clines
 239 %bufputstr(s)
 240 %       char *s;
 241 %{
 242 %  object buf;
 243 %
 244 %  resetbuf();
 245 %  prbuf(s);
 246 %}
 247 )
 248
 249 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 250 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 251 ;;;;
 252 ;;;;             Lisp Interfaces to Linear Algebra Routines
 253 ;;;;
 254 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 255 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 256
 257 ;;;;
 258 ;;;; Cholesky Decomposition
 259 ;;;;
 260
 261 (defentry chol-decomp-front (int int int) (int "chol_decomp_front"))
 262
 263 ;;;;
 264 ;;;; LU Decomposition
 265 ;;;;
 266
 267 (defentry lu-decomp-front (int int int int int) (int "lu_decomp_front"))
 268 (defentry lu-solve-front (int int int int int) (int "lu_solve_front"))
 269 (defentry lu-inverse-front (int int int int int int) (int "lu_inverse_front"))
 270
 271 ;;;;
 272 ;;;; SV Decomposition
 273 ;;;;
 274
 275 (defentry sv-decomp-front (int int int int int) (int "sv_decomp_front"))
 276
 277 ;;;;
 278 ;;;; QR Decomposition
 279 ;;;;
 280
 281 (defentry qr-decomp-front (int int int int int int) (int "qr_decomp_front"))
 282
 283 ;;;;
 284 ;;;; Estimate of Condition Number for Lower Triangular Matrix
 285 ;;;;
 286
 287 (defentry rcondest-front (int int) (double "rcondest_front"))
 288
 289 ;;;;
 290 ;;;; Make Rotation Matrix
 291 ;;;;
 292
 293 (defentry make-rotation-front
 294   (int int int int int double)
 295   (int "make_rotation_front"))
 296
 297 ;;;;
 298 ;;;; Eigenvalues and Eigenvectors
 299 ;;;;
 300
 301 (defentry eigen-front (int int int int int) (int "eigen_front"))
 302
 303 ;;;;
 304 ;;;; Spline Interpolation
 305 ;;;;
 306
 307 (defentry la-range-to-rseq (int int int int) (int "range_to_rseq"))
 308 (defentry spline-front (int int int int int int int) (int "spline_front"))
 309
 310 ;;;;
 311 ;;;; Kernel Density Estimators and Smoothers
 312 ;;;;
 313
 314 (defentry kernel-dens-front
 315   (int int double int int int int)
 316   (int "kernel_dens_front"))
 317
 318 (defentry kernel-smooth-front
 319   (int int int double int int int int)
 320   (int "kernel_smooth_front"))
 321
 322 ;;;;
 323 ;;;; Lowess Smoother Interface
 324 ;;;;
 325
 326 (defentry base-lowess-front
 327   (int int int double int double int int int)
 328   (int "base_lowess_front"))
 329
 330 ;;;;
 331 ;;;; FFT
 332 ;;;;
 333
 334 (defentry fft-front (int int int int) (int "fft_front"))
 335
 336 ;;;;
 337 ;;;; Maximization and Numerical Derivatives
 338 ;;;;
 339
 340 (defCfun "maximize_callback(n, px, pfval, pgrad, phess, pderivs)
 341      int n, pderivs, px, pfval, pgrad, phess;" 0
 342 %{
 343    (maximize-callback (int "n")
 344                       (int "px")
 345                       (int "pfval")
 346                       (int "pgrad")
 347                       (int "phess")
 348                       (int "pderivs"))
 349 %}
 350 )
 351
 352 (defentry numgrad-front (int int int double int) (int "numgrad_front"))
 353 (defentry numhess-front (int int int int int double int) (int "numhess_front"))
 354 (defentry base-minfo-maximize (int int int int int int) (int "minfo_maximize"))
 355
 356 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 357 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 358 ;;;;
 359 ;;;;                       Probability Distributions
 360 ;;;;
 361 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 362 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 363
 364 (Clines
 365 % extern double unirand(), gamma();
 366 % extern double normalcdf(), normalquant(), normaldens(), normalrand();
 367 % extern double bnormcdf();
 368 % extern double cauchycdf(), cauchyquant(), cauchydens(), cauchyrand();
 369 % extern double gammacdf(), gammaquant(), gammadens(), gammarand();
 370 % extern double chisqcdf(), chisqquant(), chisqdens(), chisqrand();
 371 % extern double betacdf(), betaquant(), betadens(), betarand();
 372 % extern double tcdf(), tquant(), tdens(), trand();
 373 % extern double fcdf(), fquant(), fdens(), frand();
 374 % extern double poissoncdf(), poissonpmf();
 375 % extern int poissonquant(), poissonrand();
 376 % extern double binomialcdf(), binomialpmf();
 377 % extern int binomialquant(), binomialrand();
 378 )
 379
 380 (defCfun "double uni()" 0
 381 %{
 382 %  double x;
 383    ((random (double "1.0")) (double "x"))
 384 %  return(x);
 385 %}
 386 )
 387
 388 (defentry one-uniform-rand () (double "unirand"))
 389 (defentry base-log-gamma (double) (double "gamma"))
 390
 391 (defentry base-normal-cdf (double) (double "normalcdf"))
 392 (defentry base-normal-quant (double) (double "normalquant"))
 393 (defentry base-normal-dens (double) (double "normaldens"))
 394 (defentry one-normal-rand () (double "normalrand"))
 395 (defentry base-bivnorm-cdf (double double double) (double "bnormcdf"))
 396
 397 (defentry base-cauchy-cdf (double) (double "cauchycdf"))
 398 (defentry base-cauchy-quant (double) (double "cauchyquant"))
 399 (defentry base-cauchy-dens (double) (double "cauchydens"))
 400 (defentry one-cauchy-rand () (double "cauchyrand"))
 401
 402 (defentry base-gamma-cdf (double double) (double "gammacdf"))
 403 (defentry base-gamma-quant (double double) (double "gammaquant"))
 404 (defentry base-gamma-dens (double double) (double "gammadens"))
 405 (defentry one-gamma-rand (double) (double "gammarand"))
 406
 407 (defentry base-chisq-cdf (double double) (double "chisqcdf"))
 408 (defentry base-chisq-quant (double double) (double "chisqquant"))
 409 (defentry base-chisq-dens (double double) (double "chisqdens"))
 410 (defentry one-chisq-rand (double) (double "chisqrand"))
 411
 412 (defentry base-beta-cdf (double double double) (double "betacdf"))
 413 (defentry base-beta-quant (double double double) (double "betaquant"))
 414 (defentry base-beta-dens (double double double) (double "betadens"))
 415 (defentry one-beta-rand (double double) (double "betarand"))
 416
 417 (defentry base-t-cdf (double double) (double "tcdf"))
 418 (defentry base-t-quant (double double) (double "tquant"))
 419 (defentry base-t-dens (double double) (double "tdens"))
 420 (defentry one-t-rand (double) (double "trand"))
 421
 422 (defentry base-f-cdf (double double double) (double "fcdf"))
 423 (defentry base-f-quant (double double double) (double "fquant"))
 424 (defentry base-f-dens (double double double) (double "fdens"))
 425 (defentry one-f-rand (double double) (double "frand"))
 426
 427 (defentry base-poisson-cdf (double double) (double "poissoncdf"))
 428 (defentry base-poisson-quant (double double) (int "poissonquant"))
 429 (defentry base-poisson-pmf (int double) (double "poissonpmf"))
 430 (defentry one-poisson-rand (double) (int "poissonrand"))
 431
 432 (defentry base-binomial-cdf (double int double) (double "binomialcdf"))
 433 (defentry base-binomial-quant (double int double) (int "binomialquant"))
 434 (defentry base-binomial-pmf (int int double) (double "binomialpmf"))
 435 (defentry one-binomial-rand (int double) (int "binomialrand"))
 436