task.c 21 KB
Newer Older
1
/*******************************************************************************
2
 * This file is part of SWIFT.
3
 * Copyright (c) 2012 Pedro Gonnet (pedro.gonnet@durham.ac.uk)
4
5
6
7
 *                    Matthieu Schaller (matthieu.schaller@durham.ac.uk)
 *               2015 Peter W. Draper (p.w.draper@durham.ac.uk)
 *               2016 John A. Regan (john.a.regan@durham.ac.uk)
 *                    Tom Theuns (tom.theuns@durham.ac.uk)
8
 *
9
10
11
12
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU Lesser General Public License as published
 * by the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
13
 *
14
15
16
17
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
18
 *
19
20
 * You should have received a copy of the GNU Lesser General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
21
 *
22
23
24
25
26
27
28
29
30
 ******************************************************************************/

/* Config parameters. */
#include "../config.h"

/* Some standard headers. */
#include <float.h>
#include <limits.h>
#include <sched.h>
31
32
33
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
34

35
36
/* MPI headers. */
#ifdef WITH_MPI
37
#include <mpi.h>
38
39
#endif

40
41
42
/* This object's header. */
#include "task.h"

43
/* Local headers. */
Pedro Gonnet's avatar
Pedro Gonnet committed
44
#include "atomic.h"
45
#include "engine.h"
46
#include "error.h"
47
#include "inline.h"
48
#include "lock.h"
49
50

/* Task type names. */
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
const char *taskID_names[task_type_count] = {"none",
                                             "sort",
                                             "self",
                                             "pair",
                                             "sub_self",
                                             "sub_pair",
                                             "init_grav",
                                             "init_grav_out",
                                             "ghost_in",
                                             "ghost",
                                             "ghost_out",
                                             "extra_ghost",
                                             "drift_part",
                                             "drift_gpart",
                                             "drift_gpart_out",
                                             "end_force",
                                             "kick1",
                                             "kick2",
                                             "timestep",
                                             "send",
                                             "recv",
                                             "grav_long_range",
                                             "grav_mm",
                                             "grav_down_in",
                                             "grav_down",
                                             "grav_mesh",
                                             "cooling",
                                             "star_formation",
                                             "sourceterms",
                                             "logger",
                                             "stars_ghost_in",
                                             "stars_ghost",
                                             "stars_ghost_out"};
84

85
/* Sub-task type names. */
86
const char *subtaskID_names[task_subtype_count] = {
87
88
89
    "none",          "density", "gradient",     "force", "grav",
    "external_grav", "tend",    "xv",           "rho",   "gpart",
    "multipole",     "spart",   "stars_density"};
90

91
92
93
94
95
#ifdef WITH_MPI
/* MPI communicators for the subtypes. */
MPI_Comm subtaskMPI_comms[task_subtype_count];
#endif

96
97
/**
 * @brief Computes the overlap between the parts array of two given cells.
98
 *
Matthieu Schaller's avatar
Matthieu Schaller committed
99
100
101
 * @param TYPE is the type of parts (e.g. #part, #gpart, #spart)
 * @param ARRAY is the array of this specific type.
 * @param COUNT is the number of elements in the array.
102
 */
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
#define TASK_CELL_OVERLAP(TYPE, ARRAY, COUNT)                               \
  __attribute__((always_inline))                                            \
      INLINE static size_t task_cell_overlap_##TYPE(                        \
          const struct cell *restrict ci, const struct cell *restrict cj) { \
                                                                            \
    if (ci == NULL || cj == NULL) return 0;                                 \
                                                                            \
    if (ci->ARRAY <= cj->ARRAY &&                                           \
        ci->ARRAY + ci->COUNT >= cj->ARRAY + cj->COUNT) {                   \
      return cj->COUNT;                                                     \
    } else if (cj->ARRAY <= ci->ARRAY &&                                    \
               cj->ARRAY + cj->COUNT >= ci->ARRAY + ci->COUNT) {            \
      return ci->COUNT;                                                     \
    }                                                                       \
                                                                            \
    return 0;                                                               \
  }
