cgit.sxemacs.org Git - sxemacs/blob - modules/ase/ase-metric.c

   1 /*
   2   ase-metric.c -- Metrical Spaces and Distances
   3   Copyright (C) 2006, 2007, 2008 Sebastian Freundt
   4
   5   Author:  Sebastian Freundt <hroptatyr@sxemacs.org>
   6
   7   * This file is part of SXEmacs.
   8   *
   9   * Redistribution and use in source and binary forms, with or without
  10   * modification, are permitted provided that the following conditions
  11   * are met:
  12   *
  13   * 1. Redistributions of source code must retain the above copyright
  14   *    notice, this list of conditions and the following disclaimer.
  15   *
  16   * 2. Redistributions in binary form must reproduce the above copyright
  17   *    notice, this list of conditions and the following disclaimer in the
  18   *    documentation and/or other materials provided with the distribution.
  19   *
  20   * 3. Neither the name of the author nor the names of any contributors
  21   *    may be used to endorse or promote products derived from this
  22   *    software without specific prior written permission.
  23   *
  24   * THIS SOFTWARE IS PROVIDED BY THE AUTHOR "AS IS" AND ANY EXPRESS OR
  25   * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
  26   * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
  27   * DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
  28   * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  29   * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  30   * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
  31   * BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
  32   * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
  33   * OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN
  34   * IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  35   */
  36
  37 /* Synched up with: Not in FSF. */
  38
  39 #include "config.h"
  40 #include "sxemacs.h"
  41 #include "ent/ent.h"
  42 #include "ase.h"
  43 #include "ase-metric.h"
  44 #include "ase-cartesian.h"
  45 #include <bytecode.h>
  46
  47 PROVIDE(ase_metric);
  48 REQUIRE(ase_metric, "ase", "ase-cartesian");
  49
  50 Lisp_Object Qase_metric, Qase_metricp;
  51 Lisp_Object Qase_euclidean_metric, Qase_euclidean_square_metric;
  52 Lisp_Object Qase_supremum_metric, Qase_trivial_metric;
  53 Lisp_Object Qase_pmetric;
  54 Lisp_Object Qmetric_distance_error;
  55
  56 \f
  57 /* stuff for the dynacat */
  58 static inline void
  59 _ase_metric_prnt(ase_metric_t n, Lisp_Object pcf)
  60 {
  61         return;
  62 }
  63
  64 static void
  65 ase_metric_prnt(Lisp_Object obj, Lisp_Object pcf, int unused)
  66 {
  67         EMOD_ASE_DEBUG_METR("m:0x%08x@0x%08x (rc:%d)\n",
  68                             (unsigned int)(XASE_METRIC(obj)),
  69                             (unsigned int)obj, 1);
  70         write_c_string("#<", pcf);
  71         print_internal(XDYNACAT_TYPE(obj), pcf, unused);
  72         {
  73                 if (NILP(XASE_METRIC_LDIST(obj))) {
  74                         write_hex_ptr(XASE_METRIC_DIST(obj),pcf);
  75                 } else {
  76                         Lisp_Object ldist = XASE_METRIC_LDIST(obj);
  77                         if (SYMBOLP(ldist)) {
  78                                 Lisp_String *name =
  79                                         symbol_name(XSYMBOL(ldist));
  80                                 write_fmt_string(pcf, " #'%s", string_data(name));
  81                         } else if (SUBRP(ldist)) {
  82                                 const char *name = subr_name(XSUBR(ldist));
  83                                 write_fmt_string(pcf, " #'%s", name);
  84                         } else {
  85                                 write_c_string(" #'(lambda ...)", pcf);
  86                         }
  87                 }
  88         }
  89         write_c_string(">", pcf);
  90         return;
  91 }
  92
  93 static void
  94 ase_metric_fini(Lisp_Object obj, int unused)
  95 {
  96         ase_metric_t m = XASE_METRIC(obj);
  97
  98         if (ase_metric_data(m)) {
  99                 xfree(ase_metric_data(m));
 100                 ase_metric_data(m) = NULL;
 101         }
 102
 103         EMOD_ASE_DEBUG_GC("m:0x%08x@0x%08x (rc:%d) shall be freed...\n",
 104                           (unsigned int)(m), (unsigned int)obj, 1);
 105
 106         return;
 107 }
 108
 109 static inline void
 110 _ase_metric_mark(ase_metric_t n)
 111 {
 112         if (n == NULL)
 113                 return;
 114
 115         mark_object(n->ldist);
 116         mark_object(n->colour);
 117         return;
 118 }
 119
 120 static void
 121 ase_metric_mark(Lisp_Object obj)
 122 {
 123         EMOD_ASE_DEBUG_METR("m:0x%08x@0x%08x (rc:%d) shall be marked...\n",
 124                             (unsigned int)(XASE_METRIC(obj)),
 125                             (unsigned int)obj, 1);
 126         _ase_metric_mark(XASE_METRIC(obj));
 127         return;
 128 }
 129
 130 \f
 131 inline Lisp_Object
 132 _ase_wrap_metric(ase_metric_t m)
 133 {
 134         Lisp_Object result;
 135
 136         result = make_dynacat(m);
 137         XDYNACAT_TYPE(result) = Qase_metric;
 138
 139         set_dynacat_printer(result, ase_metric_prnt);
 140         set_dynacat_marker(result, ase_metric_mark);
 141         set_dynacat_finaliser(result, ase_metric_fini);
 142         set_dynacat_intprinter(result, NULL);
 143
 144         EMOD_ASE_DEBUG_METR("m:0x%08x (rc:%d) shall be wrapped to 0x%08x...\n",
 145                             (unsigned int)m, 1, (unsigned int)result);
 146
 147         return result;
 148 }
 149
 150 static inline ase_metric_t
 151 _ase_make_metric(ase_distance_f fn, void *data, Lisp_Object lambda)
 152 {
 153         ase_metric_t m = NULL;
 154
 155         m = xnew(struct ase_metric_s);
 156
 157         ase_metric_dist(m) = fn;
 158         ase_metric_ldist(m) = lambda;
 159         m->colour = Qnil;
 160         ase_metric_data(m) = data;
 161
 162         EMOD_ASE_DEBUG_METR("m:0x%08x (rc:0) shall be created...\n",
 163                             (unsigned int)m);
 164         return m;
 165 }
 166
 167 Lisp_Object ase_make_metric(ase_distance_f fn, void *data, Lisp_Object lambda)
 168 {
 169         ase_metric_t m = NULL;
 170         Lisp_Object result = Qnil;
 171
 172         m = _ase_make_metric(fn, data, lambda);
 173         XSETASE_METRIC(result, m);
 174
 175         return result;
 176 }
 177
 178 \f
 179 /* some of the more common metrics */
 180 static inline Lisp_Object
 181 _ase_metric_euclidean_1dim_sq(Lisp_Object a, Lisp_Object b)
 182 {
 183         Lisp_Object tmp = ent_binop(ASE_BINARY_OP_DIFF, a, b);
 184         return ent_binop(ASE_BINARY_OP_PROD, tmp, tmp);
 185 }
 186
 187 static Lisp_Object
 188 _ase_metric_euclidean_ndim_sq(Lisp_Object a, Lisp_Object b)
 189 {
 190         int i, dim = XASE_CARTESIAN_DIMENSION(a);
 191         Lisp_Object tmp[dim];
 192         Lisp_Object *aos = XASE_CARTESIAN_OBJECTS(a);
 193         Lisp_Object *bos = XASE_CARTESIAN_OBJECTS(b);
 194
 195         for (i = 0; i < dim; i++) {
 196                 tmp[i] = _ase_metric_euclidean_1dim_sq(aos[i], bos[i]);
 197         }
 198         return Fent_binop_sum(dim, tmp);
 199 }
 200
 201 static inline Lisp_Object
 202 _ase_metric_euclidean_1dim_fast(Lisp_Object a, Lisp_Object b)
 203 {
 204         return Fabs(ent_binop(ASE_BINARY_OP_DIFF, a, b));
 205 }
 206
 207 static inline Lisp_Object
 208 _ase_metric_euclidean_1dim(Lisp_Object a, Lisp_Object b)
 209 {
 210         return Fsqrt(_ase_metric_euclidean_1dim_sq(a, b), Qnil);
 211 }
 212
 213 static inline Lisp_Object
 214 _ase_metric_euclidean_ndim(Lisp_Object a, Lisp_Object b)
 215 {
 216         return Fsqrt(_ase_metric_euclidean_ndim_sq(a, b), Qnil);
 217 }
 218
 219 static Lisp_Object
 220 ase_metric_euclidean(void *unused, Lisp_Object a, Lisp_Object b)
 221 {
 222         if (COMPARABLEP(a) && COMPARABLEP(b)) {
 223                 return _ase_metric_euclidean_1dim(a, b);
 224         } else if (ASE_CARTESIAN_INTERIOR_P(a) &&
 225                    ASE_CARTESIAN_INTERIOR_P(b) &&
 226                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(a)) &&
 227                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(b)) &&
 228                    XASE_CARTESIAN_DIMENSION(a) ==
 229                    XASE_CARTESIAN_DIMENSION(b)) {
 230                 return _ase_metric_euclidean_ndim(a, b);
 231         }
 232         ase_cartesian_embedding_error(a, b);
 233         return Qzero;
 234 }
 235
 236 static Lisp_Object
 237 ase_metric_euclidean_sq(void *unused, Lisp_Object a, Lisp_Object b)
 238 {
 239         if (COMPARABLEP(a) && COMPARABLEP(b)) {
 240                 return _ase_metric_euclidean_1dim_sq(a, b);
 241         } else if (ASE_CARTESIAN_INTERIOR_P(a) &&
 242                    ASE_CARTESIAN_INTERIOR_P(b) &&
 243                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(a)) &&
 244                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(b)) &&
 245                    XASE_CARTESIAN_DIMENSION(a) ==
 246                    XASE_CARTESIAN_DIMENSION(b)) {
 247                 return _ase_metric_euclidean_ndim_sq(a, b);
 248         }
 249         ase_cartesian_embedding_error(a, b);
 250         return Qzero;
 251 }
 252
 253 static inline Lisp_Object
 254 _ase_metric_supremum_ndim(Lisp_Object a, Lisp_Object b)
 255 {
 256         int i, dim = XASE_CARTESIAN_DIMENSION(a);
 257         Lisp_Object sup = Qzero;
 258         Lisp_Object *aos = XASE_CARTESIAN_OBJECTS(a);
 259         Lisp_Object *bos = XASE_CARTESIAN_OBJECTS(b);
 260
 261         for (i = 0; i < dim; i++) {
 262                 Lisp_Object tmp =
 263                         _ase_metric_euclidean_1dim_fast(aos[i], bos[i]);
 264                 if (ent_binrel(ASE_BINARY_REL_LESSP, sup, tmp)) {
 265                         sup = tmp;
 266                 }
 267         }
 268         return sup;
 269 }
 270
 271 static Lisp_Object
 272 ase_metric_supremum(void *unused, Lisp_Object a, Lisp_Object b)
 273 {
 274         if (COMPARABLEP(a) && COMPARABLEP(b)) {
 275                 return _ase_metric_euclidean_1dim_fast(a, b);
 276         } else if (ASE_CARTESIAN_INTERIOR_P(a) &&
 277                    ASE_CARTESIAN_INTERIOR_P(b) &&
 278                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(a)) &&
 279                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(b)) &&
 280                    XASE_CARTESIAN_DIMENSION(a) ==
 281                    XASE_CARTESIAN_DIMENSION(b)) {
 282                 return _ase_metric_supremum_ndim(a, b);
 283         }
 284         ase_cartesian_embedding_error(a, b);
 285         return Qzero;
 286 }
 287
 288 static inline Lisp_Object
 289 _ase_metric_trivial_1dim(Lisp_Object a, Lisp_Object b)
 290 {
 291         if (!ent_binrel(ASE_BINARY_REL_EQUALP, a, b))
 292                 return make_int(1);
 293         else
 294                 return Qzero;
 295 }
 296
 297 static inline Lisp_Object
 298 _ase_metric_trivial_ndim(Lisp_Object a, Lisp_Object b)
 299 {
 300         int i, dim = XASE_CARTESIAN_DIMENSION(a);
 301         Lisp_Object *aos = XASE_CARTESIAN_OBJECTS(a);
 302         Lisp_Object *bos = XASE_CARTESIAN_OBJECTS(b);
 303
 304         for (i = 0; i < dim; i++) {
 305                 if (XINT(_ase_metric_trivial_1dim(aos[i], bos[i])) == 1)
 306                         return make_int(1);
 307         }
 308         return Qzero;
 309 }
 310
 311 static Lisp_Object
 312 ase_metric_trivial(void *unused, Lisp_Object a, Lisp_Object b)
 313 {
 314         if (COMPARABLEP(a) && COMPARABLEP(b)) {
 315                 return _ase_metric_trivial_1dim(a, b);
 316         } else if (ASE_CARTESIAN_INTERIOR_P(a) &&
 317                    ASE_CARTESIAN_INTERIOR_P(b) &&
 318                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(a)) &&
 319                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(b)) &&
 320                    XASE_CARTESIAN_DIMENSION(a) ==
 321                    XASE_CARTESIAN_DIMENSION(b)) {
 322                 return _ase_metric_trivial_ndim(a, b);
 323         }
 324         ase_cartesian_embedding_error(a, b);
 325         return Qzero;
 326 }
 327
 328 static inline Lisp_Object
 329 _ase_metric_p_1dim_p(Lisp_Object a, Lisp_Object b, unsigned int p)
 330 {
 331         Lisp_Object tmp = ent_binop(ASE_BINARY_OP_DIFF, a, b);
 332         Lisp_Object result = ent_binop(ASE_BINARY_OP_POW, tmp, make_int(p));
 333         if ((p & 1) == 0)
 334                 return result;
 335         else
 336                 return Fabs(result);
 337 }
 338
 339 static Lisp_Object
 340 _ase_metric_p_ndim_p(Lisp_Object a, Lisp_Object b, unsigned int p)
 341 {
 342         int i, dim = XASE_CARTESIAN_DIMENSION(a);
 343         Lisp_Object tmp[dim];
 344         Lisp_Object *aos = XASE_CARTESIAN_OBJECTS(a);
 345         Lisp_Object *bos = XASE_CARTESIAN_OBJECTS(b);
 346
 347         for (i = 0; i < dim; i++) {
 348                 tmp[i] = _ase_metric_p_1dim_p(aos[i], bos[i], p);
 349         }
 350         return Fent_binop_sum(dim, tmp);
 351 }
 352
 353 #ifdef HAVE_MPFR
 354 static inline Lisp_Object
 355 _ase_metric_p_1dim(Lisp_Object a, Lisp_Object b, unsigned int p)
 356 {
 357         return Froot(_ase_metric_p_1dim_p(a, b, p), make_int(p), Qnil);
 358 }
 359
 360 static inline Lisp_Object
 361 _ase_metric_p_ndim(Lisp_Object a, Lisp_Object b, unsigned int p)
 362 {
 363         return Froot(_ase_metric_p_ndim_p(a, b, p), make_int(p), Qnil);
 364 }
 365
 366 static Lisp_Object
 367 ase_metric_p(void *data, Lisp_Object a, Lisp_Object b)
 368 {
 369         unsigned int p = ((ase_pmetric_data_t)data)->p;
 370         if (COMPARABLEP(a) && COMPARABLEP(b)) {
 371                 return _ase_metric_p_1dim(a, b, p);
 372         } else if (ASE_CARTESIAN_INTERIOR_P(a) &&
 373                    ASE_CARTESIAN_INTERIOR_P(b) &&
 374                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(a)) &&
 375                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(b)) &&
 376                    XASE_CARTESIAN_DIMENSION(a) ==
 377                    XASE_CARTESIAN_DIMENSION(b)) {
 378                 return _ase_metric_p_ndim(a, b, p);
 379         }
 380         ase_cartesian_embedding_error(a, b);
 381         return Qzero;
 382 }
 383 #endif
 384
 385 static Lisp_Object
 386 ase_metric_p_p(void *data, Lisp_Object a, Lisp_Object b)
 387 {
 388         unsigned int p = ((ase_pmetric_data_t)data)->p;
 389         if (COMPARABLEP(a) && COMPARABLEP(b)) {
 390                 return _ase_metric_p_1dim_p(a, b, p);
 391         } else if (ASE_CARTESIAN_INTERIOR_P(a) &&
 392                    ASE_CARTESIAN_INTERIOR_P(b) &&
 393                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(a)) &&
 394                    COMPARABLEP(XASE_CARTESIAN_FIRST_OBJECT(b)) &&
 395                    XASE_CARTESIAN_DIMENSION(a) ==
 396                    XASE_CARTESIAN_DIMENSION(b)) {
 397                 return _ase_metric_p_ndim_p(a, b, p);
 398         }
 399         ase_cartesian_embedding_error(a, b);
 400         return Qzero;
 401 }
 402
 403 \f
 404 /* ###autoload */
 405 DEFUN("ase-p-metric", Fase_p_metric, 1, 1, 0, /*
 406 Return a p-metric for some natural number P.
 407 */
 408       (p))
 409 #ifndef HAVE_MPFR
 410 {
 411         error("MPFR not available which is mandatory for p-metrics");
 412         return Qnull_pointer;
 413 }
 414 #else
 415 {
 416         ase_pmetric_data_t data;
 417         CHECK_NATNUM(p);
 418
 419         data = xnew(struct ase_pmetric_data_s);
 420         data->p = XUINT(p);
 421         return ase_make_metric(ase_metric_p, data, Qnil);
 422 }
 423 #endif
 424
 425 /* ###autoload */
 426 DEFUN("ase-p-metric*", Fase_p_metricX, 1, 1, 0, /*
 427 Return a p-metric without the final root for some natural number P.
 428 */
 429       (p))
 430 {
 431         ase_pmetric_data_t data;
 432         CHECK_NATNUM(p);
 433
 434         data = xnew(struct ase_pmetric_data_s);
 435         data->p = XUINT(p);
 436         return ase_make_metric(ase_metric_p_p, data, Qnil);
 437 }
 438
 439 /* ###autoload */
 440 DEFUN("ase-metric", Fase_metric, 1, 1, 0, /*
 441 Return a metric from a distance function FN.
 442
 443 FN should take two arguments and return the distance between those,
 444 a distance by definition lives in the reals.
 445 */
 446       (fn))
 447 {
 448         if (!SUBRP(fn) && !SYMBOLP(fn) &&
 449             !COMPILED_FUNCTIONP(fn) &&
 450             !(CONSP(fn) && EQ(XCAR(fn), Qlambda))) {
 451                 signal_invalid_function_error(fn);
 452                 return Qnil;
 453         }
 454
 455         return ase_make_metric(NULL, NULL, fn);
 456 }
 457
 458 DEFUN("ase-metric-distance", Fase_metric_distance, 3, 3, 0, /*
 459 Return the distance of P1 and P2 with respect to METRIC.
 460 */
 461       (metric, p1, p2))
 462 {
 463         ase_distance_f dist;
 464         Lisp_Object ldist;
 465         CHECK_ASE_METRIC(metric);
 466
 467         if ((dist = XASE_METRIC_DIST(metric))) {
 468                 void *data = XASE_METRIC_DATA(metric);
 469                 return dist(data, p1, p2);
 470         } else if (!NILP((ldist = XASE_METRIC_LDIST(metric)))) {
 471                 /* This portion can GC */
 472                 Lisp_Object args[3] = {ldist, p1, p2};
 473                 Lisp_Object res = Qnil;
 474                 struct gcpro ngcpro1, ngcpro2;
 475                 NGCPRO1n(res, args, countof(args));
 476                 res = Ffuncall(countof(args), args);
 477                 NUNGCPRO;
 478                 if (!NILP(Fnonnegativep(res))) {
 479                         return res;
 480                 } else {
 481                         signal_error(Qmetric_distance_error, list1(ldist));
 482                         return Qnil;
 483                 }
 484         }
 485
 486         dead_wrong_type_argument(Qase_metricp, metric);
 487         return Qnil;
 488 }
 489
 490 \f
 491 /* initialiser code */
 492 #define EMODNAME        ase_metric
 493
 494 void
 495 EMOD_PUBINIT(void)
 496 {
 497         DEFSUBR(Fase_p_metric);
 498         DEFSUBR(Fase_p_metricX);
 499         DEFSUBR(Fase_metric);
 500         DEFSUBR(Fase_metric_distance);
 501
 502         defsymbol(&Qase_metric, "ase:metric");
 503         defsymbol(&Qase_metricp, "ase:metricp");
 504
 505         DEFERROR(Qmetric_distance_error,
 506                  "Distance function must have non-negative image",
 507                  Qdomain_error);
 508
 509         DEFVAR_CONST_LISP("ase-euclidean-metric", &Qase_euclidean_metric /*
 510                                                                           */);
 511         DEFVAR_CONST_LISP("ase-euclidean-square-metric",
 512                           &Qase_euclidean_square_metric /*
 513                                                          */);
 514         DEFVAR_CONST_LISP("ase-supremum-metric", &Qase_supremum_metric /*
 515                                                                         */);
 516         DEFVAR_CONST_LISP("ase-trivial-metric", &Qase_trivial_metric /*
 517                                                                       */);
 518
 519         EMOD_PUBREINIT();
 520         Fprovide(intern("ase-metric"));
 521 }
 522
 523 void
 524 EMOD_PUBREINIT(void)
 525 {
 526         Qase_euclidean_metric =
 527                 ase_make_metric(ase_metric_euclidean, NULL, Qnil);
 528         Qase_euclidean_square_metric =
 529                 ase_make_metric(ase_metric_euclidean_sq, NULL, Qnil);
 530         Qase_supremum_metric =
 531                 ase_make_metric(ase_metric_supremum, NULL, Qnil);
 532         Qase_trivial_metric =
 533                 ase_make_metric(ase_metric_trivial, NULL, Qnil);
 534 }
 535
 536 void
 537 EMOD_PUBDEINIT(void)
 538 {
 539         Frevoke(intern("ase-metric"));
 540 }
 541
 542 /* ase-metric ends here */