| 165 |  | double  d; | 
| 166 |  | int     p[4]; | 
| 167 |  | int     i, j; | 
| 168 | – |  | 
| 169 | – | if ((x0 == x1) & (y0 == y1)) | 
| 170 | – | return(0.); | 
| 168 |  | /* check radius */ | 
| 169 |  | p[0] = x0; p[1] = y0; p[2] = x1; p[3] = y1; | 
| 170 |  | for (i = 4; i--; ) { | 
| 383 |  | #undef discount | 
| 384 |  | } | 
| 385 |  |  | 
| 386 | + | #ifdef DUMP_MATRIX | 
| 387 | + | /* Dump transport plan and corresponding price matrix to a text file */ | 
| 388 | + | static void | 
| 389 | + | dump_matrix(const MIGRATION *me, const PRICEMAT *pm) | 
| 390 | + | { | 
| 391 | + | char    fname[256]; | 
| 392 | + | FILE    *fp; | 
| 393 | + | int     i, j; | 
| 394 | + |  | 
| 395 | + | sprintf(fname, "edge_%d-%d.txt", me->rbfv[0]->ord, me->rbfv[1]->ord); | 
| 396 | + | if ((fp = fopen(fname, "w")) == NULL) | 
| 397 | + | return; | 
| 398 | + | for (j = 0; j < 2; j++) { | 
| 399 | + | fprintf(fp, "Available from %d source RBF lobes in node %d:\n", | 
| 400 | + | me->rbfv[j]->nrbf, me->rbfv[j]->ord); | 
| 401 | + | for (i = 0; i < me->rbfv[j]->nrbf; i++) | 
| 402 | + | fprintf(fp, " %.4e", rbf_volume(&me->rbfv[j]->rbfa[i]) / | 
| 403 | + | me->rbfv[j]->vtotal); | 
| 404 | + | fputc('\n', fp); | 
| 405 | + | } | 
| 406 | + | fprintf(fp, "Price (quadratic distance metric) matrix:\n"); | 
| 407 | + | for (i = 0; i < pm->nrows; i++) { | 
| 408 | + | for (j = 0; j < pm->ncols; j++) | 
| 409 | + | fprintf(fp, " %.4e", pricerow(pm,i)[j]); | 
| 410 | + | fputc('\n', fp); | 
| 411 | + | } | 
| 412 | + | fprintf(fp, "Solution matrix (transport plan):\n"); | 
| 413 | + | for (i = 0; i < mtx_nrows(me); i++) { | 
| 414 | + | for (j = 0; j < mtx_ncols(me); j++) | 
| 415 | + | fprintf(fp, " %.4e", mtx_coef(me,i,j)); | 
| 416 | + | fputc('\n', fp); | 
| 417 | + | } | 
| 418 | + | fclose(fp); | 
| 419 | + | } | 
| 420 | + | #endif | 
| 421 | + |  | 
| 422 |  | /* Compute and insert migration along directed edge (may fork child) */ | 
| 423 |  | static MIGRATION * | 
| 424 |  | create_migration(RBFNODE *from_rbf, RBFNODE *to_rbf) | 
| 474 |  | for (j = to_rbf->nrbf; j--; ) | 
| 475 |  | mtx_coef(newmig,i,j) *= nf;     /* row now sums to 1.0 */ | 
| 476 |  | } | 
| 477 | + | #ifdef DUMP_MATRIX | 
| 478 | + | dump_matrix(newmig, &pmtx); | 
| 479 | + | #endif | 
| 480 |  | end_subprocess();                       /* exit here if subprocess */ | 
| 481 |  | free_routes(&pmtx);                     /* free working arrays */ | 
| 482 |  | free(src_rem); | 
| 592 |  | static void | 
| 593 |  | check_normal_incidence(void) | 
| 594 |  | { | 
| 595 | < | static const FVECT      norm_vec = {.0, .0, 1.}; | 
| 595 | > | static FVECT            norm_vec = {.0, .0, 1.}; | 
| 596 |  | const int               saved_nprocs = nprocs; | 
| 597 |  | RBFNODE                 *near_rbf, *mir_rbf, *rbf; | 
| 598 |  | double                  bestd; | 
| 599 |  | int                     n; | 
| 600 |  |  | 
| 565 | – |  | 
| 601 |  | if (dsf_list == NULL) | 
| 602 |  | return;                         /* XXX should be error? */ | 
| 603 |  | near_rbf = dsf_list; | 
| 647 |  | nprocs = 1;                             /* compute migration matrix */ | 
| 648 |  | if (create_migration(mir_rbf, near_rbf) == NULL) | 
| 649 |  | exit(1);                        /* XXX should never happen! */ | 
| 650 | < | /* interpolate normal dist. */ | 
| 650 | > | norm_vec[2] = input_orient;             /* interpolate normal dist. */ | 
| 651 |  | rbf = e_advect_rbf(mig_list, norm_vec, 2*near_rbf->nrbf); | 
| 652 |  | nprocs = saved_nprocs;                  /* final clean-up */ | 
| 653 |  | free(mir_rbf); |