120

121
TASK_CELL_OVERLAP(part, hydro.parts, hydro.count);
122
123
TASK_CELL_OVERLAP(gpart, grav.parts, grav.count);
TASK_CELL_OVERLAP(spart, stars.parts, stars.count);
Loic Hausammann's avatar
Loic Hausammann committed
124

125
126
127
128
129
/**
 * @brief Returns the #task_actions for a given task.
 *
 * @param t The #task.
 */
130
131
__attribute__((always_inline)) INLINE static enum task_actions task_acts_on(
    const struct task *t) {
132
133
134
135
136
137
138

  switch (t->type) {

    case task_type_none:
      return task_action_none;
      break;

139
    case task_type_drift_part:
140
141
    case task_type_sort:
    case task_type_ghost:
142
    case task_type_extra_ghost:
Stefan Arridge's avatar
Stefan Arridge committed
143
    case task_type_cooling:
144
    case task_type_sourceterms:
145
146
147
      return task_action_part;
      break;

148
149
150
    case task_type_star_formation:
      return task_action_all;

151
    case task_type_stars_ghost:
Loic Hausammann's avatar
Loic Hausammann committed
152
153
154
      return task_action_spart;
      break;

155
156
157
158
159
160
161
    case task_type_self:
    case task_type_pair:
    case task_type_sub_self:
    case task_type_sub_pair:
      switch (t->subtype) {

        case task_subtype_density:
162
        case task_subtype_gradient:
163
164
165
166
        case task_subtype_force:
          return task_action_part;
          break;

167
        case task_subtype_stars_density:
168
169
          return task_action_all;
          break;
170

171
        case task_subtype_grav:
172
        case task_subtype_external_grav:
173
174
175
176
177
178
179
180
181
182
          return task_action_gpart;
          break;

        default:
          error("Unknow task_action for task");
          return task_action_none;
          break;
      }
      break;

183
    case task_type_end_force:
184
185
    case task_type_kick1:
    case task_type_kick2:
Loikki's avatar
Loikki committed
186
    case task_type_logger:
187
    case task_type_timestep:
188
189
    case task_type_send:
    case task_type_recv:
190
      if (t->ci->hydro.count > 0 && t->ci->grav.count > 0)
191
        return task_action_all;
192
      else if (t->ci->hydro.count > 0)
193
        return task_action_part;
194
      else if (t->ci->grav.count > 0)
195
196
197
        return task_action_gpart;
      else
        error("Task without particles");
198
199
      break;

200
    case task_type_init_grav:
201
    case task_type_grav_mm:
202
    case task_type_grav_long_range:
203
204
205
      return task_action_multipole;
      break;

206
    case task_type_drift_gpart:
207
    case task_type_grav_down:
208
    case task_type_grav_mesh:
209
      return task_action_gpart;
210
      break;
211

212
    default:
213
      error("Unknown task_action for task");
214
215
216
      return task_action_none;
      break;
  }
217

218
  /* Silence compiler warnings */
219
220
  error("Unknown task_action for task");
  return task_action_none;
221
222
}

223
224
225
226
227
228
229
/**
 * @brief Compute the Jaccard similarity of the data used by two
 *        different tasks.
 *
 * @param ta The first #task.
 * @param tb The second #task.
 */
