Skip to content

Commit

Permalink
offline computation of cost statistics (#1383)
Browse files Browse the repository at this point in the history
* offline computation of cost statistics

* fix bug when calling create_structure_and_set_materials with chunk_layout

* oops, forget to add back in support for symmetries

* better refactorization of choose_chunkdivision

* handle Python default for zero num_chunks
  • Loading branch information
oskooi authored Oct 7, 2020
1 parent 5ceb8ca commit 39a8f67
Show file tree
Hide file tree
Showing 5 changed files with 101 additions and 59 deletions.
12 changes: 11 additions & 1 deletion python/meep.i
Original file line number Diff line number Diff line change
Expand Up @@ -1817,7 +1817,8 @@ meep::structure *create_structure_and_set_materials(vector3 cell_size,
meep_geom::material_type_list extra_materials,
bool split_chunks_evenly,
bool set_materials,
meep::structure *existing_s) {
meep::structure *existing_s,
bool output_chunk_costs) {
// Initialize fragment_stats static members (used for creating chunks in choose_chunkdivision)
meep_geom::fragment_stats::geom = gobj_list;
meep_geom::fragment_stats::dft_data_list = dft_data_list_;
Expand All @@ -1832,6 +1833,15 @@ meep::structure *create_structure_and_set_materials(vector3 cell_size,
meep_geom::fragment_stats::split_chunks_evenly = split_chunks_evenly;
meep_geom::fragment_stats::init_libctl(_default_material, _ensure_periodicity,
&gv, cell_size, center, &gobj_list);

if (output_chunk_costs) {
meep::volume thev = gv.surroundings();
std::vector<grid_volume> chunk_vols = meep::choose_chunkdivision(gv, thev, num_chunks, sym);
for (size_t i = 0; i < chunk_vols.size(); ++i)
master_printf("CHUNK:, %2zu, %f\n",i,chunk_vols[i].get_cost());
return NULL;
}

meep::structure *s;
if (existing_s) {
s = existing_s;
Expand Down
34 changes: 19 additions & 15 deletions python/simulation.py
Original file line number Diff line number Diff line change
Expand Up @@ -1659,18 +1659,17 @@ def _init_structure(self, k=False):

if self._output_stats and isinstance(self.default_material, mp.Medium) and verbosity.meep > 0:
stats = self._compute_fragment_stats(gv)
print("STATS: aniso_eps: {}".format(stats.num_anisotropic_eps_pixels))
print("STATS: anis_mu: {}".format(stats.num_anisotropic_mu_pixels))
print("STATS: nonlinear: {}".format(stats.num_nonlinear_pixels))
print("STATS: susceptibility: {}".format(stats.num_susceptibility_pixels))
print("STATS: nonzero_cond: {}".format(stats.num_nonzero_conductivity_pixels))
print("STATS: pml_1d: {}".format(stats.num_1d_pml_pixels))
print("STATS: pml_2d: {}".format(stats.num_2d_pml_pixels))
print("STATS: pml_3d: {}".format(stats.num_3d_pml_pixels))
print("STATS: dft: {}".format(stats.num_dft_pixels))
print("STATS: total_pixels: {}".format(stats.num_pixels_in_box))
print("STATS: num_cores: {}".format(mp.count_processors()))
sys.exit(0)
print("FRAGMENT:, aniso_eps:, {}".format(stats.num_anisotropic_eps_pixels))
print("FRAGMENT:, aniso_mu:, {}".format(stats.num_anisotropic_mu_pixels))
print("FRAGMENT:, nonlinear:, {}".format(stats.num_nonlinear_pixels))
print("FRAGMENT:, susceptibility:, {}".format(stats.num_susceptibility_pixels))
print("FRAGMENT:, conductivity:, {}".format(stats.num_nonzero_conductivity_pixels))
print("FRAGMENT:, pml_1d:, {}".format(stats.num_1d_pml_pixels))
print("FRAGMENT:, pml_2d:, {}".format(stats.num_2d_pml_pixels))
print("FRAGMENT:, pml_3d:, {}".format(stats.num_3d_pml_pixels))
print("FRAGMENT:, dft:, {}".format(stats.num_dft_pixels))
print("FRAGMENT:, total_pixels:, {}".format(stats.num_pixels_in_box))
print("FRAGMENT:, procs:, {}".format(mp.count_processors()))

fragment_vols = self._make_fragment_lists(gv)
self.dft_data_list = fragment_vols[0]
Expand Down Expand Up @@ -1703,8 +1702,12 @@ def _init_structure(self, k=False):
self.extra_materials,
self.split_chunks_evenly,
False if self.chunk_layout else True,
None
)
None,
True if self._output_stats is not None else False
)

if self._output_stats is not None:
sys.exit(0)

if self.chunk_layout:
self.load_chunk_layout(br, self.chunk_layout)
Expand Down Expand Up @@ -1850,7 +1853,8 @@ def set_materials(self, geometry=None, default_material=None):
self.extra_materials,
self.split_chunks_evenly,
True,
self.structure
self.structure,
False
)

def dump_structure(self, fname):
Expand Down
6 changes: 6 additions & 0 deletions src/meep.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -807,6 +807,12 @@ class structure {
void write_susceptibility_params(h5file *file, const char *dname, int EorH);
};

// defined in structure.cpp
std::vector<grid_volume> choose_chunkdivision(grid_volume &gv,
volume &v,
int num_chunks,
const symmetry &s);

class src_vol;
class fields;
class fields_chunk;
Expand Down
20 changes: 10 additions & 10 deletions src/meepgeom.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2159,16 +2159,16 @@ double fragment_stats::cost() const {

void fragment_stats::print_stats() const {
master_printf("Fragment stats\n");
master_printf(" num_anisotropic_eps_pixels: %zd\n", num_anisotropic_eps_pixels);
master_printf(" num_anisotropic_mu_pixels: %zd\n", num_anisotropic_mu_pixels);
master_printf(" num_nonlinear_pixels: %zd\n", num_nonlinear_pixels);
master_printf(" num_susceptibility_pixels: %zd\n", num_susceptibility_pixels);
master_printf(" num_nonzero_conductivity_pixels: %zd\n", num_nonzero_conductivity_pixels);
master_printf(" num_1d_pml_pixels: %zd\n", num_1d_pml_pixels);
master_printf(" num_2d_pml_pixels: %zd\n", num_2d_pml_pixels);
master_printf(" num_3d_pml_pixels: %zd\n", num_3d_pml_pixels);
master_printf(" num_dft_pixels: %zd\n", num_dft_pixels);
master_printf(" num_pixels_in_box: %zd\n", num_pixels_in_box);
master_printf(" anisotropic_eps: %zu\n", num_anisotropic_eps_pixels);
master_printf(" anisotropic_mu: %zu\n", num_anisotropic_mu_pixels);
master_printf(" nonlinear: %zu\n", num_nonlinear_pixels);
master_printf(" susceptibility: %zu\n", num_susceptibility_pixels);
master_printf(" conductivity: %zu\n", num_nonzero_conductivity_pixels);
master_printf(" pml_1d: %zu\n", num_1d_pml_pixels);
master_printf(" pml_2d: %zu\n", num_2d_pml_pixels);
master_printf(" pml_3d: %zu\n", num_3d_pml_pixels);
master_printf(" dft: %zu\n", num_dft_pixels);
master_printf(" pixels_in_box: %zu\n", num_pixels_in_box);
master_printf(" box.low: {%f, %f, %f}\n", box.low.x, box.low.y, box.low.z);
master_printf(" box.high: {%f, %f, %f}\n\n", box.high.x, box.high.y, box.high.z);
}
Expand Down
88 changes: 55 additions & 33 deletions src/structure.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -122,25 +122,72 @@ static void split_by_cost(std::vector<int> factors, grid_volume gvol,

void structure::choose_chunkdivision(const grid_volume &thegv, int desired_num_chunks,
const boundary_region &br, const symmetry &s) {
user_volume = thegv;
if (desired_num_chunks == 0) desired_num_chunks = count_processors();

if (thegv.dim == Dcyl && thegv.get_origin().r() < 0) abort("r < 0 origins are not supported");

user_volume = thegv;
gv = thegv;
v = gv.surroundings();
S = s;
a = gv.a;
dt = Courant / a;

// create the chunks:
std::vector<grid_volume> chunk_volumes = meep::choose_chunkdivision(gv, v, desired_num_chunks, s);

int my_num_chunks = chunk_volumes.size();

// initialize effort volumes
num_effort_volumes = 1;
effort_volumes = new grid_volume[num_effort_volumes];
effort_volumes[0] = gv;
effort = new double[num_effort_volumes];
effort[0] = 1.0;

// Next, add effort volumes for PML boundary regions:
br.apply(this);

// Break off PML regions into their own chunks
num_chunks = 0;
chunks = new structure_chunk_ptr[my_num_chunks * num_effort_volumes];
for (size_t i = 0, stop = chunk_volumes.size(); i < stop; ++i) {
const int proc = i * count_processors() / my_num_chunks;
for (int j = 0; j < num_effort_volumes; ++j) {
grid_volume vc;
if (chunk_volumes[i].intersect_with(effort_volumes[j], &vc)) {
chunks[num_chunks] = new structure_chunk(vc, v, Courant, proc);
br.apply(this, chunks[num_chunks++]);
}
}
}

check_chunks();

if (meep_geom::fragment_stats::resolution != 0) {
// Save cost of each chunk's grid_volume
for (int i = 0; i < num_chunks; ++i) {
chunks[i]->cost = chunks[i]->gv.get_cost();
}
}

}

std::vector<grid_volume> choose_chunkdivision(grid_volume &gv, volume &v, int desired_num_chunks,
const symmetry &S) {

if (desired_num_chunks == 0) desired_num_chunks = count_processors();
if (gv.dim == Dcyl && gv.get_origin().r() < 0) abort("r < 0 origins are not supported");

// First, reduce overall grid_volume gv by symmetries:
if (S.multiplicity() > 1) {
bool break_this[3];
for (int dd = 0; dd < 3; dd++) {
const direction d = (direction)dd;
break_this[dd] = false;
for (int n = 0; n < S.multiplicity(); n++)
if (has_direction(thegv.dim, d) &&
if (has_direction(gv.dim, d) &&
(S.transform(d, n).d != d || S.transform(d, n).flipped)) {
if (thegv.num_direction(d) & 1 && !break_this[d] && verbosity > 0)
if (gv.num_direction(d) & 1 && !break_this[d] && verbosity > 0)
master_printf("Padding %s to even number of grid points.\n", direction_name(d));
break_this[dd] = true;
}
Expand All @@ -164,15 +211,12 @@ void structure::choose_chunkdivision(const grid_volume &thegv, int desired_num_c
}

// initialize effort volumes
num_effort_volumes = 1;
effort_volumes = new grid_volume[num_effort_volumes];
int num_effort_volumes = 1;
grid_volume *effort_volumes = new grid_volume[num_effort_volumes];
effort_volumes[0] = gv;
effort = new double[num_effort_volumes];
double *effort = new double[num_effort_volumes];
effort[0] = 1.0;

// Next, add effort volumes for PML boundary regions:
br.apply(this);

std::vector<int> prime_factors = get_prime_factors(desired_num_chunks);

// We may have to use a different number of chunks than the user requested
Expand All @@ -184,11 +228,8 @@ void structure::choose_chunkdivision(const grid_volume &thegv, int desired_num_c
meep_geom::fragment_stats::split_chunks_evenly ? desired_num_chunks : adjusted_num_chunks;

// Finally, create the chunks:
num_chunks = 0;
chunks = new structure_chunk_ptr[my_num_chunks * num_effort_volumes];
std::vector<grid_volume> chunk_volumes;

bool by_cost = false;
if (meep_geom::fragment_stats::resolution == 0 ||
meep_geom::fragment_stats::split_chunks_evenly) {
if (verbosity > 0 && my_num_chunks > 1)
Expand All @@ -203,28 +244,9 @@ void structure::choose_chunkdivision(const grid_volume &thegv, int desired_num_c
if (verbosity > 0 && my_num_chunks > 1)
master_printf("Splitting into %d chunks by cost\n", my_num_chunks);
split_by_cost(prime_factors, gv, chunk_volumes);
by_cost = true;
}

// Break off PML regions into their own chunks
for (size_t i = 0, stop = chunk_volumes.size(); i < stop; ++i) {
const int proc = i * count_processors() / my_num_chunks;
for (int j = 0; j < num_effort_volumes; ++j) {
grid_volume vc;
if (chunk_volumes[i].intersect_with(effort_volumes[j], &vc)) {
chunks[num_chunks] = new structure_chunk(vc, v, Courant, proc);
br.apply(this, chunks[num_chunks++]);
}
}
}
check_chunks();

if (by_cost) {
// Save cost of each chunk's grid_volume
for (int i = 0; i < num_chunks; ++i) {
chunks[i]->cost = chunks[i]->gv.get_cost();
}
}
return chunk_volumes;
}

double structure::estimated_cost(int process) {
Expand Down

0 comments on commit 39a8f67

Please sign in to comment.