mpi_complete_tree_debug/
mpi_complete_tree_debug.rs

1//! Test the computation of a complete octree.
2
3use bempp_octree::{
4    morton::MortonKey,
5    octree::{is_complete_linear_and_balanced, KeyType, Octree},
6    tools::{gather_to_all, generate_random_points},
7};
8use itertools::Itertools;
9use mpi::traits::Communicator;
10use rand::prelude::*;
11use rand_chacha::ChaCha8Rng;
12
13pub fn main() {
14    // Initialise MPI
15    let universe = mpi::initialize().unwrap();
16
17    // Get the world communicator
18    let comm = universe.world();
19
20    // Initialise a seeded Rng.
21    let mut rng = ChaCha8Rng::seed_from_u64(comm.rank() as u64);
22
23    // Create `npoints` per rank.
24    let npoints = 10000;
25
26    // Generate random points.
27
28    let mut points = generate_random_points(npoints, &mut rng, &comm);
29    // Make sure that the points live on the unit sphere.
30    for point in points.iter_mut() {
31        let len = point.coords()[0] * point.coords()[0]
32            + point.coords()[1] * point.coords()[1]
33            + point.coords()[2] * point.coords()[2];
34        let len = len.sqrt();
35        point.coords_mut()[0] /= len;
36        point.coords_mut()[1] /= len;
37        point.coords_mut()[2] /= len;
38    }
39
40    let tree = Octree::new(&points, 15, 50, &comm);
41
42    // We now check that each node of the tree has all its neighbors available.
43
44    let leaf_tree = tree.leaf_keys();
45    let all_keys = tree.all_keys();
46
47    assert!(is_complete_linear_and_balanced(leaf_tree, &comm));
48    for &key in leaf_tree {
49        // We only check interior keys. Leaf keys may not have a neighbor
50        // on the same level.
51        let mut parent = key.parent();
52        while parent.level() > 0 {
53            // Check that the key itself is there.
54            assert!(all_keys.contains_key(&key));
55            // Check that all its neighbours are there.
56            for neighbor in parent.neighbours().iter().filter(|&key| key.is_valid()) {
57                assert!(all_keys.contains_key(neighbor));
58            }
59            parent = parent.parent();
60            // Check that the parent is there.
61            assert!(all_keys.contains_key(&parent));
62        }
63    }
64
65    // At the end check that the root of the tree is also contained.
66    assert!(all_keys.contains_key(&MortonKey::root()));
67
68    // Count the number of ghosts on each rank
69    // Count the number of global keys on each rank.
70
71    // Assert that all ghosts are from a different rank and count them.
72
73    let nghosts = all_keys
74        .iter()
75        .filter_map(|(_, &value)| {
76            if let KeyType::Ghost(rank) = value {
77                assert!(rank != comm.size() as usize);
78                Some(rank)
79            } else {
80                None
81            }
82        })
83        .count();
84
85    if comm.size() == 1 {
86        assert_eq!(nghosts, 0);
87    } else {
88        assert!(nghosts > 0);
89    }
90
91    let nglobal = all_keys
92        .iter()
93        .filter(|(_, &value)| matches!(value, KeyType::Global))
94        .count();
95
96    // Assert that all globals across all ranks have the same count.
97
98    let nglobals = gather_to_all(std::slice::from_ref(&nglobal), &comm);
99
100    assert_eq!(nglobals.iter().unique().count(), 1);
101
102    // Check that the points are associated with the correct leaf keys.
103    let mut npoints = 0;
104    let leaf_point_map = tree.leaf_keys_to_local_point_indices();
105
106    for (key, point_indices) in leaf_point_map {
107        for &index in point_indices {
108            assert!(key.is_ancestor(tree.point_keys()[index]));
109        }
110        npoints += point_indices.len();
111    }
112
113    // Make sure that the number of points and point keys lines up
114    // with the points stored for each leaf key.
115    assert_eq!(npoints, tree.points().len());
116    assert_eq!(npoints, tree.point_keys().len());
117
118    // Check the neighbour relationships.
119
120    let all_neighbours = tree.neighbour_map();
121    let all_keys = tree.all_keys();
122
123    for (key, key_type) in all_keys {
124        // Ghost keys should not be in the neighbour map.
125        match key_type {
126            KeyType::Ghost(_) => assert!(!all_neighbours.contains_key(key)),
127            _ => {
128                // If it is not a ghost the key should be in the neighbour map.
129                assert!(all_neighbours.contains_key(key));
130            }
131        }
132    }
133
134    if comm.rank() == 0 {
135        println!("No errors were found in setting up tree.");
136    }
137}
mpi_complete_tree_debug/mpi_complete_tree_debug.rs

mpi_complete_tree_debug/
mpi_complete_tree_debug.rs
