WIP: Openacc

Makefile

@@ -22,10 +22,12 @@ LDFLAGS                         = -lm $(EXTRA_LDFLAGS)
 RM                              = rm -rf
 MKDIR                           = mkdir -p
 
+# Final binary dependencies.
 COMMON_OBJECT_FILES             = $(BUILD_DIR)/config.o $(BUILD_DIR)/csv.o $(BUILD_DIR)/functions.o $(BUILD_DIR)/initialization.o $(BUILD_DIR)/main.o
 EXTRA_OBJECT_FILES              =
 OBJECT_FILES                    = $(COMMON_OBJECT_FILES) $(EXTRA_OBJECT_FILES)
 
+# main.cpp dependencies.
 COMMON_MAIN_DEPENDENCIES        = $(SRC_CXX_DIR)/main.cpp $(INC_DIR)/config.h $(INC_DIR)/csv.h $(INC_DIR)/data.h $(INC_DIR)/functions.h $(INC_DIR)/initialization.h
 EXTRA_MAIN_DEPENDENCIES         =
 MAIN_O_DEPENDENCIES             = $(COMMON_MAIN_DEPENDENCIES) $(EXTRA_MAIN_DEPENDENCIES)
@@ -45,6 +47,11 @@ ifdef PROFILING
 EXTRA_LDFLAGS                   = -lprofiler
 endif
 
+# If compiling with OpenACC directives
+ifdef OPENACC
+EXTRA_FLAGS                     += -fopenacc
+endif
+
 ###############################################################################
 # Main program compilation
 

src/c/functions.c

@@ -54,6 +54,7 @@ inline void compute_acceleration(const size_t particle_index,
 inline void apply_gravity(const size_t size,
                           const ParticleProperties *particles_properties,
                           Vector *forces) {
+#pragma acc parallel loop vector
   for (size_t i = 0; i < size; ++i) {
     forces[i].y_component -= (particles_properties[i].mass * 9.81d);
   }
@@ -112,7 +113,7 @@ inline void compute_forces(const double dt, const size_t particles_size,
                            const ParticleProperties *properties, const Contact *contacts,
                            const Vector *velocities, double *normal_forces,
                            double *tangent_forces, Vector *forces) {
-
+#pragma acc parallel loop vector
   for (size_t i = 0; i < contacts_size; ++i) {
     const size_t p1_idx = contacts[i].p1_idx;
     const size_t p2_idx = contacts[i].p2_idx;

src/cpp/main.cpp

@@ -58,6 +58,7 @@ void simulation_step(const size_t particles_size, const double dt) {
   compute_forces(dt, particles_size, contacts_size, particles, properties,
                  contacts_buffer, velocities, normal_forces, tangent_forces, forces);
 
+#pragma acc parallel loop vector
   for (size_t part = 0; part < particles_size; ++part) {
     compute_acceleration(part, properties, forces, accelerations);
     compute_velocity(dt, part, accelerations, velocities);
@@ -139,6 +140,11 @@ int main(int argc, char *argv[]) {
   }
 #endif
 
+#pragma acc data copyin(particles[0:num_particles], properties[0:num_particles], \
+                        contacts_buffer[0:size_triangular_matrix(num_particles)], \
+                        normal_forces[0:num_particles], tangent_forces[0:num_particles], \
+                        forces[0:num_particles], accelerations[0:num_particles], \
+                        velocities[0:num_particles], displacements[0:num_particles])
   for (unsigned long step = 1; step <= max_steps; ++step) {
 #ifdef DEBUG_STEP
     current_step = step;