230
231
float task_overlap(const struct task *restrict ta,
                   const struct task *restrict tb) {
232
233
234
235
236
237

  if (ta == NULL || tb == NULL) return 0.f;

  const enum task_actions ta_act = task_acts_on(ta);
  const enum task_actions tb_act = task_acts_on(tb);

238
239
  /* First check if any of the two tasks are of a type that don't
     use cells. */
240
241
242
243
244
  if (ta_act == task_action_none || tb_act == task_action_none) return 0.f;

  const int ta_part = (ta_act == task_action_part || ta_act == task_action_all);
  const int ta_gpart =
      (ta_act == task_action_gpart || ta_act == task_action_all);
245
246
  const int ta_spart =
      (ta_act == task_action_spart || ta_act == task_action_all);
247
248
249
  const int tb_part = (tb_act == task_action_part || tb_act == task_action_all);
  const int tb_gpart =
      (tb_act == task_action_gpart || tb_act == task_action_all);
250
251
  const int tb_spart =
      (tb_act == task_action_spart || tb_act == task_action_all);
252
253
254
255
256
257

  /* In the case where both tasks act on parts */
  if (ta_part && tb_part) {

    /* Compute the union of the cell data. */
    size_t size_union = 0;
258
259
260
261
    if (ta->ci != NULL) size_union += ta->ci->hydro.count;
    if (ta->cj != NULL) size_union += ta->cj->hydro.count;
    if (tb->ci != NULL) size_union += tb->ci->hydro.count;
    if (tb->cj != NULL) size_union += tb->cj->hydro.count;
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276

    /* Compute the intersection of the cell data. */
    const size_t size_intersect = task_cell_overlap_part(ta->ci, tb->ci) +
                                  task_cell_overlap_part(ta->ci, tb->cj) +
                                  task_cell_overlap_part(ta->cj, tb->ci) +
                                  task_cell_overlap_part(ta->cj, tb->cj);

    return ((float)size_intersect) / (size_union - size_intersect);
  }

  /* In the case where both tasks act on gparts */
  else if (ta_gpart && tb_gpart) {

    /* Compute the union of the cell data. */
    size_t size_union = 0;
277
278
279
280
    if (ta->ci != NULL) size_union += ta->ci->grav.count;
    if (ta->cj != NULL) size_union += ta->cj->grav.count;
    if (tb->ci != NULL) size_union += tb->ci->grav.count;
    if (tb->cj != NULL) size_union += tb->cj->grav.count;
281
282
283
284
285
286
287
288
289

    /* Compute the intersection of the cell data. */
    const size_t size_intersect = task_cell_overlap_gpart(ta->ci, tb->ci) +
                                  task_cell_overlap_gpart(ta->ci, tb->cj) +
                                  task_cell_overlap_gpart(ta->cj, tb->ci) +
                                  task_cell_overlap_gpart(ta->cj, tb->cj);

    return ((float)size_intersect) / (size_union - size_intersect);
  }
290

Loic Hausammann's avatar
Loic Hausammann committed
291
292
293
294
295
  /* In the case where both tasks act on sparts */
  else if (ta_spart && tb_spart) {

    /* Compute the union of the cell data. */
    size_t size_union = 0;
296
297
298
299
    if (ta->ci != NULL) size_union += ta->ci->stars.count;
    if (ta->cj != NULL) size_union += ta->cj->stars.count;
    if (tb->ci != NULL) size_union += tb->ci->stars.count;
    if (tb->cj != NULL) size_union += tb->cj->stars.count;
Loic Hausammann's avatar
Loic Hausammann committed
300
301
302
303
304
305
306
307
308
309

    /* Compute the intersection of the cell data. */
    const size_t size_intersect = task_cell_overlap_spart(ta->ci, tb->ci) +
                                  task_cell_overlap_spart(ta->ci, tb->cj) +
                                  task_cell_overlap_spart(ta->cj, tb->ci) +
                                  task_cell_overlap_spart(ta->cj, tb->cj);

    return ((float)size_intersect) / (size_union - size_intersect);
  }

310
311
  /* Else, no overlap */
  return 0.f;
312
}
313

314
315
/**
 * @brief Unlock the cell held by this task.
316
 *
317
318
 * @param t The #task.
 */
