Skip to content
Snippets Groups Projects
Commit 827b23e2 authored by Armin Damon Riess's avatar Armin Damon Riess
Browse files

change particles_ to positions_

parent a76e120c
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
lib/node.cpp
+ 20
− 22
View file @ 827b23e2
......@@ -18,9 +18,10 @@ Node::Node() {
}
*/
Node::Node(double* masses, double* particles, unsigned nParticles, unsigned depth,
Node::Node(double* masses, double* positions, unsigned nParticles, unsigned depth,
double size, Vector center, List localParticles, unsigned nLocalParticles) {
masses_ = masses;
positions_ = positions;
nParticles_ = nParticles;
depth_ = depth;
......@@ -44,9 +45,9 @@ Node::Node(double* masses, double* particles, unsigned nParticles, unsigned dept
} else {
isLeaf_ = false;
// array of Lists is empty -> segmentation fault ###########################################
// ##################### segmentation fault in line 52 #####################
// for each local particle, determine which octant it is in
List* particleDivision = new List[nLocalParticles];
List particleDivision[8] = {List{}};
for (unsigned i=0; i<nLocalParticles; ++i) {
particleDivision[getOctant(localParticles_[i])].push_back(i);
}
......@@ -61,11 +62,10 @@ Node::Node(double* masses, double* particles, unsigned nParticles, unsigned dept
for (unsigned j=0; j<3; ++j) {
childCenter[j] = center_[j] + (i & (1 << j) ? childSize : -childSize);
}
children_[i] = new Node(masses, particles, nParticles, depth+1,
children_[i] = new Node(masses, positions, nParticles, depth+1,
childSize, childCenter, particleDivision[i], particleDivision[i].size());
}
}
delete[] particleDivision;
}
}
......@@ -76,7 +76,7 @@ Node::~Node() {
}
Node& Node::operator=(const Node& other) {
particles_ = other.particles_;
positions_ = other.positions_;
nParticles_ = other.nParticles_;
depth_ = other.depth_;
size_ = other.size_;
......@@ -96,7 +96,7 @@ void Node::update(List& localParticles, unsigned nLocalParticles) {
// for each local particle, determine which octant it is in
// get the mass in this node while we're at it
mass_ = 0;
List* particleDivision = new List[nLocalParticles_];
List particleDivision[8] = {List{}};
for (unsigned i=0; i<nLocalParticles_; ++i) {
particleDivision[getOctant(localParticles_[i])].push_back(i);
mass_ += masses_[localParticles_[i]];
......@@ -122,23 +122,21 @@ void Node::update(List& localParticles, unsigned nLocalParticles) {
for (unsigned j=0; j<3; ++j) {
childCenter[j] = center_[j] + (i & (1 << j) ? childSize : -childSize);
}
children_[i] = new Node(masses_, particles_, nParticles_, depth_+1,
children_[i] = new Node(masses_, positions_, nParticles_, depth_+1,
childSize, childCenter, particleDivision[i], particleDivision[i].size());
isLeaf_ = false;
}
}
}
delete[] particleDivision;
}
Vector Node::calculateForce(unsigned particle) {
Vector force(3);
for (unsigned i=0; i<3; ++i) force[i] = 0;
// calculate theta
double dx = center_[0] - particles_[3*particle + 0];
double dy = center_[1] - particles_[3*particle + 1];
double dz = center_[2] - particles_[3*particle + 2];
double dx = center_[0] - positions_[3*particle + 0];
double dy = center_[1] - positions_[3*particle + 1];
double dz = center_[2] - positions_[3*particle + 2];
double lambda = std::sqrt(dx*dx + dy*dy + dz*dz);
double theta = size_/lambda;
......@@ -152,16 +150,16 @@ Vector Node::calculateForce(unsigned particle) {
for (unsigned i : localParticles_) {
if (i != particle) {
double dx = particles_[3*i + 0] - particles_[3*particle + 0];
double dy = particles_[3*i + 1] - particles_[3*particle + 1];
double dz = particles_[3*i + 2] - particles_[3*particle + 2];
double dx = positions_[3*i + 0] - positions_[3*particle + 0];
double dy = positions_[3*i + 1] - positions_[3*particle + 1];
double dz = positions_[3*i + 2] - positions_[3*particle + 2];
double r = std::sqrt(dx*dx + dy*dy + dz*dz);
double r3 = r*r*r;
double m = masses_[i];
for (unsigned j=0; j<3; ++j) {
force[j] += m * (particles_[3*i + j] - particles_[3*particle + j]) / r3;
force[j] += m * (positions_[3*i + j] - positions_[3*particle + j]) / r3;
}
}
}
......@@ -192,7 +190,7 @@ Vector Node::calculateForce(unsigned particle) {
double m = mass_;
for (unsigned j=0; j<3; ++j) {
force[j] += m * (centerOfMass_[j] - particles_[3*particle + j]) / r3;
force[j] += m * (centerOfMass_[j] - positions_[3*particle + j]) / r3;
}
}
......@@ -206,13 +204,13 @@ void Node::calculateCenterOfMass() {
for (unsigned i : localParticles_)
for (unsigned j=0; j<3; ++j)
centerOfMass_[j] += particles_[3*i + j] / nLocalParticles_;
centerOfMass_[j] += positions_[3*i + j] / nLocalParticles_;
}
unsigned Node::getOctant(unsigned particle) {
double px = particles_[localParticles_[0]];
double py = particles_[localParticles_[1]];
double pz = particles_[localParticles_[2]];
double px = positions_[3*particle + 0];
double py = positions_[3*particle + 1];
double pz = positions_[3*particle + 2];
if (px < center_[0]) {
if (py < center_[1]) {
if (pz < center_[2]) {
......
lib/node.hpp
+ 4
− 3
View file @ 827b23e2
......@@ -12,7 +12,7 @@ public:
Node() = default;
// Whoever creates the node is responsible for deciding which particles are in which octant and for allocating the localParticles array
// Node will write to localParticles array
Node(double* masses, double* particles, unsigned nParticles, unsigned depth,
Node(double* masses, double* positions, unsigned nParticles, unsigned depth,
double size, Vector center, List localParticles, unsigned nLocalParticles);
// Destuctor
~Node();
......@@ -26,7 +26,8 @@ public:
Vector calculateForce(unsigned particle);
// Calculate center of mass
void calculateCenterOfMass();
// Function that takes a particle and returns the octant it is in. This function doesn't check if the particle is actually in the node.
// Function that takes a particle index in the global particles_ array and returns the octant it is in relative to the center fo the current node
// this function doesn't check if the particle is actually in the node.
unsigned getOctant(unsigned particle);
private:
......@@ -35,7 +36,7 @@ private:
Node* children_[8];
double* masses_;
double* particles_;
double* positions_;
unsigned nParticles_;
unsigned depth_;
......
lib/tree.cpp
+ 3
− 3
View file @ 827b23e2
......@@ -4,8 +4,8 @@
#include <vector>
#include <iostream>
Tree::Tree(double* masses, double* particles, double* forces, unsigned nParticles, double size, Vector center) {
particles_ = particles;
Tree::Tree(double* masses, double* positions, double* forces, unsigned nParticles, double size, Vector center) {
positions_ = positions;
forces_ = forces;
nParticles_ = nParticles;
size_ = size;
......@@ -17,7 +17,7 @@ Tree::Tree(double* masses, double* particles, double* forces, unsigned nParticle
allParticles_ = List(nParticles_);
for (unsigned i=0; i < nParticles_; ++i) allParticles_[i] = i;
root_ = new Node(masses_, particles_, nParticles_, 0, size_, center_, allParticles_, nParticles_);
root_ = new Node(masses_, positions_, nParticles_, 0, size_, center_, allParticles_, nParticles_);
}
Tree::~Tree() {
......
lib/tree.hpp
+ 2
− 2
View file @ 827b23e2
......@@ -9,7 +9,7 @@ class Tree {
public:
// Constructor
Tree() = delete;
Tree(double* masses, double* particles, double* forces_, unsigned nParticles, double size, Vector center);
Tree(double* masses, double* positions, double* forces_, unsigned nParticles, double size, Vector center);
// Destructor
~Tree();
// Calculate force on a particle
......@@ -21,7 +21,7 @@ private:
Node* root_;
double* masses_;
double* particles_;
double* positions_;
double* forces_;
unsigned nParticles_;
......
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