main.c 52.4 KB
Newer Older
1
/*******************************************************************************
2
 * This file is part of SWIFT.
3
 * Copyright (c) 2012 Pedro Gonnet (pedro.gonnet@durham.ac.uk),
4
 *                    Matthieu Schaller (matthieu.schaller@durham.ac.uk)
5
 *               2015 Peter W. Draper (p.w.draper@durham.ac.uk)
6
7
8
 *                    Angus Lepper (angus.lepper@ed.ac.uk)
 *               2016 John A. Regan (john.a.regan@durham.ac.uk)
 *                    Tom Theuns (tom.theuns@durham.ac.uk)
9
 *
10
11
12
13
 * 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.
14
 *
15
16
17
18
 * 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.
19
 *
20
21
 * 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/>.
22
 *
23
24
 ******************************************************************************/

Pedro Gonnet's avatar
Pedro Gonnet committed
25
26
/* Config parameters. */
#include "../config.h"
Pedro Gonnet's avatar
Pedro Gonnet committed
27
28

/* Some standard headers. */
29
#include <errno.h>
30
#include <fenv.h>
31
#include <libgen.h>
Pedro Gonnet's avatar
Pedro Gonnet committed
32
33
34
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
35
#include <sys/stat.h>
36
#include <unistd.h>
Pedro Gonnet's avatar
Pedro Gonnet committed
37

38
39
/* MPI headers. */
#ifdef WITH_MPI
40
#include <mpi.h>
41
42
#endif

Pedro Gonnet's avatar
Pedro Gonnet committed
43
/* Local headers. */
44
#include "argparse.h"
45
#include "swift.h"
Pedro Gonnet's avatar
Pedro Gonnet committed
46

47
48
/* Engine policy flags. */
#ifndef ENGINE_POLICY
49
#define ENGINE_POLICY engine_policy_none
50
51
#endif

James Willis's avatar
James Willis committed
52
53
54
/* Global profiler. */
struct profiler prof;

55
/*  Usage string. */
56
static const char *const swift_usage[] = {
Peter W. Draper's avatar
Peter W. Draper committed
57
58
59
60
61
    "swift [options] [[--] param-file]",
    "swift [options] param-file",
    "swift_mpi [options] [[--] param-file]",
    "swift_mpi [options] param-file",
    NULL,
62
63
};

64
/* Function to handle multiple -P arguments. */
65
66
67
68
69
70
71
72
73
74
75
struct cmdparams {
  const char *param[PARSER_MAX_NO_OF_PARAMS];
  int nparam;
};

static int handle_cmdparam(struct argparse *self,
                           const struct argparse_option *opt) {
  struct cmdparams *cmdps = (struct cmdparams *)opt->data;
  cmdps->param[cmdps->nparam] = *(char **)opt->value;
  cmdps->nparam++;
  return 1;
76
77
}

Pedro Gonnet's avatar
Pedro Gonnet committed
78
79
80
81
/**
 * @brief Main routine that loads a few particles and generates some output.
 *
 */
82
83
int main(int argc, char *argv[]) {

84
  struct clocks_time tic, toc;
85
  struct engine e;
86

87
88
  /* Structs used by the engine. Declare now to make sure these are always in
   * scope.  */
lhausamm's avatar
lhausamm committed
89
  struct chemistry_global_data chemistry;
90
  struct cooling_function_data cooling_func;
91
  struct cosmology cosmo;
92
  struct external_potential potential;
93
  struct star_formation starform;
94
  struct pm_mesh mesh;
95
96
97
  struct gpart *gparts = NULL;
  struct gravity_props gravity_properties;
  struct hydro_props hydro_properties;
Loic Hausammann's avatar
Loic Hausammann committed
98
  struct stars_props stars_properties;
99
  struct feedback_props feedback_properties;
100
  struct entropy_floor_properties entropy_floor;
101
  struct black_holes_props black_holes_properties;
102
  struct fof_props fof_properties;
103
104
105
106
  struct part *parts = NULL;
  struct phys_const prog_const;
  struct space s;
  struct spart *sparts = NULL;
107
  struct bpart *bparts = NULL;
108
  struct unit_system us;
109

110
  int nr_nodes = 1, myrank = 0;
111

112
#ifdef WITH_MPI
113
  /* Start by initializing MPI. */
114
  int res = 0, prov = 0;
115
116
  if ((res = MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &prov)) !=
      MPI_SUCCESS)
117
118
    error("Call to MPI_Init failed with error %i.", res);
  if (prov != MPI_THREAD_MULTIPLE)
119
    error(
120
121
        "MPI does not provide the level of threading"
        " required (MPI_THREAD_MULTIPLE).");
122
  if ((res = MPI_Comm_size(MPI_COMM_WORLD, &nr_nodes)) != MPI_SUCCESS)
123
124
125
    error("MPI_Comm_size failed with error %i.", res);
  if ((res = MPI_Comm_rank(MPI_COMM_WORLD, &myrank)) != MPI_SUCCESS)
    error("Call to MPI_Comm_rank failed with error %i.", res);
126

127
  /* Make sure messages are stamped with the correct rank and step. */
128
  engine_rank = myrank;
129
  engine_current_step = 0;
130

131
132
133
  if ((res = MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_RETURN)) !=
      MPI_SUCCESS)
    error("Call to MPI_Comm_set_errhandler failed with error %i.", res);
134
135
  if (myrank == 0)
    printf("[0000] [00000.0] main: MPI is up and running with %i node(s).\n\n",
Peter W. Draper's avatar
Peter W. Draper committed
136
           nr_nodes);
137
138
139
140
  if (nr_nodes == 1) {
    message("WARNING: you are running with one MPI rank.");
    message("WARNING: you should use the non-MPI version of this program.");
  }
141
  fflush(stdout);
142

143
#endif
144

145
  /* Welcome to SWIFT, you made the right choice */
146
  if (myrank == 0) greetings(/*fof=*/0);
147

148
  int with_aff = 0;
149
  int dry_run = 0;
150
  int dump_tasks = 0;
Loic Hausammann's avatar
Loic Hausammann committed
151
  int dump_cells = 0;
152
  int dump_threadpool = 0;
153
  int nsteps = -2;
154
  int restart = 0;
155
156
  int with_cosmology = 0;
  int with_external_gravity = 0;
157
  int with_temperature = 0;
158
  int with_cooling = 0;
159
160
  int with_self_gravity = 0;
  int with_hydro = 0;
161
  int with_stars = 0;
162
  int with_fof = 0;
163
  int with_star_formation = 0;
Loic Hausammann's avatar
Loic Hausammann committed
164
  int with_feedback = 0;
165
  int with_black_holes = 0;
166
167
  int with_timestep_limiter = 0;
  int with_timestep_sync = 0;
168
  int with_fp_exceptions = 0;
169
  int with_drift_all = 0;
170
  int with_mpole_reconstruction = 0;
171
  int with_structure_finding = 0;
Loic Hausammann's avatar
Loic Hausammann committed
172
  int with_logger = 0;
173
  int with_eagle = 0;
174
  int verbose = 0;
175
  int nr_threads = 1;
176
  int with_verbose_timers = 0;
177
  char *output_parameters_filename = NULL;
178
  char *cpufreqarg = NULL;
179
  char *param_filename = NULL;
180
  char restart_file[200] = "";
181
  unsigned long long cpufreq = 0;
182
183
184
185
  struct cmdparams cmdps;
  cmdps.nparam = 0;
  cmdps.param[0] = NULL;
  char *buffer = NULL;
186

187
188
189
  /* Parse the command-line parameters. */
  struct argparse_option options[] = {
      OPT_HELP(),
190

191
192
      OPT_GROUP("  Simulation options:\n"),
      OPT_BOOLEAN('b', "feedback", &with_feedback, "Run with stars feedback.",
193
194
195
                  NULL, 0, 0),
      OPT_BOOLEAN('c', "cosmology", &with_cosmology,
                  "Run with cosmological time integration.", NULL, 0, 0),
196
197
198
      OPT_BOOLEAN(0, "temperature", &with_temperature,
                  "Run with temperature calculation.", NULL, 0, 0),
      OPT_BOOLEAN('C', "cooling", &with_cooling,
199
200
                  "Run with cooling (also switches on --temperature).", NULL, 0,
                  0),
201
202
203
204
205
      OPT_BOOLEAN('D', "drift-all", &with_drift_all,
                  "Always drift all particles even the ones far from active "
                  "particles. This emulates Gadget-[23] and GIZMO's default "
                  "behaviours.",
                  NULL, 0, 0),
206
      OPT_BOOLEAN('F', "star-formation", &with_star_formation,
207
                  "Run with star formation.", NULL, 0, 0),
208
209
210
211
212
213
      OPT_BOOLEAN('g', "external-gravity", &with_external_gravity,
                  "Run with an external gravitational potential.", NULL, 0, 0),
      OPT_BOOLEAN('G', "self-gravity", &with_self_gravity,
                  "Run with self-gravity.", NULL, 0, 0),
      OPT_BOOLEAN('M', "multipole-reconstruction", &with_mpole_reconstruction,
                  "Reconstruct the multipoles every time-step.", NULL, 0, 0),
214
215
      OPT_BOOLEAN('s', "hydro", &with_hydro, "Run with hydrodynamics.", NULL, 0,
                  0),
216
      OPT_BOOLEAN('S', "stars", &with_stars, "Run with stars.", NULL, 0, 0),
217
218
      OPT_BOOLEAN('B', "black-holes", &with_black_holes,
                  "Run with black holes.", NULL, 0, 0),
219
220
221
222
      OPT_BOOLEAN(
          'u', "fof", &with_fof,
          "Run Friends-of-Friends algorithm to perform black hole seeding.",
          NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
223
      OPT_BOOLEAN('x', "velociraptor", &with_structure_finding,
224
                  "Run with structure finding.", NULL, 0, 0),
225
226
      OPT_BOOLEAN(0, "limiter", &with_timestep_limiter,
                  "Run with time-step limiter.", NULL, 0, 0),
227
228
229
      OPT_BOOLEAN(0, "sync", &with_timestep_sync,
                  "Run with time-step synchronization of particles hit by "
                  "feedback events.",
230
                  NULL, 0, 0),
231

232
233
234
235
236
237
238
239
      OPT_GROUP("  Simulation meta-options:\n"),
      OPT_BOOLEAN(
          0, "eagle", &with_eagle,
          "Run with all the options needed for the EAGLE model. This is "
          "equivalent to --hydro --limiter --sync --self-gravity --stars "
          "--star-formation --cooling --feedback --black-holes --fof.",
          NULL, 0, 0),

240
      OPT_GROUP("  Control options:\n"),
241
      OPT_BOOLEAN('a', "pin", &with_aff,
242
                  "Pin runners using processor affinity.", NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
243
244
245
246
247
      OPT_BOOLEAN('d', "dry-run", &dry_run,
                  "Dry run. Read the parameter file, allocates memory but does "
                  "not read the particles from ICs. Exits before the start of "
                  "time integration. Checks the validity of parameters and IC "
                  "files as well as memory limits.",
248
                  NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
249
250
251
252
253
254
255
      OPT_BOOLEAN('e', "fpe", &with_fp_exceptions,
                  "Enable floating-point exceptions (debugging mode).", NULL, 0,
                  0),
      OPT_STRING('f', "cpu-frequency", &cpufreqarg,
                 "Overwrite the CPU "
                 "frequency (Hz) to be used for time measurements.",
                 NULL, 0, 0),
Loic Hausammann's avatar
Loic Hausammann committed
256
      OPT_BOOLEAN(0, "logger", &with_logger,
Loic Hausammann's avatar
Loic Hausammann committed
257
                  "Run with the particle logger.", NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
258
      OPT_INTEGER('n', "steps", &nsteps,
259
260
                  "Execute a fixed number of time steps. When unset use the "
                  "time_end parameter to stop.",
261
                  NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
262
      OPT_STRING('o', "output-params", &output_parameters_filename,
263
264
265
                 "Generate a parameter file with the options for selecting the "
                 "output fields.",
                 NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
266
      OPT_STRING('P', "param", &buffer,
267
268
                 "Set parameter value, overiding the value read from the "
                 "parameter file. Can be used more than once {sec:par:value}.",
269
270
271
                 handle_cmdparam, (intptr_t)&cmdps, 0),
      OPT_BOOLEAN('r', "restart", &restart, "Continue using restart files.",
                  NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
272
      OPT_INTEGER('t', "threads", &nr_threads,
273
274
                  "The number of threads to use on each MPI rank. Defaults to "
                  "1 if not specified.",
275
276
                  NULL, 0, 0),
      OPT_INTEGER('T', "timers", &with_verbose_timers,
Peter W. Draper's avatar
Peter W. Draper committed
277
278
279
                  "Print timers every time-step.", NULL, 0, 0),
      OPT_INTEGER('v', "verbose", &verbose,
                  "Run in verbose mode, in MPI mode 2 outputs from all ranks.",
280
                  NULL, 0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
281
282
283
      OPT_INTEGER('y', "task-dumps", &dump_tasks,
                  "Time-step frequency at which task graphs are dumped.", NULL,
                  0, 0),
Loic Hausammann's avatar
Loic Hausammann committed
284
285
286
      OPT_INTEGER(0, "cell-dumps", &dump_cells,
                  "Time-step frequency at which cell graphs are dumped.", NULL,
                  0, 0),
Peter W. Draper's avatar
Peter W. Draper committed
287
      OPT_INTEGER('Y', "threadpool-dumps", &dump_threadpool,
288
289
                  "Time-step frequency at which threadpool tasks are dumped.",
                  NULL, 0, 0),
290
      OPT_END(),
291
292
  };
  struct argparse argparse;
293
  argparse_init(&argparse, options, swift_usage, 0);
Peter W. Draper's avatar
Peter W. Draper committed
294
295
296
  argparse_describe(&argparse, "\nParameters:",
                    "\nSee the file examples/parameter_example.yml for an "
                    "example of parameter file.");
297
298
  int nargs = argparse_parse(&argparse, argc, (const char **)argv);

299
300
301
302
303
304
305
306
307
308
309
310
311
312
  /* Deal with meta options */
  if (with_eagle) {
    with_hydro = 1;
    with_timestep_limiter = 1;
    with_timestep_sync = 1;
    with_self_gravity = 1;
    with_stars = 1;
    with_star_formation = 1;
    with_cooling = 1;
    with_feedback = 1;
    with_black_holes = 1;
    with_fof = 1;
  }

313
314
315
316
317
318
319
  /* Write output parameter file */
  if (myrank == 0 && output_parameters_filename != NULL) {
    io_write_output_field_parameter(output_parameters_filename);
    printf("End of run.\n");
    return 0;
  }

320
321
  /* Need a parameter file. */
  if (nargs != 1) {
Peter W. Draper's avatar
Peter W. Draper committed
322
323
324
    if (myrank == 0) argparse_usage(&argparse);
    printf("\nError: no parameter file was supplied.\n");
    return 1;
325
326
  }
  param_filename = argv[0];
327

328
  /* Checks of options. */
Peter W. Draper's avatar
Peter W. Draper committed
329
#if !defined(HAVE_SETAFFINITY) || !defined(HAVE_LIBNUMA)
330
  if (with_aff) {
Peter W. Draper's avatar
Peter W. Draper committed
331
332
    printf("Error: no NUMA support for thread affinity\n");
    return 1;
333
  }
334
#endif
335

Loic Hausammann's avatar
Loic Hausammann committed
336
337
#if !defined(WITH_LOGGER)
  if (with_logger) {
338
    printf("Error: the particle logger is not available, please compile with --enable-logger.");
Loic Hausammann's avatar
Loic Hausammann committed
339
340
341
342
    return 1;
  }
#endif

343
#ifndef HAVE_FE_ENABLE_EXCEPT
344
  if (with_fp_exceptions) {
Peter W. Draper's avatar
Peter W. Draper committed
345
346
    printf("Error: no support for floating point exceptions\n");
    return 1;
347
  }
348
#endif
349
350

#ifndef HAVE_VELOCIRAPTOR
351
  if (with_structure_finding) {
Peter W. Draper's avatar
Peter W. Draper committed
352
353
    printf("Error: VELOCIraptor is not available\n");
    return 1;
354
  }
355
#endif
356

357
#ifndef SWIFT_DEBUG_TASKS
358
  if (dump_tasks) {
359
    if (myrank == 0) {
Peter W. Draper's avatar
Peter W. Draper committed
360
361
362
      message(
          "WARNING: complete task dumps are only created when "
          "configured with --enable-task-debugging.");
363
364
      message("         Basic task statistics will be output.");
    }
365
  }
366
#endif
367

Loic Hausammann's avatar
Loic Hausammann committed
368
369
370
371
372
373
374
375
376
377
#ifndef SWIFT_CELL_GRAPH
  if (dump_cells) {
    if (myrank == 0) {
      error(
          "complete cell dumps are only created when "
          "configured with --enable-cell-graph.");
    }
  }
#endif

378
#ifndef SWIFT_DEBUG_THREADPOOL
379
  if (dump_threadpool) {
Peter W. Draper's avatar
Peter W. Draper committed
380
381
382
383
    printf(
        "Error: threadpool dumping is only possible if SWIFT was "
        "configured with the --enable-threadpool-debugging option.\n");
    return 1;
384
  }
385
#endif
386

387
388
389
  /* The CPU frequency is a long long, so we need to parse that ourselves. */
  if (cpufreqarg != NULL) {
    if (sscanf(cpufreqarg, "%llu", &cpufreq) != 1) {
Peter W. Draper's avatar
Peter W. Draper committed
390
391
      if (myrank == 0)
        printf("Error parsing CPU frequency (%s).\n", cpufreqarg);
392
      return 1;
393
    }
394
  }
395

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
396
  if (!with_self_gravity && !with_hydro && !with_external_gravity) {
Peter W. Draper's avatar
Peter W. Draper committed
397
398
399
400
    if (myrank == 0) {
      argparse_usage(&argparse);
      printf("\nError: At least one of -s, -g or -G must be chosen.\n");
    }
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
401
402
    return 1;
  }
403
  if (with_stars && !with_external_gravity && !with_self_gravity) {
Peter W. Draper's avatar
Peter W. Draper committed
404
405
    if (myrank == 0) {
      argparse_usage(&argparse);
406
      printf(
Peter W. Draper's avatar
Peter W. Draper committed
407
408
409
          "\nError: Cannot process stars without gravity, -g or -G "
          "must be chosen.\n");
    }
410
411
    return 1;
  }
412

413
  if (with_black_holes && !with_self_gravity) {
414
415
416
    if (myrank == 0) {
      argparse_usage(&argparse);
      printf(
417
418
          "\nError: Cannot process black holes without self-gravity, -G must "
          "be chosen.\n");
419
420
421
422
    }
    return 1;
  }

423
424
425
426
427
428
  if (with_fof) {
#ifndef WITH_FOF
    error("Running with FOF but compiled without it!");
#endif
  }

429
  if (with_fof && !with_self_gravity) {
430
431
432
433
434
435
    if (myrank == 0)
      printf(
          "Error: Cannot perform FOF search without gravity, -g or -G must be "
          "chosen.\n");
    return 1;
  }
436

437
438
439
440
441
442
443
444
  if (with_fof && !with_black_holes) {
    if (myrank == 0)
      printf(
          "Error: Cannot perform FOF seeding without black holes being in use, "
          "-B must be chosen.\n");
    return 1;
  }

445
446
447
448
449
450
451
452
453
454
  if (!with_stars && with_star_formation) {
    if (myrank == 0) {
      argparse_usage(&argparse);
      printf(
          "\nError: Cannot process star formation without stars, --stars must "
          "be chosen.\n");
    }
    return 1;
  }

Loic Hausammann's avatar
Loic Hausammann committed
455
  if (!with_stars && with_feedback) {
Peter W. Draper's avatar
Peter W. Draper committed
456
457
    if (myrank == 0) {
      argparse_usage(&argparse);
Loic Hausammann's avatar
Loic Hausammann committed
458
      printf(
459
460
461
462
463
464
465
466
467
468
469
          "\nError: Cannot process feedback without stars, --stars must be "
          "chosen.\n");
    }
    return 1;
  }

  if (!with_hydro && with_feedback) {
    if (myrank == 0) {
      argparse_usage(&argparse);
      printf(
          "\nError: Cannot process feedback without gas, --hydro must be "
Loic Hausammann's avatar
Loic Hausammann committed
470
          "chosen.\n");
Peter W. Draper's avatar
Peter W. Draper committed
471
    }
Loic Hausammann's avatar
Loic Hausammann committed
472
473
474
    return 1;
  }

475
476
477
478
479
480
481
482
483
484
  if (!with_hydro && with_black_holes) {
    if (myrank == 0) {
      argparse_usage(&argparse);
      printf(
          "\nError: Cannot process black holes without gas, --hydro must be "
          "chosen.\n");
    }
    return 1;
  }

485
486
/* Let's pin the main thread, now we know if affinity will be used. */
#if defined(HAVE_SETAFFINITY) && defined(HAVE_LIBNUMA) && defined(_GNU_SOURCE)
487
488
  if (with_aff &&
      ((ENGINE_POLICY)&engine_policy_setaffinity) == engine_policy_setaffinity)
489
490
491
    engine_pin();
#endif

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
492
493
  /* Genesis 1.1: And then, there was time ! */
  clocks_set_cpufreq(cpufreq);
494

495
496
497
  /* Are we running with gravity */
  const int with_gravity = (with_self_gravity || with_external_gravity);

498
499
500
  /* How vocal are we ? */
  const int talking = (verbose == 1 && myrank == 0) || (verbose == 2);

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
501
502
503
  if (myrank == 0 && dry_run)
    message(
        "Executing a dry run. No i/o or time integration will be performed.");
504

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
505
506
507
508
509
  /* Report CPU frequency.*/
  cpufreq = clocks_get_cpufreq();
  if (myrank == 0) {
    message("CPU frequency used for tick conversion: %llu Hz", cpufreq);
  }
510

Matthieu Schaller's avatar
Matthieu Schaller committed
511
/* Report host name(s). */
512
#ifdef WITH_MPI
513
  if (talking) {
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
514
515
    message("Rank %d running on: %s", myrank, hostname());
  }
516
#else
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
517
  message("Running on: %s", hostname());
518
519
#endif

520
521
/* Do we have debugging checks ? */
#ifdef SWIFT_DEBUG_CHECKS
522
523
  if (myrank == 0)
    message("WARNING: Debugging checks activated. Code will be slower !");
524
525
#endif

526
527
528
529
530
531
532
/* Do we have debugging checks ? */
#ifdef SWIFT_USE_NAIVE_INTERACTIONS
  if (myrank == 0)
    message(
        "WARNING: Naive cell interactions activated. Code will be slower !");
#endif

533
534
535
536
537
538
539
540
541
/* Do we have gravity accuracy checks ? */
#ifdef SWIFT_GRAVITY_FORCE_CHECKS
  if (myrank == 0)
    message(
        "WARNING: Checking 1/%d of all gpart for gravity accuracy. Code will "
        "be slower !",
        SWIFT_GRAVITY_FORCE_CHECKS);
#endif

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
542
543
  /* Do we choke on FP-exceptions ? */
  if (with_fp_exceptions) {
544
#ifdef HAVE_FE_ENABLE_EXCEPT
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
545
    feenableexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW);
546
#endif
547
548
    if (myrank == 0)
      message("WARNING: Floating point exceptions will be reported.");
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
549
  }
550

551
/* Do we have slow barriers? */
552
#ifndef HAVE_PTHREAD_BARRIERS
553
  if (myrank == 0)
554
    message("WARNING: Non-optimal thread barriers are being used.");
555
556
#endif

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
557
558
  /* How large are the parts? */
  if (myrank == 0) {
559
560
561
    message("sizeof(part)        is %4zi bytes.", sizeof(struct part));
    message("sizeof(xpart)       is %4zi bytes.", sizeof(struct xpart));
    message("sizeof(spart)       is %4zi bytes.", sizeof(struct spart));
562
    message("sizeof(bpart)       is %4zi bytes.", sizeof(struct bpart));
563
564
565
566
567
    message("sizeof(gpart)       is %4zi bytes.", sizeof(struct gpart));
    message("sizeof(multipole)   is %4zi bytes.", sizeof(struct multipole));
    message("sizeof(grav_tensor) is %4zi bytes.", sizeof(struct grav_tensor));
    message("sizeof(task)        is %4zi bytes.", sizeof(struct task));
    message("sizeof(cell)        is %4zi bytes.", sizeof(struct cell));
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
568
  }
569

Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
570
  /* Read the parameter file */
Matthieu Schaller's avatar
Matthieu Schaller committed
571
572
  struct swift_params *params =
      (struct swift_params *)malloc(sizeof(struct swift_params));
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
573
574
  if (params == NULL) error("Error allocating memory for the parameter file.");
  if (myrank == 0) {
575
576
    message("Reading runtime parameters from file '%s'", param_filename);
    parser_read_file(param_filename, params);
577
578

    /* Handle any command-line overrides. */
579
    if (cmdps.nparam > 0) {
580
581
582
      message(
          "Overwriting values read from the YAML file with command-line "
          "values.");
Peter W. Draper's avatar
Peter W. Draper committed
583
584
      for (int k = 0; k < cmdps.nparam; k++)
        parser_set_param(params, cmdps.param[k]);
585
    }
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
586
  }
587
#ifdef WITH_MPI
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
588
589
  /* Broadcast the parameter file */
  MPI_Bcast(params, sizeof(struct swift_params), MPI_BYTE, 0, MPI_COMM_WORLD);
590
#endif
591

592
593
  /* Temporary early aborts for modes not supported over MPI. */
#ifdef WITH_MPI
Matthieu Schaller's avatar
Matthieu Schaller committed
594
595
  if (with_mpole_reconstruction && nr_nodes > 1)
    error("Cannot reconstruct m-poles every step over MPI (yet).");
596
597
598
599
  if (with_timestep_limiter)
    error("Can't run with time-step limiter over MPI (yet)");
  if (with_timestep_sync)
    error("Can't run with time-step synchronization over MPI (yet)");
600
601
#endif

602
  /* Temporary early aborts for modes not supported with hand-vec. */
603
604
#if defined(WITH_VECTORIZATION) && defined(GADGET2_SPH) && \
    !defined(CHEMISTRY_NONE)
605
606
607
  error(
      "Cannot run with chemistry and hand-vectorization (yet). "
      "Use --disable-hand-vec at configure time.");
608
609
#endif

610
  /* Check that we can write the snapshots by testing if the output
611
   * directory exists and is searchable and writable. */
612
613
614
  char basename[PARSER_MAX_LINE_SIZE];
  parser_get_param_string(params, "Snapshots:basename", basename);
  const char *dirp = dirname(basename);
615
616
  if (access(dirp, W_OK | X_OK) != 0) {
    error("Cannot write snapshots in directory %s (%s)", dirp, strerror(errno));
617
618
  }

Matthieu Schaller's avatar
Matthieu Schaller committed
619
  /* Check that we can write the structure finding catalogues by testing if the
620
   * output directory exists and is searchable and writable. */
Matthieu Schaller's avatar
Matthieu Schaller committed
621
  if (with_structure_finding) {
622
    char stfbasename[PARSER_MAX_LINE_SIZE];
623
    parser_get_param_string(params, "StructureFinding:basename", stfbasename);
624
625
    const char *stfdirp = dirname(stfbasename);
    if (access(stfdirp, W_OK | X_OK) != 0) {
Matthieu Schaller's avatar
Matthieu Schaller committed
626
627
      error("Cannot write stf catalogues in directory %s (%s)", stfdirp,
            strerror(errno));
628
629
    }
  }
630

631
  /* Prepare the domain decomposition scheme */
632
  struct repartition reparttype;
633
#ifdef WITH_MPI
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
634
635
636
637
638
  struct partition initial_partition;
  partition_init(&initial_partition, &reparttype, params, nr_nodes);

  /* Let's report what we did */
  if (myrank == 0) {
639
640
#if defined(HAVE_PARMETIS)
    if (reparttype.usemetis)
Matthieu Schaller's avatar
Matthieu Schaller committed
641
      message("Using METIS serial partitioning:");
642
    else
Matthieu Schaller's avatar
Matthieu Schaller committed
643
      message("Using ParMETIS partitioning:");
644
#elif defined(HAVE_METIS)
645
    message("Using METIS serial partitioning:");
646
647
#else
    message("Non-METIS partitioning:");
648
649
#endif
    message("  initial partitioning: %s",
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
650
651
            initial_partition_name[initial_partition.type]);
    if (initial_partition.type == INITPART_GRID)
652
      message("    grid set to [ %i %i %i ].", initial_partition.grid[0],
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
653
              initial_partition.grid[1], initial_partition.grid[2]);
654
    message("  repartitioning: %s", repartition_name[reparttype.type]);
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
655
  }
656
#endif
657

658
  /* Common variables for restart and IC sections. */
659
  int clean_smoothing_length_values = 0;
660
661
662
  int flag_entropy_ICs = 0;

  /* Work out where we will read and write restart files. */
663
  char restart_dir[PARSER_MAX_LINE_SIZE];
Peter W. Draper's avatar
Peter W. Draper committed
664
665
  parser_get_opt_param_string(params, "Restarts:subdir", restart_dir,
                              "restart");
666
667

  /* The directory must exist. */
668
  if (myrank == 0) {
669
    if (access(restart_dir, W_OK | X_OK) != 0) {
670
      if (restart) {
671
        error("Cannot restart as no restart subdirectory: %s (%s)", restart_dir,
672
              strerror(errno));
673
      } else {
Peter W. Draper's avatar
Peter W. Draper committed
674
675
676
        if (mkdir(restart_dir, 0777) != 0)
          error("Failed to create restart directory: %s (%s)", restart_dir,
                strerror(errno));
677
      }
678
    }
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
679
680
  }

681
  /* Basename for any restart files. */
682
683
  char restart_name[PARSER_MAX_LINE_SIZE];
  parser_get_opt_param_string(params, "Restarts:basename", restart_name,
684
                              "swift");
685

686
687
  /* How often to check for the stop file and dump restarts and exit the
   * application. */
688
  const int restart_stop_steps =
Peter W. Draper's avatar
Peter W. Draper committed
689
      parser_get_opt_param_int(params, "Restarts:stop_steps", 100);
690

691
692
693
694
695
696
697
698
699
700
701
  /* Get the maximal wall-clock time of this run */
  const float restart_max_hours_runtime =
      parser_get_opt_param_float(params, "Restarts:max_run_time", FLT_MAX);

  /* Do we want to resubmit when we hit the limit? */
  const int resubmit_after_max_hours =
      parser_get_opt_param_int(params, "Restarts:resubmit_on_exit", 0);

  /* What command should we run to resubmit at the end? */
  char resubmit_command[PARSER_MAX_LINE_SIZE];
  if (resubmit_after_max_hours)
702
703
    parser_get_param_string(params, "Restarts:resubmit_command",
                            resubmit_command);
704

705
  /* If restarting, look for the restart files. */
706
  if (restart) {
707

708
    /* Attempting a restart. */
709
710
    char **restart_files = NULL;
    int restart_nfiles = 0;
711
712

    if (myrank == 0) {
713
      message("Restarting SWIFT");
714

715
      /* Locate the restart files. */
Peter W. Draper's avatar
Peter W. Draper committed
716
717
      restart_files =
          restart_locate(restart_dir, restart_name, &restart_nfiles);
718
719
      if (restart_nfiles == 0)
        error("Failed to locate any restart files in %s", restart_dir);
720
721

      /* We need one file per rank. */
722
      if (restart_nfiles != nr_nodes)
723
        error("Incorrect number of restart files, expected %d found %d",
724
              nr_nodes, restart_nfiles);
725
726

      if (verbose > 0)
727
728
        for (int i = 0; i < restart_nfiles; i++)
          message("found restart file: %s", restart_files[i]);
729
730
731
    }

#ifdef WITH_MPI
732
    /* Distribute the restart files, need one for each rank. */
733
734
735
    if (myrank == 0) {

      for (int i = 1; i < nr_nodes; i++) {
736
737
        strcpy(restart_file, restart_files[i]);
        MPI_Send(restart_file, 200, MPI_BYTE, i, 0, MPI_COMM_WORLD);
738
739
740
      }

      /* Keep local file. */
741
      strcpy(restart_file, restart_files[0]);
742
743

      /* Finished with the list. */
744
      restart_locate_free(restart_nfiles, restart_files);
745
746

    } else {
747
      MPI_Recv(restart_file, 200, MPI_BYTE, 0, 0, MPI_COMM_WORLD,
748
749
               MPI_STATUS_IGNORE);
    }
750
    if (verbose > 1) message("local restart file = %s", restart_file);
751
#else
752

753
    /* Just one restart file. */
754
    strcpy(restart_file, restart_files[0]);
755
756
#endif

757
    /* Now read it. */
758
    restart_read(&e, restart_file);
759
760
761

    /* And initialize the engine with the space and policies. */
    if (myrank == 0) clocks_gettime(&tic);
762
763
    engine_config(/*restart=*/1, /*fof=*/0, &e, params, nr_nodes, myrank,
                  nr_threads, with_aff, talking, restart_file);
764
765
    if (myrank == 0) {
      clocks_gettime(&toc);
766
      message("engine_config took %.3f %s.", clocks_diff(&tic, &toc),
767
768
769
770
              clocks_getunit());
      fflush(stdout);
    }

771
772
    /* Check if we are already done when given steps on the command-line. */
    if (e.step >= nsteps && nsteps > 0)
773
      error("Not restarting, already completed %d steps", e.step);
774

775
  } else {
776

777
    /* Not restarting so look for the ICs. */
778
    /* Initialize unit system and constants */
779
    units_init_from_params(&us, params, "InternalUnitSystem");
780
    phys_const_init(&us, params, &prog_const);
781
    if (myrank == 0) {
782
783
784
785
786
787
788
      message("Internal unit system: U_M = %e g.", us.UnitMass_in_cgs);
      message("Internal unit system: U_L = %e cm.", us.UnitLength_in_cgs);
      message("Internal unit system: U_t = %e s.", us.UnitTime_in_cgs);
      message("Internal unit system: U_I = %e A.", us.UnitCurrent_in_cgs);
      message("Internal unit system: U_T = %e K.", us.UnitTemperature_in_cgs);
      phys_const_print(&prog_const);
    }
Pedro Gonnet's avatar
oops.    
Pedro Gonnet committed
789

790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
    /* Read particles and space information from ICs */
    char ICfileName[200] = "";
    parser_get_param_string(params, "InitialConditions:file_name", ICfileName);
    const int periodic =
        parser_get_param_int(params, "InitialConditions:periodic");
    const int replicate =
        parser_get_opt_param_int(params, "InitialConditions:replicate", 1);
    clean_smoothing_length_values = parser_get_opt_param_int(
        params, "InitialConditions:cleanup_smoothing_lengths", 0);
    const int cleanup_h = parser_get_opt_param_int(
        params, "InitialConditions:cleanup_h_factors", 0);
    const int cleanup_sqrt_a = parser_get_opt_param_int(
        params, "InitialConditions:cleanup_velocity_factors", 0);
    const int generate_gas_in_ics = parser_get_opt_param_int(
        params, "InitialConditions:generate_gas_in_ics", 0);

806
    /* Initialise the cosmology */
807
808
809
    if (with_cosmology)
      cosmology_init(params, &us, &prog_const, &cosmo);
    else
810
      cosmology_init_no_cosmo(&cosmo);
811
    if (myrank == 0 && with_cosmology) cosmology_print(&cosmo);
812

813
    /* Initialise the hydro properties */
814
815
    if (with_hydro)
      hydro_props_init(&hydro_properties, &prog_const, &us, params);
816
817
818
819
820
821
822
823
    else
      bzero(&hydro_properties, sizeof(struct hydro_props));

    /* Initialise the equation of state */
    if (with_hydro)
      eos_init(&eos, &prog_const, &us, params);
    else
      bzero(&eos, sizeof(struct eos_parameters));
824

825
826
827
828
829
830
831
    /* Initialise the entropy floor */
    if (with_hydro)
      entropy_floor_init(&entropy_floor, &prog_const, &us, &hydro_properties,
                         params);
    else
      bzero(&entropy_floor, sizeof(struct entropy_floor_properties));

Loic Hausammann's avatar
Loic Hausammann committed
832
833
834
835
836
837
838
    /* Initialise the pressure floor */
    if (with_hydro)
      pressure_floor_init(&pressure_floor_props, &prog_const, &us,
                          &hydro_properties, params);
    else
      bzero(&pressure_floor_props, sizeof(struct pressure_floor_properties));

Loic Hausammann's avatar
Loic Hausammann committed
839
840
    /* Initialise the stars properties */
    if (with_stars)
841
      stars_props_init(&stars_properties, &prog_const, &us, params,
842
                       &hydro_properties, &cosmo);
Loic Hausammann's avatar
Loic Hausammann committed
843
844
845
    else
      bzero(&stars_properties, sizeof(struct stars_props));

846
847
848
849
850
851
    /* Initialise the feedback properties */
    if (with_feedback) {
#ifdef FEEDBACK_NONE
      error("ERROR: Running with feedback but compiled without it.");
#endif
      feedback_props_init(&feedback_properties, &prog_const, &us, params,
852
853
                          &hydro_properties, &cosmo);
    } else
854
      bzero(&feedback_properties, sizeof(struct feedback_props));
855

856
857
858
859
860
861
862
863
864
865
    /* Initialise the black holes properties */
    if (with_black_holes) {
#ifdef BLACK_HOLES_NONE
      error("ERROR: Running with black_holes but compiled without it.");
#endif
      black_holes_props_init(&black_holes_properties, &prog_const, &us, params,
                             &hydro_properties, &cosmo);
    } else
      bzero(&black_holes_properties, sizeof(struct black_holes_props));

Peter W. Draper's avatar
Peter W. Draper committed
866
      /* Initialise the cooling function properties */
867
#ifdef COOLING_NONE
868
    if (with_cooling) {
Peter W. Draper's avatar
Peter W. Draper committed
869
      error(
870
          "ERROR: Running with cooling calculation"
Peter W. Draper's avatar
Peter W. Draper committed
871
          " but compiled without it.");
872
873
874
    }
#else
    if (!with_cooling && !with_temperature) {
Peter W. Draper's avatar
Peter W. Draper committed
875
876
877
      error(
          "ERROR: Compiled with cooling but running without it. "
          "Did you forget the --cooling or --temperature flags?");
878
    }
879
#endif
880
881
    bzero(&cooling_func, sizeof(struct cooling_function_data));
    if (with_cooling || with_temperature) {
882
      cooling_init(params, &us, &prog_const, &hydro_properties, &cooling_func);
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
    }
    if (myrank == 0) cooling_print(&cooling_func);

    /* Initialise the star formation law and its properties */
    bzero(&starform, sizeof(struct star_formation));
    if (with_star_formation) {
#ifdef STAR_FORMATION_NONE
      error("ERROR: Running with star formation but compiled without it!");
#endif
      starformation_init(params, &prog_const, &us, &hydro_properties,
                         &starform);
    }
    if (with_star_formation && myrank == 0) starformation_print(&starform);

    /* Initialise the chemistry */
    bzero(&chemistry, sizeof(struct chemistry_global_data));
    chemistry_init(params, &us, &prog_const, &chemistry);
    if (myrank == 0) chemistry_print(&chemistry);