Skip to content
Snippets Groups Projects
Commit dc4238ea authored by Paul Heinzlreiter's avatar Paul Heinzlreiter
Browse files

* debugging

parent f0037313
Branches
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 95 additions and 26 deletions
bh_tree_mpi/datastructures/BarnesHutTree.cpp
+ 10
− 0
View file @ dc4238ea
......@@ -132,6 +132,11 @@ namespace nbody {
//put all new bodies into fitting leaves, walk through tree and split
Node* current;
for (vector<Body>::iterator it = bodies.begin(); it != bodies.end(); it++) {
if (!this->getRootBB().contained(it->position)) {
cout << it - bodies.begin() << " Error: not contained: " << it->position[0] << " " << it->position[1] << " " << it->position[2] << endl;
}
}
for (vector<Body>::iterator it = bodies.begin(); it != bodies.end(); it++) {
current = this->nodes->next;
while (!current->leaf) {
......@@ -140,11 +145,16 @@ namespace nbody {
while (child != NULL && !child->getBB().contained(it->position)) {
child = child->nextSibling;
}
if (child != NULL) {
current = child;
}
/*
if (child == NULL) {
cout << it - bodies.begin() << " Error: not contained: " << it->position[0] << " " << it->position[1] << " " << it->position[2] << endl;
} else {
current = child;
}
*/
}
current->bodies.push_back(*it);
current->bodies.back().refinement = true;
......
bh_tree_mpi/datastructures/Body.cpp
+ 2
− 0
View file @ dc4238ea
......@@ -119,9 +119,11 @@ namespace nbody {
dvdt[i] = 1.0 / 6.0 * (a.dv[i] + 2.0f * (b.dv[i] + c.dv[i]) + d.dv[i]);
this->position[i] = this->position[i] + dxdt[i] * timestep;
this->velocity[i] = this->velocity[i] + dvdt[i] * timestep;
/*
if (fabs(dxdt[i] * timestep) > FLT_EPSILON) {
cout << "M: " << dxdt[i] * timestep << endl;
}
*/
}
//this->print();
......
bh_tree_mpi/datastructures/Box.cpp
+ 39
− 1
View file @ dc4238ea
......@@ -34,7 +34,9 @@ namespace nbody {
}
void Box::extendForBodies(vector<Body> bodies) {
if (bodies.empty()) return;
if (bodies.empty()) {
return;
}
for (vector<Body>::iterator it = bodies.begin(); it != bodies.end(); it++) {
for (int i = 0; i < 3; i++) {
if (it->position[i] < this->min[i]) {
......@@ -120,6 +122,9 @@ namespace nbody {
}
double Box::volume() {
if (!this->isValid()) {
return -1.0;
}
double result = 1.0;
for (int i = 0; i < 3; i++) {
......@@ -131,6 +136,9 @@ namespace nbody {
bool Box::overlapsSphere(double* sphereCenter, double sphereRadius) {
double dmin = 0.0;
if (!this->isValid()) {
return false;
}
for (int i = 0; i < 3; i++) {
if (sphereCenter[i] < this->min[i]) {
double dist = sphereCenter[i] - this->min[i];
......@@ -148,6 +156,9 @@ namespace nbody {
double Box::maxSidelength() {
double max = 0.0;
if (!this->isValid()) {
return -1.0;
}
for (int i = 0; i < 3; i++) {
if (this->max[i] - this->min[i] > max) {
max = this->max[i] - this->min[i];
......@@ -160,6 +171,9 @@ namespace nbody {
int inside = 0;
double nextPosition[3] = {position[0], position[1], position[2]};
if (!this->isValid()) {
return DBL_MAX;
}
for (int i = 0; i < 3; i++) {
if (nextPosition[i] < this->min[i]) {
nextPosition[i] = this->min[i];
......@@ -184,6 +198,9 @@ namespace nbody {
double Box::distanceToBox(Box box) {
double length = 0.0;
if (!this->isValid()) {
return DBL_MAX;
}
for (int i = 0; i < 3; i++) {
double elem;
......@@ -202,6 +219,9 @@ namespace nbody {
vector<Box> Box::octreeSplit() {
vector<Box> result;
if (!this->isValid()) {
return result;
}
for (unsigned int i = 0; i < 8; i++) {
Box current = *this;
......@@ -224,6 +244,9 @@ namespace nbody {
int longestIndex = -1;
double longestSide = -1.0;
if (!this->isValid()) {
return result;
}
for (int i = 0; i < 3; i++) {
if (this->max[i] - this->min[i] > longestSide) {
longestSide = this->max[i] - this->min[i];
......@@ -248,6 +271,15 @@ namespace nbody {
return true;
}
bool Box::isValid() {
for (int i = 0; i < 3; i++) {
if (this->max[i] < this->min[i]) {
return false;
}
}
return true;
}
void Box::print() {
cout << " min ";
for (int i = 0; i < 3; i++) {
......@@ -261,6 +293,9 @@ namespace nbody {
}
bool Box::contained(double* position) {
if (!this->isValid()) {
return false;
}
for (int i = 0; i < 3; i++) {
if (position[i] < this->min[i] || position[i] > this->max[i]) {
return false;
......@@ -270,6 +305,9 @@ namespace nbody {
}
void Box::extend(Box extender) {
if (!extender.isValid()) {
return;
}
for (int i = 0; i < 3; i++) {
if (this->min[i] > extender.min[i]) {
this->min[i] = extender.min[i];
......
bh_tree_mpi/datastructures/Box.hpp
+ 1
− 0
View file @ dc4238ea
......@@ -29,6 +29,7 @@ namespace nbody {
virtual double volume();
virtual double maxSidelength();
virtual bool isCorrect();
virtual bool isValid();
virtual void print();
virtual bool overlapsSphere(double* sphereCenter, double sphereRadius);
virtual double distanceToPosition(double* position);
......
bh_tree_mpi/datastructures/Tree.cpp
+ 16
− 4
View file @ dc4238ea
......@@ -20,7 +20,7 @@ namespace nbody {
}
void Tree::clean() {
//remove all nodes; refresh root node
//remove all nodes; refresh dummy first node
while (this->nodes->next != this->nodes) {
Node* node = this->nodes->next;
......@@ -57,9 +57,10 @@ namespace nbody {
body.resetAcceleration();
while (n != this->nodes) {
//cout << "n" << endl;
//cout << "n " << n << " " << n->leaf << endl;
if (n->sufficientForBody(body)) {
n->representative.accumulateForceOnto(body);
n = n->afterSubtree;
//cout << "rb" << endl;
} else if (n->leaf) {
//cout << "l" << endl;
......@@ -67,8 +68,11 @@ namespace nbody {
//cout << "lb" << endl;
it->accumulateForceOnto(body);
}
n = n->afterSubtree;
} else {
n = n->next;
}
n = n->afterSubtree;
//n = n->afterSubtree;
}
//body.print();
}
......@@ -121,6 +125,10 @@ namespace nbody {
vector<Body> result;
Node* current = this->nodes->next;
if (!current->bb.isValid()) {
//empty tree means no refinements
return result;
}
while (current != this->nodes) {
bool sufficient = current->sufficientForBox(domain);
......@@ -149,16 +157,20 @@ namespace nbody {
this->build(this->extractLocalBodies());
}
void Tree::advance() {
Box Tree::advance() {
Box bb;
for (Node* n = this->nodes->next; n != this->nodes; n = n->next) {
if (n->leaf) {
for (vector<Body>::iterator it = n->bodies.begin(); it != n->bodies.end(); it++) {
if (!it->refinement) {
it->integrate();
bb.extend(*it);
}
}
}
}
return bb;
}
Box Tree::getLocalBB() {
......
bh_tree_mpi/datastructures/Tree.hpp
+ 1
− 1
View file @ dc4238ea
......@@ -38,7 +38,7 @@ namespace nbody {
virtual Box getRootBB();
virtual void print();
virtual bool isCorrect();
virtual void advance();
virtual Box advance();
virtual Box getLocalBB();
};
......
bh_tree_mpi/mpimain.cpp
+ 2
− 1
View file @ dc4238ea
......@@ -11,7 +11,8 @@ int main(int argc, char* argv[]) {
MpiSimulation simulation(argc, argv);
simulation.distributeBodies();
for (int i = 0; i < 6; i++) {
for (int i = 0; i < 1; i++) {
//local tree is built and correct domains are distributed
simulation.runStep();
}
simulation.cleanup();
......
bh_tree_mpi/simulation/MpiSimulation.cpp
+ 24
− 19
View file @ dc4238ea
......@@ -4,6 +4,7 @@
#include <Box.hpp>
#include <Tree.hpp>
#include <Node.hpp>
#include <iostream>
#include "MpiSimulation.hpp"
namespace nbody {
......@@ -112,7 +113,7 @@ namespace nbody {
while (nodes.size() < (unsigned int) this->parallelSize) {
int mostBodiesIndex = 0;
for (unsigned int i = 0; i < nodes.size(); i++) {
for (unsigned int i = 1; i < nodes.size(); i++) {
if (nodes[i].getBodies().size() > nodes[mostBodiesIndex].getBodies().size()) {
mostBodiesIndex = i;
}
......@@ -120,16 +121,13 @@ namespace nbody {
vector<Box> subdomains = nodes[mostBodiesIndex].getBB().splitLongestSide();
vector<Body> buf = nodes[mostBodiesIndex].getBodies();
Node n(NULL);
Box bb;
n.setBodies(subdomains[0].extractBodies(buf));
bb.extendForBodies(n.getBodies());
n.setBB(bb);
n.setBB(subdomains[0]);
nodes.insert(nodes.begin() + mostBodiesIndex, n);
n = Node(NULL);
n.setBodies(subdomains[1].extractBodies(buf));
bb.extendForBodies(n.getBodies());
n.setBB(bb);
n.setBB(subdomains[1]);
nodes.insert(nodes.begin() + mostBodiesIndex, n);
nodes.erase(nodes.begin() + mostBodiesIndex + 2);
}
......@@ -143,13 +141,24 @@ namespace nbody {
this->comms.push_back(MpiBodyComm(&this->bodyType));
this->comms[0].recvBlocking(0, this->bodies);
}
this->tree.build(this->bodies);
this->domains[this->parallelRank] = this->tree.getRootBB();
//this->tree.isCorrect();
this->distributeDomains();
Box overallDomain;
for (int i = 0; i < this->parallelSize; i++) {
overallDomain.extend(this->domains[i]);
}
this->tree.build(this->bodies, overallDomain);
this->tree.isCorrect();
}
void MpiSimulation::distributeDomains() {
MPI_Allgather(&this->domains[this->parallelRank], 1, this->boxType, &this->domains[0], 1, this->boxType, MPI_COMM_WORLD);
/*
for (int i = 0; i < this->parallelSize; i++) {
this->domains[i].print();
cout << "-" << endl;
}
*/
}
void MpiSimulation::distributeLETs() {
......@@ -187,25 +196,21 @@ namespace nbody {
void MpiSimulation::runStep() {
Box overallDomain;
this->distributeLETs();
this->tree.computeForces();
this->domains[this->parallelRank] = this->tree.advance();
this->distributeDomains();
//this->domains[this->mpiRank].print();
for (int i = 0; i < this->parallelSize; i++) {
overallDomain.extend(this->domains[i]);
}
/*
if (this->mpiRank == 0) {
overallDomain.print();
}
*/
this->tree.rebuild(overallDomain);
this->distributeLETs();
this->tree.computeForces();
this->tree.advance();
this->domains[this->parallelRank] = this->tree.getLocalBB();
/*
if (!this->tree.isCorrect()) {
cout << "tree not correct" << endl;
cout << "Tree wrong" << endl;
}
*/
}
}
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment