/*******************************************************************************
* This file is part of SWIFT.
* Copyright (c) 2019 Stuart McAlpine (stuart.mcalpine@helsinki.fi)
* 2024 Will J. Roper (w.roper@sussex.ac.uk)
*
* 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.
*
* 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.
*
* You should have received a copy of the GNU Lesser General Public License
* along with this program. If not, see .
*
******************************************************************************/
/* Config header */
#include
/* Standard includes */
#include
/* Local includes */
#include "cell.h"
#include "engine.h"
#include "multipole.h"
#include "space.h"
#include "zoom.h"
/**
* @brief Find which background or buffer cells contain the zoom region.
*
* A void cell is a low resolution cell above the zoom region (or part
* of it).
*
* A void cell is always in the cell grid directly above the zoom cells, i.e. if
* there are buffer cells, the void cells are the buffer cells, if there are no
* buffer cells, the void cells are the background cells.
*
* @param s The space.
* @param verbose Are we talking?
*/
void zoom_find_void_cells(struct space *s, const int verbose) {
/* Get the zoom properties */
struct zoom_region_properties *zoom_props = s->zoom_props;
/* Get the cell pointers. */
struct cell *cells = s->cells_top;
/* Zero the void cell count. */
zoom_props->nr_void_cells = 0;
/* Get the offset and the number of cells we're dealing with. */
int offset = zoom_props->bkg_cell_offset;
int ncells = zoom_props->nr_bkg_cells;
/* Work out how many void cells we should have. */
int void_cdim = s->zoom_props->void_dim[0] * s->iwidth[0] * 1.0001;
int target_void_count = void_cdim * void_cdim * void_cdim;
/* Allocate the indices of void cells */
if (swift_memalign(
"void_cell_indices", (void **)&s->zoom_props->void_cell_indices,
SWIFT_STRUCT_ALIGNMENT, target_void_count * sizeof(int)) != 0)
error("Failed to allocate indices of local top-level background cells.");
bzero(s->zoom_props->void_cell_indices, target_void_count * sizeof(int));
/* Loop over the background cells and find cells containing
* the void region. */
for (int cid = offset; cid < offset + ncells; cid++) {
/* Get the cell */
struct cell *c = &cells[cid];
/* Label this cell if it contains the zoom region. */
if (zoom_cell_overlaps_zoom_region(c, s)) {
c->subtype = cell_subtype_void;
zoom_props->void_cell_indices[zoom_props->nr_void_cells++] = cid;
}
}
if (target_void_count != zoom_props->nr_void_cells)
error(
"Not all void cells were found and labelled! (target_void_count=%d, "
"found_void_count=%d)",
target_void_count, zoom_props->nr_void_cells);
if (verbose) message("Found %d void cells", zoom_props->nr_void_cells);
#ifdef SWIFT_DEBUG_CHECKS
/* Check the void cells are in the right place. */
for (int i = 0; i < zoom_props->nr_void_cells; i++) {
const int cid = zoom_props->void_cell_indices[i];
if (cid < offset || cid >= offset + ncells)
error("Void cell index is out of range (cid=%d, offset=%d, ncells=%d)",
cid, offset, ncells);
if (zoom_cell_overlaps_zoom_region(&cells[cid], s) == 0)
error("Void cell is not inside the zoom region (cid=%d)", cid);
}
#endif
/* We also need to label the buffer cells as void cells if they
* are within the zoom region. */
int nr_buffer_void = 0;
if (zoom_props->with_buffer_cells) {
offset = zoom_props->buffer_cell_offset;
ncells = zoom_props->nr_buffer_cells;
/* Loop over the buffer cells and find cells containing
* the zoom region. */
for (int cid = offset; cid < offset + ncells; cid++) {
/* Get the cell */
struct cell *c = &cells[cid];
/* Label this cell if it contains the zoom region. */
if (zoom_cell_overlaps_zoom_region(c, s)) {
c->subtype = cell_subtype_void;
nr_buffer_void++;
}
}
if (verbose)
message("Found %d void cells in the buffer region", nr_buffer_void);
}
}
/**
* @brief Find what background or buffer cells surround the zoom region.
*
* When interacting background cells and zoom TL cells, it helps to know
* which background cells are within the mesh distance of the zoom region.
* These cells are the cells which will interact via pair/grav or long_range
* gravity tasks, and get tagged as "neighbour" cells.
*
* @param s The space.
* @param verbose Are we talking?
*/
void zoom_find_neighbouring_cells(struct space *s, const int verbose) {
/* Get the zoom properties */
struct zoom_region_properties *zoom_props = s->zoom_props;
/* Get the right cell cdim. */
int cdim[3];
double iwidth[3];
if (zoom_props->with_buffer_cells) {
cdim[0] = zoom_props->buffer_cdim[0];
cdim[1] = zoom_props->buffer_cdim[1];
cdim[2] = zoom_props->buffer_cdim[2];
iwidth[0] = zoom_props->buffer_iwidth[0];
iwidth[1] = zoom_props->buffer_iwidth[1];
iwidth[2] = zoom_props->buffer_iwidth[2];
} else {
cdim[0] = s->cdim[0];
cdim[1] = s->cdim[1];
cdim[2] = s->cdim[2];
iwidth[0] = s->iwidth[0];
iwidth[1] = s->iwidth[1];
iwidth[2] = s->iwidth[2];
}
/* Get the cell pointers. */
struct cell *cells = s->cells_top;
/* Get the right offset and the number of cells we're dealing with. */
int offset;
int ncells;
if (zoom_props->with_buffer_cells) {
offset = zoom_props->buffer_cell_offset;
ncells = zoom_props->nr_buffer_cells;
} else {
offset = zoom_props->bkg_cell_offset;
ncells = zoom_props->nr_bkg_cells;
}
/* At this point we can only define neighbour cells by cell properties,
* leaving the fancy gravity distance criterion for task creation later.
* Here we just make sure all possible neighbour cells are flagged
* as such. */
/* Maximal distance any interaction can take place
* before the mesh kicks in, rounded up to the next integer */
const int delta_cells = ceil(zoom_props->neighbour_distance *
max3(iwidth[0], iwidth[1], iwidth[2])) +
1;
/* Turn this into upper and lower bounds for loops */
int delta_m = delta_cells;
int delta_p = delta_cells;
/* Special case where every cell is in range of every other one */
if (delta_cells > cdim[0]) {
delta_m = cdim[0];
delta_p = cdim[0];
}
/* Let's be verbose about this choice */
if (verbose)
message(
"Looking for neighbouring background cells up to %d background "
"top-level "
"cells away from the zoom region (delta_m=%d "
"delta_p=%d)",
delta_cells, delta_m, delta_p);
/* Allocate the indices of neighbour background cells. */
/* We don't know how many we will need at this point so have
* to allocate assuming we will have all cells in range. */
if (swift_memalign("neighbour_cells_top",
(void **)&s->zoom_props->neighbour_cells_top,
SWIFT_STRUCT_ALIGNMENT, ncells * sizeof(int)) != 0)
error("Failed to allocate indices of local top-level background cells.");
bzero(s->zoom_props->neighbour_cells_top, ncells * sizeof(int));
/* Get a pointer to the neighbour cells index array. */
int *neighbour_cells_top = s->zoom_props->neighbour_cells_top;
/* Loop over background/buffer cells. We'll talk out delta_cells cells from
* any void cells. */
for (int i = 0; i < cdim[0]; i++) {
for (int j = 0; j < cdim[1]; j++) {
for (int k = 0; k < cdim[2]; k++) {
/* Get this cell. */
const size_t cid = cell_getid_offset(cdim, offset, i, j, k);
struct cell *ci = &cells[cid];
/* Skip non-void cells, we only want to walk out. */
if (ci->subtype != cell_subtype_void) continue;
/* Loop over every other cell within (Manhattan) range delta */
for (int ii = i - delta_m; ii <= i + delta_p; ii++) {
/* Escape beyond range */
if (ii < 0 || ii >= cdim[0]) continue;
for (int jj = j - delta_m; jj <= j + delta_p; jj++) {
/* Escape beyond range */
if (jj < 0 || jj >= cdim[1]) continue;
for (int kk = k - delta_m; kk <= k + delta_p; kk++) {
/* Escape beyond range */
if (kk < 0 || kk >= cdim[2]) continue;
/* Get this cell. */
const int cjd = cell_getid_offset(cdim, offset, ii, jj, kk);
struct cell *cj = &cells[cjd];
/* Ensure this neighbour isn't a void cell or an already
* counted neighbour. */
if (cj->subtype != cell_subtype_void &&
cj->subtype != cell_subtype_neighbour) {
/* Record that we've found a neighbour. */
cj->subtype = cell_subtype_neighbour;
neighbour_cells_top[zoom_props->nr_neighbour_cells++] = cjd;
}
} /* neighbour k loop */
} /* neighbour j loop */
} /* neighbour i loop */
} /* k loop */
} /* j loop */
} /* i loop */
if (verbose)
message("%i cells neighbour the zoom region",
zoom_props->nr_neighbour_cells);
}
/**
* @brief Run through all cells and ensure they have the correct cell type and
* width for their position in s->cells_top.
*
* @param s The space.
*/
void zoom_verify_cell_type(struct space *s) {
#ifdef SWIFT_DEBUG_CHECKS
/* Get the cells array and cell properties */
struct cell *cells = s->cells_top;
const int bkg_cell_offset = s->zoom_props->bkg_cell_offset;
const int buffer_offset = s->zoom_props->buffer_cell_offset;
const double *zoom_width = s->zoom_props->width;
const double *width = s->width;
/* Loop over all cells */
for (int cid = 0; cid < s->nr_cells; cid++) {
/* Check cell type */
if (cid < bkg_cell_offset && cells[cid].type != cell_type_zoom)
error(
"Cell has the wrong cell type for it's array position (cid=%d, "
"c->type=%s, "
"s->zoom_props->bkg_cell_offset=%d)",
cid, cellID_names[cells[cid].type], bkg_cell_offset);
if (cid >= bkg_cell_offset && cells[cid].type == cell_type_zoom)
error(
"Cell has the wrong cell type for it's array position (cid=%d, "
"c->type=%s, "
"s->zoom_props->bkg_cell_offset=%d)",
cid, cellID_names[cells[cid].type], bkg_cell_offset);
/* Check cell widths */
for (int i = 0; i < 3; i++) {
if (cid < bkg_cell_offset && cells[cid].width[i] != zoom_width[i])
error(
"Cell has the wrong cell width for it's array position (cid=%d, "
"c->type=%s, "
"s->zoom_props->bkg_cell_offset=%d, c->width=[%f %f %f], "
"s->zoom_props->width=[%f %f %f])",
cid, cellID_names[cells[cid].type], bkg_cell_offset,
cells[cid].width[0], cells[cid].width[1], cells[cid].width[2],
s->zoom_props->width[0], s->zoom_props->width[1],
s->zoom_props->width[2]);
if ((cid >= bkg_cell_offset && cid < buffer_offset) &&
cells[cid].width[i] != width[i])
error(
"Cell has the wrong cell width for it's array position (cid=%d, "
"c->type=%s, "
"s->zoom_props->bkg_cell_offset=%d, c->width=[%f %f %f], "
"s->zoom_props->width=[%f %f %f])",
cid, cellID_names[cells[cid].type], bkg_cell_offset,
cells[cid].width[0], cells[cid].width[1], cells[cid].width[2],
s->zoom_props->width[0], s->zoom_props->width[1],
s->zoom_props->width[2]);
}
}
/* Ensure the cell and region boundaries align. */
if (s->zoom_props->with_buffer_cells) {
int found_bkg_bufferi_low = 0;
int found_bkg_bufferj_low = 0;
int found_bkg_bufferk_low = 0;
int found_bkg_bufferi_up = 0;
int found_bkg_bufferj_up = 0;
int found_bkg_bufferk_up = 0;
double tol = 0.001 * s->width[0]; /* Tolerence on edge matching */
for (int i = 0; i < s->cdim[0]; i++) {
for (int j = 0; j < s->cdim[1]; j++) {
for (int k = 0; k < s->cdim[2]; k++) {
/* Define the cell position. */
const double pos[3] = {s->width[0] * i, s->width[1] * j,
s->width[2] * k};
/* Test the lower boundary. */
if (fabs(pos[0] - s->zoom_props->buffer_lower_bounds[0]) < tol)
found_bkg_bufferi_low = 1;
if (fabs(pos[1] - s->zoom_props->buffer_lower_bounds[1]) < tol)
found_bkg_bufferj_low = 1;
if (fabs(pos[2] - s->zoom_props->buffer_lower_bounds[2]) < tol)
found_bkg_bufferk_low = 1;
/* Test the upper boundary. */
if (fabs(pos[0] - s->zoom_props->buffer_upper_bounds[0]) < tol)
found_bkg_bufferi_up = 1;
if (fabs(pos[1] - s->zoom_props->buffer_upper_bounds[1]) < tol)
found_bkg_bufferj_up = 1;
if (fabs(pos[2] - s->zoom_props->buffer_upper_bounds[2]) < tol)
found_bkg_bufferk_up = 1;
}
}
}
/* Did we find the boundaries?. */
if (!found_bkg_bufferi_low || !found_bkg_bufferj_low ||
!found_bkg_bufferk_low || !found_bkg_bufferi_up ||
!found_bkg_bufferj_up || !found_bkg_bufferk_up)
error("The background cell and buffer region edges don't match!");
/* And for the zoom cells. */
int found_buffer_zoomi_low = 0;
int found_buffer_zoomj_low = 0;
int found_buffer_zoomk_low = 0;
int found_buffer_zoomi_up = 0;
int found_buffer_zoomj_up = 0;
int found_buffer_zoomk_up = 0;
tol = 0.001 * s->zoom_props->buffer_width[0];
for (int i = 0; i < s->zoom_props->buffer_cdim[0]; i++) {
for (int j = 0; j < s->zoom_props->buffer_cdim[1]; j++) {
for (int k = 0; k < s->zoom_props->buffer_cdim[2]; k++) {
/* Define the cell position. */
const double pos[3] = {s->zoom_props->buffer_lower_bounds[0] +
s->zoom_props->buffer_width[0] * i,
s->zoom_props->buffer_lower_bounds[1] +
s->zoom_props->buffer_width[1] * j,
s->zoom_props->buffer_lower_bounds[2] +
s->zoom_props->buffer_width[2] * k};
/* Test the lower boundary. */
if (fabs(pos[0] - s->zoom_props->region_lower_bounds[0]) < tol)
found_buffer_zoomi_low = 1;
if (fabs(pos[1] - s->zoom_props->region_lower_bounds[1]) < tol)
found_buffer_zoomj_low = 1;
if (fabs(pos[2] - s->zoom_props->region_lower_bounds[2]) < tol)
found_buffer_zoomk_low = 1;
/* Test the upper boundary. */
if (fabs(pos[0] - s->zoom_props->region_upper_bounds[0]) < tol)
found_buffer_zoomi_up = 1;
if (fabs(pos[1] - s->zoom_props->region_upper_bounds[1]) < tol)
found_buffer_zoomj_up = 1;
if (fabs(pos[2] - s->zoom_props->region_upper_bounds[2]) < tol)
found_buffer_zoomk_up = 1;
}
}
}
/* Did we find the boundaries?. */
if (!found_buffer_zoomi_low || !found_buffer_zoomj_low ||
!found_buffer_zoomk_low || !found_buffer_zoomi_up ||
!found_buffer_zoomj_up || !found_buffer_zoomk_up)
error("The buffer cell and zoom region edges don't match!");
} else {
/* Check the background and zoom cells align. */
int found_bkg_zoomi_low = 0;
int found_bkg_zoomj_low = 0;
int found_bkg_zoomk_low = 0;
int found_bkg_zoomi_up = 0;
int found_bkg_zoomj_up = 0;
int found_bkg_zoomk_up = 0;
const double tol = 0.001 * s->width[0];
for (int i = 0; i < s->cdim[0]; i++) {
for (int j = 0; j < s->cdim[1]; j++) {
for (int k = 0; k < s->cdim[2]; k++) {
/* Define the cell position. */
const double pos[3] = {s->width[0] * i, s->width[1] * j,
s->width[2] * k};
/* Test the lower boundary. */
if (fabs(pos[0] - s->zoom_props->region_lower_bounds[0]) < tol)
found_bkg_zoomi_low = 1;
if (fabs(pos[1] - s->zoom_props->region_lower_bounds[1]) < tol)
found_bkg_zoomj_low = 1;
if (fabs(pos[2] - s->zoom_props->region_lower_bounds[2]) < tol)
found_bkg_zoomk_low = 1;
/* Test the upper boundary. */
if (fabs(pos[0] - s->zoom_props->region_upper_bounds[0]) < tol)
found_bkg_zoomi_up = 1;
if (fabs(pos[1] - s->zoom_props->region_upper_bounds[1]) < tol)
found_bkg_zoomj_up = 1;
if (fabs(pos[2] - s->zoom_props->region_upper_bounds[2]) < tol)
found_bkg_zoomk_up = 1;
}
}
}
/* Did we find the boundaries?. */
if (!found_bkg_zoomi_low || !found_bkg_zoomj_low || !found_bkg_zoomk_low ||
!found_bkg_zoomi_up || !found_bkg_zoomj_up || !found_bkg_zoomk_up)
error("The background cell and zoom region edges don't match!");
}
#else
error("zoom_verify_cell_type() called without SWIFT_DEBUG_CHECKS defined");
#endif
}
/**
* @brief Build the TL cells, with a zoom region.
*
* This replaces the loop in space_regrid when running with a zoom region. It
* constructs all zoom, background and buffer cells (if required) and sets their
* initial values.
*
* Zoom cells occupy the first cells in the cells array, followed by background
* cells and then buffer cells (if required).
*
* @param s The space.
* @param ti_current The current time.
* @param verbose Are we talking?
*/
void zoom_construct_tl_cells(struct space *s, const integertime_t ti_current,
int verbose) {
const ticks tic = getticks();
/* Get the zoom properties */
struct zoom_region_properties *zoom_props = s->zoom_props;
/* Get some zoom region properties */
const float dmin_zoom =
min3(zoom_props->width[0], zoom_props->width[1], zoom_props->width[2]);
const int bkg_cell_offset = zoom_props->bkg_cell_offset;
const double zoom_bounds[3] = {zoom_props->region_lower_bounds[0],
zoom_props->region_lower_bounds[1],
zoom_props->region_lower_bounds[2]};
/* Loop over zoom cells and set locations and initial values */
for (int i = 0; i < zoom_props->cdim[0]; i++) {
for (int j = 0; j < zoom_props->cdim[1]; j++) {
for (int k = 0; k < zoom_props->cdim[2]; k++) {
const size_t cid = cell_getid(zoom_props->cdim, i, j, k);
/* Create the zoom cell and it's multipoles */
struct cell *c = &s->cells_top[cid];
c->loc[0] = (i * zoom_props->width[0]) + zoom_bounds[0];
c->loc[1] = (j * zoom_props->width[1]) + zoom_bounds[1];
c->loc[2] = (k * zoom_props->width[2]) + zoom_bounds[2];
c->width[0] = zoom_props->width[0];
c->width[1] = zoom_props->width[1];
c->width[2] = zoom_props->width[2];
if (s->with_self_gravity) c->grav.multipole = &s->multipoles_top[cid];
c->type = cell_type_zoom;
c->subtype = cell_subtype_regular;
c->dmin = dmin_zoom;
c->h_min_allowed =
c->dmin * 0.5 * (1. / kernel_gamma); /* Only zoom cells do hydro */
c->h_max_allowed =
c->dmin * (1. / kernel_gamma); /* Only zoom cells do hydro */
c->depth = 0;
c->split = 0;
c->hydro.count = 0;
c->grav.count = 0;
c->stars.count = 0;
c->sinks.count = 0;
c->top = c;
c->super = c;
c->hydro.super = c;
c->grav.super = c;
c->hydro.ti_old_part = ti_current;
c->grav.ti_old_part = ti_current;
c->stars.ti_old_part = ti_current;
c->sinks.ti_old_part = ti_current;
c->black_holes.ti_old_part = ti_current;
c->grav.ti_old_multipole = ti_current;
#ifdef WITH_MPI
c->mpi.tag = -1;
c->mpi.recv = NULL;
c->mpi.send = NULL;
#endif
#if defined(SWIFT_DEBUG_CHECKS) || defined(SWIFT_CELL_GRAPH)
cell_assign_top_level_cell_index(c, s);
#endif
}
}
}
/* Attach zoom cells to their cell array (they are the first cells in
* cells_top so nothing fancy needed here). */
zoom_props->zoom_cells_top = s->cells_top;
if (verbose)
message("Set zoom cell dimensions to [ %i %i %i ].", zoom_props->cdim[0],
zoom_props->cdim[1], zoom_props->cdim[2]);
/* Get the minimum cell size. */
const double dmin = min3(s->width[0], s->width[1], s->width[2]);
/* Loop over background cells and set locations and initial values */
for (int i = 0; i < s->cdim[0]; i++) {
for (int j = 0; j < s->cdim[1]; j++) {
for (int k = 0; k < s->cdim[2]; k++) {
const size_t cid = cell_getid_offset(s->cdim, bkg_cell_offset, i, j, k);
/* Background top level cells. */
struct cell *c = &s->cells_top[cid];
c->loc[0] = i * s->width[0];
c->loc[1] = j * s->width[1];
c->loc[2] = k * s->width[2];
c->width[0] = s->width[0];
c->width[1] = s->width[1];
c->width[2] = s->width[2];
c->dmin = dmin;
if (s->with_self_gravity) c->grav.multipole = &s->multipoles_top[cid];
c->type = cell_type_bkg;
c->subtype = cell_subtype_regular;
c->depth = 0;
c->split = 0;
c->hydro.count = 0;
c->grav.count = 0;
c->stars.count = 0;
c->sinks.count = 0;
c->top = c;
c->super = c;
c->hydro.super = c;
c->grav.super = c;
c->hydro.ti_old_part = ti_current;
c->grav.ti_old_part = ti_current;
c->stars.ti_old_part = ti_current;
c->sinks.ti_old_part = ti_current;
c->black_holes.ti_old_part = ti_current;
c->grav.ti_old_multipole = ti_current;
#ifdef WITH_MPI
c->mpi.tag = -1;
c->mpi.recv = NULL;
c->mpi.send = NULL;
#endif
#if defined(SWIFT_DEBUG_CHECKS) || defined(SWIFT_CELL_GRAPH)
cell_assign_top_level_cell_index(c, s);
#endif
}
}
}
/* Attach background cells to their cell array. */
zoom_props->bkg_cells_top = &s->cells_top[bkg_cell_offset];
if (verbose)
message("Set background cell dimensions to [ %i %i %i ].", s->cdim[0],
s->cdim[1], s->cdim[2]);
/* If we have a buffer region create buffer cells. */
if (zoom_props->with_buffer_cells) {
/* Get relevant information. */
const float dmin_buffer =
min3(zoom_props->buffer_width[0], zoom_props->buffer_width[1],
zoom_props->buffer_width[2]);
const int buffer_offset = zoom_props->buffer_cell_offset;
const double buffer_bounds[3] = {zoom_props->buffer_lower_bounds[0],
zoom_props->buffer_lower_bounds[1],
zoom_props->buffer_lower_bounds[2]};
/* Loop over buffer cells and set locations and initial values */
for (int i = 0; i < zoom_props->buffer_cdim[0]; i++) {
for (int j = 0; j < zoom_props->buffer_cdim[1]; j++) {
for (int k = 0; k < zoom_props->buffer_cdim[2]; k++) {
const size_t cid = cell_getid_offset(zoom_props->buffer_cdim,
buffer_offset, i, j, k);
/* Buffer top level cells. */
struct cell *c = &s->cells_top[cid];
c->loc[0] = (i * zoom_props->buffer_width[0]) + buffer_bounds[0];
c->loc[1] = (j * zoom_props->buffer_width[1]) + buffer_bounds[1];
c->loc[2] = (k * zoom_props->buffer_width[2]) + buffer_bounds[2];
c->width[0] = zoom_props->buffer_width[0];
c->width[1] = zoom_props->buffer_width[1];
c->width[2] = zoom_props->buffer_width[2];
c->dmin = dmin_buffer;
if (s->with_self_gravity) c->grav.multipole = &s->multipoles_top[cid];
c->type = cell_type_buffer;
c->subtype = cell_subtype_regular;
c->depth = 0;
c->split = 0;
c->hydro.count = 0;
c->grav.count = 0;
c->stars.count = 0;
c->sinks.count = 0;
c->top = c;
c->super = c;
c->hydro.super = c;
c->grav.super = c;
c->hydro.ti_old_part = ti_current;
c->grav.ti_old_part = ti_current;
c->stars.ti_old_part = ti_current;
c->sinks.ti_old_part = ti_current;
c->black_holes.ti_old_part = ti_current;
c->grav.ti_old_multipole = ti_current;
#ifdef WITH_MPI
c->mpi.tag = -1;
c->mpi.recv = NULL;
c->mpi.send = NULL;
#endif
#if defined(SWIFT_DEBUG_CHECKS) || defined(SWIFT_CELL_GRAPH)
cell_assign_top_level_cell_index(c, s);
#endif
}
}
}
/* Attach background cells to their cell array. */
zoom_props->buffer_cells_top = &s->cells_top[buffer_offset];
if (verbose)
message("Set buffer cell dimensions to [ %i %i %i ].",
zoom_props->buffer_cdim[0], zoom_props->buffer_cdim[1],
zoom_props->buffer_cdim[2]);
}
/* Now find and label what cells contain the zoom region. */
zoom_find_void_cells(s, verbose);
/* Find neighbours of the zoom region (cells within the distance past which
* the gravity mesh takes over). */
zoom_find_neighbouring_cells(s, verbose);
#ifdef SWIFT_DEBUG_CHECKS
/* Lets check all the cells are in the right place with the correct widths */
zoom_verify_cell_type(s);
#endif
if (verbose)
message("took %.3f %s.", clocks_from_ticks(getticks() - tic),
clocks_getunit());
}
/**
* @brief Collect the timestep from the zoom region to the void level.
*
* This will recurse to the void leaves and grab the timestep from the zoom
* cells, populating the void cell tree as it goes.
*
* @param c The #cell.
*/
static void zoom_void_timestep_collect_recursive(struct cell *c) {
/* Define the timestep info we'll be collecting. */
integertime_t ti_hydro_end_min = max_nr_timesteps, ti_hydro_beg_max = 0;
integertime_t ti_rt_end_min = max_nr_timesteps, ti_rt_beg_max = 0;
integertime_t ti_grav_end_min = max_nr_timesteps, ti_grav_beg_max = 0;
integertime_t ti_stars_end_min = max_nr_timesteps, ti_stars_beg_max = 0;
integertime_t ti_black_holes_end_min = max_nr_timesteps,
ti_black_holes_beg_max = 0;
integertime_t ti_sinks_end_min = max_nr_timesteps, ti_sinks_beg_max = 0;
#ifdef SWIFT_DEBUG_CHECKS
/* Ensure we have the right kind of cell. */
if (c->subtype != cell_subtype_void) {
error(
"Trying to update timesteps on cell which isn't a void cell! "
"(c->type=%s, c->subtype=%s)",
cellID_names[c->type], subcellID_names[c->subtype]);
}
#endif
/* Collect the values from the progeny. */
for (int k = 0; k < 8; k++) {
struct cell *cp = c->progeny[k];
if (cp != NULL) {
/* Recurse, if we're still in the void cell tree. */
if (cp->subtype == cell_subtype_void) {
zoom_void_timestep_collect_recursive(cp);
}
/* And update */
ti_hydro_end_min = min(ti_hydro_end_min, cp->hydro.ti_end_min);
ti_hydro_beg_max = max(ti_hydro_beg_max, cp->hydro.ti_beg_max);
ti_rt_end_min = min(cp->rt.ti_rt_end_min, ti_rt_end_min);
ti_rt_beg_max = max(cp->rt.ti_rt_beg_max, ti_rt_beg_max);
ti_grav_end_min = min(ti_grav_end_min, cp->grav.ti_end_min);
ti_grav_beg_max = max(ti_grav_beg_max, cp->grav.ti_beg_max);
ti_stars_end_min = min(ti_stars_end_min, cp->stars.ti_end_min);
ti_stars_beg_max = max(ti_stars_beg_max, cp->stars.ti_beg_max);
ti_black_holes_end_min =
min(ti_black_holes_end_min, cp->black_holes.ti_end_min);
ti_black_holes_beg_max =
max(ti_black_holes_beg_max, cp->black_holes.ti_beg_max);
ti_sinks_end_min = min(ti_sinks_end_min, cp->sinks.ti_end_min);
ti_sinks_beg_max = max(ti_sinks_beg_max, cp->sinks.ti_beg_max);
}
}
/* Store the collected values in the cell. */
c->hydro.ti_end_min = ti_hydro_end_min;
c->hydro.ti_beg_max = ti_hydro_beg_max;
c->rt.ti_rt_end_min = ti_rt_end_min;
c->rt.ti_rt_beg_max = ti_rt_beg_max;
c->grav.ti_end_min = ti_grav_end_min;
c->grav.ti_beg_max = ti_grav_beg_max;
c->stars.ti_end_min = ti_stars_end_min;
c->stars.ti_beg_max = ti_stars_beg_max;
c->black_holes.ti_end_min = ti_black_holes_end_min;
c->black_holes.ti_beg_max = ti_black_holes_beg_max;
c->sinks.ti_end_min = ti_sinks_end_min;
c->sinks.ti_beg_max = ti_sinks_beg_max;
}
/**
* @brief Collect the timestep from the zoom region to the void level.
*
* @param e The engine.
*/
void zoom_void_timestep_collect(struct engine *e) {
ticks tic = getticks();
/* Get the space and void cells. */
struct space *s = e->s;
struct zoom_region_properties *zoom_props = s->zoom_props;
const int nr_void_cells = zoom_props->nr_void_cells;
const int *void_cell_indices = zoom_props->void_cell_indices;
/* Loop over all void cells and collect the timesteps. */
for (int i = 0; i < nr_void_cells; i++) {
struct cell *c = &s->cells_top[void_cell_indices[i]];
zoom_void_timestep_collect_recursive(c);
}
if (e->verbose)
message("took %.3f %s.", clocks_from_ticks(getticks() - tic),
clocks_getunit());
}