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

* debugging

parent 61198dd3
Branches
No related merge requests found
......@@ -5,8 +5,7 @@
namespace nbody {
using namespace std;
BarnesHutTree::BarnesHutTree() {
this->maxLeafBodies = 16;
BarnesHutTree::BarnesHutTree(int parallelId) : Tree(parallelId) {
}
BarnesHutTree::~BarnesHutTree() {
......
......@@ -16,7 +16,7 @@ namespace nbody {
virtual void update();
static void split(Node* current);
public:
BarnesHutTree();
BarnesHutTree(int parallelId);
virtual ~BarnesHutTree();
virtual void build(vector<Body> bodies);
virtual void build(vector<Body> bodies, Box domain);
......
......@@ -3,6 +3,9 @@
#include <iostream>
#include <cmath>
#include <cfloat>
#include <climits>
#include <cstring>
#include <cstdio>
namespace nbody {
using namespace std;
......@@ -24,6 +27,7 @@ namespace nbody {
};
Body::Body() {
this->id = ULONG_MAX;
for (int i = 0; i < 3; i++) {
this->position[i] = 0.0;
this->velocity[i] = 0.0;
......@@ -34,6 +38,7 @@ namespace nbody {
}
Body::Body(double positionX, double positionY, double positionZ) {
this->id = ULONG_MAX;
for (int i = 0; i < 3; i++) {
this->velocity[i] = 0.0;
this->acceleration[i] = 0.0;
......@@ -46,6 +51,7 @@ namespace nbody {
}
Body::Body(double positionX, double positionY, double positionZ, double velocityX, double velocityY, double velocityZ) {
this->id = ULONG_MAX;
for (int i = 0; i < 3; i++) {
this->acceleration[i] = 0.0;
}
......@@ -60,6 +66,7 @@ namespace nbody {
}
Body::Body(double positionX, double positionY, double positionZ, double weight) {
this->id = ULONG_MAX;
for (int i = 0; i < 3; i++) {
this->velocity[i] = 0.0;
this->acceleration[i] = 0.0;
......@@ -72,6 +79,7 @@ namespace nbody {
}
Body::Body(double positionX, double positionY, double positionZ, double velocityX, double velocityY, double velocityZ, double weight) {
this->id = ULONG_MAX;
for (int i = 0; i < 3; i++) {
this->acceleration[i] = 0.0;
}
......@@ -89,6 +97,10 @@ namespace nbody {
}
void Body::setId(unsigned long id) {
this->id = id;
}
Derivative Body::evaluate(double dt, Derivative d) {
double velocity[3];
Derivative output;
......@@ -153,21 +165,49 @@ namespace nbody {
}
}
bool Body::valid() {
bool validDouble(double val) {
char buf[128];
for (int i = 0; i < 3; i++) {
if (this->position[i] < -FLT_MAX || this->position[i] > FLT_MAX) {
sprintf(buf, "%f", val);
for (int i = 0; i < strlen(buf); i++) {
if (!((buf[i] >= '0' && buf[i] <= '9') || buf[i] == '.' || buf[i] == '-')) {
return false;
}
}
return true;
}
bool Body::valid() {
for (int i = 0; i < 3; i++) {
if (!validDouble(this->position[i])) return false;
if (!validDouble(this->velocity[i])) return false;
if (!validDouble(this->acceleration[i])) return false;
}
if (!validDouble(this->mass)) return false;
return this->mass >= 0.0;
}
void Body::print() {
cout << "Position: " << this->position[0] << " " << this->position[1] << " " << this->position[2] << endl;
cout << "Velocity: " << this->velocity[0] << " " << this->velocity[1] << " " << this->velocity[2] << endl;
cout << "Acceleration: " << this->acceleration[0] << " " << this->acceleration[1] << " " << this->acceleration[2] << endl;
cout << "Mass: " << this->mass << endl;
cout << "Refinement: " << this->refinement << endl << endl;
bool Body::valid(vector<Body> bodies) {
for (vector<Body>::iterator it = bodies.begin(); it != bodies.end(); it++) {
if (!it->valid()) {
return false;
}
}
return true;
}
void Body::print(int parallelId) {
cout << parallelId << " Id: " << this->id << endl;
cout << parallelId << " Position: " << this->position[0] << " " << this->position[1] << " " << this->position[2] << endl;
cout << parallelId << " Velocity: " << this->velocity[0] << " " << this->velocity[1] << " " << this->velocity[2] << endl;
cout << parallelId << " Acceleration: " << this->acceleration[0] << " " << this->acceleration[1] << " " << this->acceleration[2] << endl;
cout << parallelId << " Mass: " << this->mass << endl;
cout << parallelId << " Refinement: " << this->refinement << endl << endl;
}
void Body::print(int parallelId, vector<Body> bodies) {
for (vector<Body>::iterator it = bodies.begin(); it != bodies.end(); it++) {
it->print(parallelId);
}
}
}
#ifndef BODY_HPP
#define BODY_HPP
#include <vector>
namespace nbody {
using namespace std;
static const double timestep = 1.0;
class Derivative;
......@@ -13,6 +17,7 @@ namespace nbody {
Derivative evaluate(double dt, Derivative d);
public:
unsigned long id;
double position[3];
double velocity[3];
double acceleration[3];
......@@ -29,7 +34,10 @@ namespace nbody {
virtual void integrate();
virtual void accumulateForceOnto(Body& from);
virtual bool valid();
virtual void print();
virtual void print(int parallelId);
static void print(int parallelId, vector<Body> bodies);
virtual void setId(unsigned long id);
static bool valid(vector<Body> bodies);
};
}
......
......@@ -165,7 +165,7 @@ namespace nbody {
this->bb.print(parallelId);
for (vector<Body>::iterator it = this->bodies.begin(); it != this->bodies.end(); it++) {
cout << " ";
it->print();
it->print(parallelId);
}
}
......
......@@ -8,10 +8,11 @@
namespace nbody {
using namespace std;
Tree::Tree() {
Tree::Tree(int parallelId) {
//insert dummy root node
this->nodes = new Node(this);
this->maxLeafBodies = 16;
this->parallelId = parallelId;
}
Tree::~Tree() {
......@@ -93,6 +94,7 @@ namespace nbody {
string line;
ifstream infile(filename.c_str(), ifstream::in);
double mass, px, py, pz, vx, vy, vz;
unsigned long id = 0;
while (std::getline(infile, line)) {
istringstream iss(line);
......@@ -100,6 +102,7 @@ namespace nbody {
vx = vy = vz = 0.0;
iss >> mass >> px >> py >> pz >> vx >> vy >> vz;
Body b(px, py, pz, vx, vy, vz, mass);
b.setId(id++);
result.push_back(b);
}
infile.close();
......@@ -175,7 +178,7 @@ namespace nbody {
for (vector<Body>::iterator it = n->bodies.begin(); it != n->bodies.end(); it++) {
if (!it->refinement) {
it->integrate();
it->print();
//it->print();
bb.extend(*it);
}
}
......
......@@ -19,10 +19,11 @@ namespace nbody {
protected:
Node* nodes;
unsigned int maxLeafBodies;
int parallelId;
public:
static vector<Body> dubinskiParse(string filename);
Tree();
Tree(int parallelId);
virtual ~Tree();
virtual void clean();
virtual void build(vector<Body> bodies) = 0;
......
#include <BarnesHutTree.hpp>
#include <MpiSimulation.hpp>
#include <iostream>
#include <cstdio>
#include <cstddef>
#include <vector>
#include <mpi.h>
using namespace nbody;
using namespace std;
class Car {
public:
int shifts;
int topSpeed;
double pos[3];
};
int main(int argc, char* argv[]) {
MPI_Init(&argc, &argv);
MpiSimulation simulation(argc, argv);
simulation.distributeBodies();
for (int i = 0; i < 6; i++) {
simulation.runStep();
}
simulation.cleanup();
MPI_Finalize();
const int tag = 13;
int size, rank;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
if (size < 2) {
fprintf(stderr,"Requires at least two processes.\n");
}
/* create a type for struct car */
const int nitems = 3;
int blocklengths[3] = {1, 1, 3};
MPI_Datatype types[3] = {MPI_INT, MPI_INT, MPI_DOUBLE};
MPI_Datatype mpi_car_type;
MPI_Aint offsets[3];
offsets[0] = offsetof(Car, shifts);
offsets[1] = offsetof(Car, topSpeed);
offsets[2] = offsetof(Car, pos);
MPI_Type_create_struct(nitems, blocklengths, offsets, types, &mpi_car_type);
MPI_Type_commit(&mpi_car_type);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
if (rank == 0) {
vector<Car> send;
send.push_back(Car());
send.push_back(Car());
send[0].shifts = 4;
send[0].topSpeed = 100;
send[0].pos[0] = 0.5;
send[0].pos[1] = 0.6;
send[0].pos[2] = 0.7;
send[1].shifts = 5;
send[1].topSpeed = 200;
send[1].pos[0] = 1.5;
send[1].pos[1] = 1.6;
send[1].pos[2] = 1.7;
const int dest = 1;
MPI_Send(&send[0], 2, mpi_car_type, dest, tag, MPI_COMM_WORLD);
printf("Rank %d: sent structure car\n", rank);
}
if (rank == 1) {
MPI_Status status;
const int src=0;
vector<Car> recv;
recv.reserve(2);
MPI_Recv(&recv[0], 2, mpi_car_type, src, tag, MPI_COMM_WORLD, &status);
printf("Rank %d: Received: shifts = %d topSpeed = %d | %f | %f | %f\n", rank,
recv[0].shifts, recv[0].topSpeed, recv[0].pos[0], recv[0].pos[1], recv[0].pos[2]);
printf("Rank %d: Received: shifts = %d topSpeed = %d | %f | %f | %f\n", rank,
recv[1].shifts, recv[1].topSpeed, recv[1].pos[0], recv[1].pos[1], recv[1].pos[2]);
}
MPI_Type_free(&mpi_car_type);
MPI_Finalize();
return 0;
}
......@@ -2,6 +2,7 @@
#include <MpiSimulation.hpp>
#include <iostream>
#include <mpi.h>
#include <cstdlib>
#include <Tree.hpp>
using namespace nbody;
......@@ -9,17 +10,77 @@ using namespace std;
int main(int argc, char* argv[]) {
MPI_Init(&argc, &argv);
MpiSimulation simulation(argc, argv);
//MpiSimulation simulation(argc, argv);
Body data[4];
int mysize, myrank;
Body b[2];
int bodyBlocklengths[8] = {1, 1, 3, 3, 3, 1, 1, 1};
MPI_Datatype bodyDatatypes[8] = {MPI_LB, MPI_UNSIGNED_LONG, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_UB};
MPI_Aint bodyDisplacements[8];
MPI_Aint baseAdress;
MPI_Datatype mpi_body;
simulation.distributeBodies();
simulation.distributeDomains();
simulation.buildTree();
MPI_Get_address(&(b[0]), &baseAdress);
bodyDisplacements[0] = baseAdress - baseAdress;
MPI_Get_address(&(b[0].id), &bodyDisplacements[1]);
bodyDisplacements[1] -= baseAdress;
MPI_Get_address(&(b[0].position[0]), &bodyDisplacements[2]);
bodyDisplacements[2] -= baseAdress;
MPI_Get_address(&(b[0].velocity[0]), &bodyDisplacements[3]);
bodyDisplacements[3] -= baseAdress;
MPI_Get_address(&(b[0].acceleration[0]), &bodyDisplacements[4]);
bodyDisplacements[4] -= baseAdress;
MPI_Get_address(&(b[0].mass), &bodyDisplacements[5]);
bodyDisplacements[5] -= baseAdress;
MPI_Get_address(&(b[0].refinement), &bodyDisplacements[6]);
bodyDisplacements[6] -= baseAdress;
MPI_Get_address(&(b[7]), &bodyDisplacements[7]);
bodyDisplacements[7] -= baseAdress;
MPI_Type_create_struct(8, bodyBlocklengths, bodyDisplacements, bodyDatatypes, &mpi_body);
MPI_Type_commit(&mpi_body);
MPI_Comm_size(MPI_COMM_WORLD, &mysize);
MPI_Comm_rank(MPI_COMM_WORLD, &myrank);
if (myrank == 0) {
srandom(42);
for (int i = 0; i < 4; i++) {
Body b(random() / random(), random() / random(), random() / random(), random() / random(), random() / random(), random() / random(), random() / random());
b.setId(i);
data[i] = b;
}
MPI_Send(&data[0], 4, mpi_body, 1, 0, MPI_COMM_WORLD);
} else if (myrank == 1) {
MPI_Status status;
MPI_Recv(&data[0], 4, mpi_body, 0, 0, MPI_COMM_WORLD, &status);
}
MPI_Finalize();
/*
if (myrank == 0) {
srandom(42);
for (int i = 0; i < 4; i++) {
Body b(random() / random(), random() / random(), random() / random(), random() / random(), random() / random(), random() / random(), random() / random());
b.setId(i);
data[i] = b;
}
//MPI_Send(&data[0], 4, *simulation.getDatatype(), 1, 0, MPI_COMM_WORLD);
} else if (myrank == 1) {
MPI_Status status;
//MPI_Recv(&data[0], 4, *simulation.getDatatype(), 0, 0, MPI_COMM_WORLD, &status);
}
*/
//simulation.distributeBodies();
//simulation.distributeDomains();
//simulation.buildTree();
/*
for (int i = 0; i < 1; i++) {
//local tree is built and correct domains are distributed
simulation.runStep();
}
simulation.cleanup();
MPI_Finalize();
*/
//simulation.cleanup();
//MPI_Finalize();
return 0;
}
......@@ -7,6 +7,7 @@
#include <Tree.hpp>
#include <Node.hpp>
#include <iostream>
#include <BarnesHutTree.hpp>
#include "MpiSimulation.hpp"
namespace nbody {
......@@ -14,9 +15,9 @@ namespace nbody {
MpiSimulation::MpiSimulation(int& argc, char**& argv) {
Body b[2];
int bodyBlocklengths[7] = {1, 3, 3, 3, 1, 1, 1};
MPI_Datatype bodyDatatypes[7] = {MPI_LB, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_UB};
MPI_Aint bodyDisplacements[7];
int bodyBlocklengths[8] = {1, 1, 3, 3, 3, 1, 1, 1};
MPI_Datatype bodyDatatypes[8] = {MPI_LB, MPI_UNSIGNED_LONG, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_UB};
MPI_Aint bodyDisplacements[8];
Box box[2];
int boxBlocklengths[4] = {1, 3, 3, 1};
MPI_Datatype boxDatatypes[4] = {MPI_LB, MPI_DOUBLE, MPI_DOUBLE, MPI_UB};
......@@ -27,19 +28,21 @@ namespace nbody {
this->boxType = MPI_DATATYPE_NULL;
MPI_Get_address(&(b[0]), &baseAdress);
bodyDisplacements[0] = baseAdress - baseAdress;
MPI_Get_address(&(b[0].position[0]), &bodyDisplacements[1]);
MPI_Get_address(&(b[0].id), &bodyDisplacements[1]);
bodyDisplacements[1] -= baseAdress;
MPI_Get_address(&(b[0].velocity[0]), &bodyDisplacements[2]);
MPI_Get_address(&(b[0].position[0]), &bodyDisplacements[2]);
bodyDisplacements[2] -= baseAdress;
MPI_Get_address(&(b[0].acceleration[0]), &bodyDisplacements[3]);
MPI_Get_address(&(b[0].velocity[0]), &bodyDisplacements[3]);
bodyDisplacements[3] -= baseAdress;
MPI_Get_address(&(b[0].mass), &bodyDisplacements[4]);
MPI_Get_address(&(b[0].acceleration[0]), &bodyDisplacements[4]);
bodyDisplacements[4] -= baseAdress;
MPI_Get_address(&(b[0].refinement), &bodyDisplacements[5]);
MPI_Get_address(&(b[0].mass), &bodyDisplacements[5]);
bodyDisplacements[5] -= baseAdress;
MPI_Get_address(&(b[1]), &bodyDisplacements[6]);
MPI_Get_address(&(b[0].refinement), &bodyDisplacements[6]);
bodyDisplacements[6] -= baseAdress;
MPI_Type_create_struct(7, bodyBlocklengths, bodyDisplacements, bodyDatatypes, &this->bodyType);
MPI_Get_address(&(b[7]), &bodyDisplacements[7]);
bodyDisplacements[7] -= baseAdress;
MPI_Type_create_struct(8, bodyBlocklengths, bodyDisplacements, bodyDatatypes, &this->bodyType);
MPI_Type_commit(&this->bodyType);
MPI_Get_address(&(box[0]), &baseAdress);
......@@ -57,6 +60,7 @@ namespace nbody {
MPI_Comm_rank(MPI_COMM_WORLD, &this->parallelRank);
this->domains = new Box[this->parallelSize];
this->correctState = true;
this->tree = new BarnesHutTree(this->parallelRank);
if (argc == 2) {
this->correctState = this->readInputData(string(argv[1]));
......@@ -80,6 +84,8 @@ namespace nbody {
}
MPI_Type_free(&this->bodyType);
MPI_Type_free(&this->boxType);
delete this->tree;
this->tree = NULL;
}
bool MpiSimulation::stateCorrect() {
......@@ -88,6 +94,7 @@ namespace nbody {
MpiSimulation::~MpiSimulation() {
delete[] this->domains;
delete this->tree;
}
int MpiSimulation::getNumberOfProcesses() {
......@@ -112,6 +119,7 @@ namespace nbody {
bb = nodes.front().getBB();
bb.extendForBodies(this->bodies);
nodes.front().setBB(bb);
bb.print(0);
while (nodes.size() < (unsigned int) this->parallelSize) {
int mostBodiesIndex = 0;
......@@ -143,6 +151,14 @@ namespace nbody {
this->comms.push_back(MpiBodyComm(&this->bodyType));
this->comms[0].recvBlocking(0, this->bodies);
}
Body::print(this->parallelRank, this->bodies);
/*
if (!Body::valid(this->bodies)) {
cout << "wrong distributed bodies" << endl;
} else {
cout << "bodies distributed ok" << endl;
}
*/
}
void MpiSimulation::distributeDomains(vector<Body> localBodies) {
......@@ -181,7 +197,14 @@ namespace nbody {
for (int i = 0; i < this->parallelSize; i++) {
if (i != this->parallelRank) {
this->domains[i].print(this->parallelRank);
vector<Body> refinements = this->tree.copyRefinements(this->domains[i]);
vector<Body> refinements = this->tree->copyRefinements(this->domains[i]);
/*
for (vector<Body>::iterator bit = refinements.begin(); bit != refinements.end(); bit++) {
bit->print();
}
*/
//cout << "ref " << refinements.size() << endl;
vector<MpiBodyComm>::iterator it = this->comms.begin();
......@@ -202,28 +225,28 @@ namespace nbody {
//integrate bodies in order of arrival to do communication/computation overlapping
this->comms[0].recvBlocking(MPI_ANY_SOURCE, refinements);
//cout << "recv: " << refinements.size() << endl;
this->tree.mergeLET(refinements);
this->tree->mergeLET(refinements);
//this->tree.getRootBB().print();
received++;
}
if (!this->tree.isCorrect()) {
if (!this->tree->isCorrect()) {
cout << "WRONG" << endl;
}
}
void MpiSimulation::buildTree() {
if (this->bodies.size() > 0 && this->overallDomain.isValid()) {
this->tree.build(this->bodies, this->overallDomain);
this->tree->build(this->bodies, this->overallDomain);
} else if (this->bodies.size() > 0) {
this->tree.build(this->bodies);
this->tree->build(this->bodies);
}
if (!this->tree.isCorrect()) {
if (!this->tree->isCorrect()) {
cout << "wrong tree" << endl;
}
}
void MpiSimulation::rebuildTree() {
this->tree.rebuild(this->overallDomain);
this->tree->rebuild(this->overallDomain);
//this->tree.getRootBB().print(this->parallelRank);
}
......@@ -235,7 +258,7 @@ namespace nbody {
this->distributeDomains(this->tree.advance());
this->rebuildTree();
*/
if (!this->tree.isCorrect()) {
if (!this->tree->isCorrect()) {
cout << "WRONG" << endl;
}/* else {
this->overallDomain.print(this->parallelRank);
......
......@@ -12,7 +12,7 @@ namespace nbody {
} else {
this->correctState = false;
}
this->tree.build(this->bodies);
this->tree->build(this->bodies);
}
PthreadSimulation::~PthreadSimulation() {
......@@ -45,9 +45,9 @@ namespace nbody {
//Box domain = this->tree.getRootBB();
//cout << this->tree.getRootBB().volume() << endl;
this->tree.rebuild();
this->tree.computeForces();
this->tree.advance();
this->tree->rebuild();
this->tree->computeForces();
this->tree->advance();
}
......
......@@ -8,6 +8,7 @@ namespace nbody {
this->parallelRank = -1;
this->parallelSize = -1;
this->correctState = 0;
this->tree = NULL;
}
Simulation::~Simulation() {
......
......@@ -14,7 +14,7 @@ namespace nbody {
int parallelRank;
int correctState;
vector<Body> bodies;
BarnesHutTree tree;
BarnesHutTree* tree;
public:
Simulation();
virtual ~Simulation();
......
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