319
320
void task_unlock(struct task *t) {

321
322
  const enum task_types type = t->type;
  const enum task_subtypes subtype = t->subtype;
323
324
  struct cell *ci = t->ci, *cj = t->cj;

325
  /* Act based on task type. */
326
327
  switch (type) {

328
    case task_type_end_force:
329
330
    case task_type_kick1:
    case task_type_kick2:
331
    case task_type_logger:
332
    case task_type_timestep:
333
334
335
      cell_unlocktree(ci);
      cell_gunlocktree(ci);
      break;
Matthieu Schaller's avatar
Matthieu Schaller committed
336

337
    case task_type_drift_part:
338
    case task_type_sort:
339
340
341
      cell_unlocktree(ci);
      break;

342
    case task_type_drift_gpart:
343
    case task_type_grav_mesh:
344
345
346
      cell_gunlocktree(ci);
      break;

347
    case task_type_self:
348
    case task_type_sub_self:
349
350
      if (subtype == task_subtype_grav) {
        cell_gunlocktree(ci);
351
        cell_munlocktree(ci);
352
353
354
      } else {
        cell_unlocktree(ci);
      }
355
      break;
356

357
    case task_type_pair:
358
    case task_type_sub_pair:
359
360
361
      if (subtype == task_subtype_grav) {
        cell_gunlocktree(ci);
        cell_gunlocktree(cj);
362
363
        cell_munlocktree(ci);
        cell_munlocktree(cj);
364
365
366
367
368
369
      } else {
        cell_unlocktree(ci);
        cell_unlocktree(cj);
      }
      break;

370
    case task_type_grav_down:
371
      cell_gunlocktree(ci);
372
373
374
      cell_munlocktree(ci);
      break;

375
    case task_type_grav_long_range:
376
      cell_munlocktree(ci);
377
      break;
378

379
380
381
382
383
    case task_type_grav_mm:
      cell_munlocktree(ci);
      cell_munlocktree(cj);
      break;

384
385
386
387
    default:
      break;
  }
}
388
389
390
391
392
393

/**
 * @brief Try to lock the cells associated with this task.
 *
 * @param t the #task.
 */
394
395
int task_lock(struct task *t) {

396
397
  const enum task_types type = t->type;
  const enum task_subtypes subtype = t->subtype;
398
  struct cell *ci = t->ci, *cj = t->cj;
399
400
401
402
#ifdef WITH_MPI
  int res = 0, err = 0;
  MPI_Status stat;
#endif
403

404
  switch (type) {
405

406
407
408
    /* Communication task? */
    case task_type_recv:
    case task_type_send:
409
#ifdef WITH_MPI
410
411
412
413
414
      /* Check the status of the MPI request. */
      if ((err = MPI_Test(&t->req, &res, &stat)) != MPI_SUCCESS) {
        char buff[MPI_MAX_ERROR_STRING];
        int len;
        MPI_Error_string(err, buff, &len);
415
416
417
418
        error(
            "Failed to test request on send/recv task (type=%s/%s tag=%lld, "
            "%s).",
            taskID_names[t->type], subtaskID_names[t->subtype], t->flags, buff);
419
420
      }
      return res;
421
#else
422
      error("SWIFT was not compiled with MPI support.");
423
#endif
424
      break;
425

426
    case task_type_end_force:
427
428
    case task_type_kick1:
    case task_type_kick2:
Loikki's avatar
Loikki committed
429
    case task_type_logger:
430
    case task_type_timestep:
431
      if (ci->hydro.hold || ci->grav.phold) return 0;
432
433
      if (cell_locktree(ci) != 0) return 0;
      if (cell_glocktree(ci) != 0) {
Matthieu Schaller's avatar
Matthieu Schaller committed
434
435
        cell_unlocktree(ci);
        return 0;
436
437
438
      }
      break;

439
    case task_type_drift_part:
440
    case task_type_sort:
441
      if (ci->hydro.hold) return 0;
442
443
      if (cell_locktree(ci) != 0) return 0;
      break;
444

445
    case task_type_drift_gpart:
446
    case task_type_grav_mesh:
447
      if (ci->grav.phold) return 0;
448
449
450
      if (cell_glocktree(ci) != 0) return 0;
      break;

451
    case task_type_self:
452
    case task_type_sub_self:
453
      if (subtype == task_subtype_grav) {
454
        /* Lock the gparts and the m-pole */
455
        if (ci->grav.phold || ci->grav.mhold) return 0;
456
457
458
459
460
461
        if (cell_glocktree(ci) != 0)
          return 0;
        else if (cell_mlocktree(ci) != 0) {
          cell_gunlocktree(ci);
          return 0;
        }
462
463
464
465
      } else {
        if (cell_locktree(ci) != 0) return 0;
      }
      break;
466

467
    case task_type_pair:
468
    case task_type_sub_pair:
469
      if (subtype == task_subtype_grav) {
470
        /* Lock the gparts and the m-pole in both cells */
471
        if (ci->grav.phold || cj->grav.phold) return 0;
472
473
474
475
        if (cell_glocktree(ci) != 0) return 0;
        if (cell_glocktree(cj) != 0) {
          cell_gunlocktree(ci);
          return 0;
476
477
478
479
480
481
482
483
484
        } else if (cell_mlocktree(ci) != 0) {
          cell_gunlocktree(ci);
          cell_gunlocktree(cj);
          return 0;
        } else if (cell_mlocktree(cj) != 0) {
          cell_gunlocktree(ci);
          cell_gunlocktree(cj);
          cell_munlocktree(ci);
          return 0;
485
486
        }
      } else {
487
        /* Lock the parts in both cells */
488
        if (ci->hydro.hold || cj->hydro.hold) return 0;
489
490
491
492
493
494
495
        if (cell_locktree(ci) != 0) return 0;
        if (cell_locktree(cj) != 0) {
          cell_unlocktree(ci);
          return 0;
        }
      }
      break;
496

497
498
    case task_type_grav_down:
      /* Lock the gparts and the m-poles */
499
      if (ci->grav.phold || ci->grav.mhold) return 0;
500
501
502
503
504
505
506
507
      if (cell_glocktree(ci) != 0)
        return 0;
      else if (cell_mlocktree(ci) != 0) {
        cell_gunlocktree(ci);
        return 0;
      }
      break;

508
    case task_type_grav_long_range:
509
      /* Lock the m-poles */
510
      if (ci->grav.mhold) return 0;
511
      if (cell_mlocktree(ci) != 0) return 0;
Matthieu Schaller's avatar
Matthieu Schaller committed
512
513
      break;

514
515
    case task_type_grav_mm:
      /* Lock both m-poles */
516
      if (ci->grav.mhold || cj->grav.mhold) return 0;
517
518
519
520
521
522
      if (cell_mlocktree(ci) != 0) return 0;
      if (cell_mlocktree(cj) != 0) {
        cell_munlocktree(ci);
        return 0;
      }

523
524
    default:
      break;
525
526
527
528
529
  }

  /* If we made it this far, we've got a lock. */
  return 1;
}
530

531
532
533
534
535
536
537
538
539
540
541
/**
 * @brief Print basic information about a task.
 *
 * @param t The #task.
 */
void task_print(const struct task *t) {

  message("Type:'%s' sub_type:'%s' wait=%d nr_unlocks=%d skip=%d",
          taskID_names[t->type], subtaskID_names[t->subtype], t->wait,
          t->nr_unlock_tasks, t->skip);
}
542
543
544
545
546
547
548
549
550
551
552

#ifdef WITH_MPI
/**
 * @brief Create global communicators for each of the subtasks.
 */
void task_create_mpi_comms(void) {
  for (int i = 0; i < task_subtype_count; i++) {
    MPI_Comm_dup(MPI_COMM_WORLD, &subtaskMPI_comms[i]);
  }
}
#endif
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658

/**
 * @brief dump all the tasks of all the known engines into a file for postprocessing.
 *
 * Dumps the information to a file "thread_info-stepn.dat" where n is the
 * given step value, or "thread_info_MPI-stepn.dat", if we are running
 * under MPI. Note if running under MPIU all the ranks are dumped into this
 * one file, which has an additional field to identify the rank.
 *
 * @param e the #engine
 * @param step the current step.
 */
void task_dump_all(struct engine *e, int step) {

#ifdef SWIFT_DEBUG_TASKS

  /* Need this to convert ticks to seconds. */
  unsigned long long cpufreq = clocks_get_cpufreq();

#ifdef WITH_MPI
  /* Make sure output file is empty, only on one rank. */
  char dumpfile[35];
  snprintf(dumpfile, sizeof(dumpfile), "thread_info_MPI-step%d.dat", step);
  FILE *file_thread;
  if (engine_rank == 0) {
    file_thread = fopen(dumpfile, "w");
    fclose(file_thread);
  }
  MPI_Barrier(MPI_COMM_WORLD);

  for (int i = 0; i < e->nr_nodes; i++) {

    /* Rank 0 decides the index of the writing node, this happens
     * one-by-one. */
    int kk = i;
    MPI_Bcast(&kk, 1, MPI_INT, 0, MPI_COMM_WORLD);

    if (i == engine_rank) {

      /* Open file and position at end. */
      file_thread = fopen(dumpfile, "a");

      /* Add some information to help with the plots and conversion of ticks to
       * seconds. */
      fprintf(file_thread, " %03d 0 0 0 0 %lld %lld %lld %lld %lld 0 0 %lld\n",
              engine_rank, e->tic_step, e->toc_step, e->updates, e->g_updates,
              e->s_updates, cpufreq);
      int count = 0;
      for (int l = 0; l < e->sched.nr_tasks; l++) {
        if (!e->sched.tasks[l].implicit && e->sched.tasks[l].toc != 0) {
          fprintf(
              file_thread, " %03i %i %i %i %i %lli %lli %i %i %i %i %lli %i\n",
              engine_rank, e->sched.tasks[l].rid, e->sched.tasks[l].type,
              e->sched.tasks[l].subtype, (e->sched.tasks[l].cj == NULL),
              e->sched.tasks[l].tic, e->sched.tasks[l].toc,
              (e->sched.tasks[l].ci != NULL) ? e->sched.tasks[l].ci->hydro.count
                                             : 0,
              (e->sched.tasks[l].cj != NULL) ? e->sched.tasks[l].cj->hydro.count
                                             : 0,
              (e->sched.tasks[l].ci != NULL) ? e->sched.tasks[l].ci->grav.count
                                             : 0,
              (e->sched.tasks[l].cj != NULL) ? e->sched.tasks[l].cj->grav.count
                                             : 0,
              e->sched.tasks[l].flags, e->sched.tasks[l].sid);
        }
        count++;
      }
      fclose(file_thread);
    }

    /* And we wait for all to synchronize. */
    MPI_Barrier(MPI_COMM_WORLD);
  }

#else
  /* Non-MPI, so just a single engine's worth of tasks to dump. */
  char dumpfile[32];
  snprintf(dumpfile, sizeof(dumpfile), "thread_info-step%d.dat", step);
  FILE *file_thread;
  file_thread = fopen(dumpfile, "w");

  /* Add some information to help with the plots and conversion of ticks to
   * seconds. */
  fprintf(file_thread, " %d %d %d %d %lld %lld %lld %lld %lld %d %lld\n", -2,
          -1, -1, 1, e->tic_step, e->toc_step, e->updates, e->g_updates,
          e->s_updates, 0, cpufreq);
  for (int l = 0; l < e->sched.nr_tasks; l++) {
    if (!e->sched.tasks[l].implicit && e->sched.tasks[l].toc != 0) {
      fprintf(
          file_thread, " %i %i %i %i %lli %lli %i %i %i %i %i\n",
          e->sched.tasks[l].rid, e->sched.tasks[l].type,
          e->sched.tasks[l].subtype, (e->sched.tasks[l].cj == NULL),
          e->sched.tasks[l].tic, e->sched.tasks[l].toc,
          (e->sched.tasks[l].ci == NULL) ? 0
                                         : e->sched.tasks[l].ci->hydro.count,
          (e->sched.tasks[l].cj == NULL) ? 0
                                         : e->sched.tasks[l].cj->hydro.count,
          (e->sched.tasks[l].ci == NULL) ? 0 : e->sched.tasks[l].ci->grav.count,
          (e->sched.tasks[l].cj == NULL) ? 0 : e->sched.tasks[l].cj->grav.count,
          e->sched.tasks[l].sid);
    }
  }
  fclose(file_thread);
#endif  // WITH_MPI
#endif // SWIFT_DEBUG_TASKS
}