Commit 155f54cf authored by gsaxena's avatar gsaxena

Parallelized the /examples/performance_test_voxels.cpp example - now every...

Parallelized the /examples/performance_test_voxels.cpp example - now every process updates its central voxel. A bit unfair for the parallel version
but the major work done is in comparison and not writing a single voxel substrate density. To let a single process (the appropriate process which conatins the
central voxel would be unfair as ALL processes must execute i.e. the ENTIRE DOMAIN must be covered).
Further, it is very difficult to change the signature of the supply function and target function pointers - this will need a change in the core files of BioFVM which would break the other
examples. This is also UNDESIRABLE as this is "just" a performance test and not a "real" example.
The output is produced in a scaling_*.txt file.
parent 6489af67
......@@ -114,6 +114,8 @@ main_experiment: ./examples/main_experiment.cpp $(BioFVM_OBJECTS) $(pugixml_OBJE
$(COMPILE_COMMAND) -o ./examples/tutorial3 ./examples/main_experiment.cpp $(BioFVM_OBJECTS) $(pugixml_OBJECTS)
3DtumorX: ./examples/3D_tumor_example_X.cpp $(BioFVM_OBJECTS) $(pugixml_OBJECTS)
$(COMPILE_COMMAND) -o ./examples/3DtumorX ./examples/3D_tumor_example_X.cpp $(BioFVM_OBJECTS) $(pugixml_OBJECTS)
perfVoxels: ./examples/performance_test_voxels.cpp $(BioFVM_OBJECTS) $(pugixml_OBJECTS)
$(COMPILE_COMMAND) -o ./examples/perfVoxels ./examples/performance_test_voxels.cpp $(BioFVM_OBJECTS) $(pugixml_OBJECTS)
clean:
rm -f *.o
......
......@@ -56,10 +56,15 @@
#include "../BioFVM.h"
/*----------------Include MPI Header-------------*/
#include <mpi.h>
/*----------------------------------------------*/
using namespace BioFVM;
int omp_num_threads = 8; // set number of threads for parallel computing
// set this to # of CPU cores x 2 (for hyperthreading)
/*------------------*/
/* Global variables */
/*------------------*/
int num_substrates;
int center_voxel_index=0;
......@@ -110,6 +115,59 @@ void process_output(double t, double dt, double mesh_resolution)
int main( int argc, char* argv[] )
{
unsigned int mpi_Error;
int mpi_Requested_level = MPI_THREAD_FUNNELED;
int mpi_Provided_level;
MPI_Comm mpi_Comm = MPI_COMM_WORLD;
MPI_Comm mpi_Cart_comm;
int mpi_Size, mpi_Rank;
int mpi_Dims[3], mpi_Is_periodic[3], mpi_Coords[3];
int mpi_Reorder;
mpi_Error = MPI_Init_thread(NULL, NULL, mpi_Requested_level, &mpi_Provided_level);
MPI_Barrier(mpi_Comm);
if(mpi_Error != MPI_SUCCESS)
MPI_Abort(mpi_Comm,-1);
MPI_Comm_size(mpi_Comm, &mpi_Size);
MPI_Comm_rank(mpi_Comm, &mpi_Rank);
if(mpi_Provided_level != mpi_Requested_level && mpi_Rank == 0)
{
std::cout<<"The MPI Requested Level is not available, please lower the level and try again"<<std::endl;
MPI_Abort(mpi_Comm, -1);
}
/*------------------------------------------------------------*/
/* MPI Comm Size and Rank */
/* Right now it should be 1-D X decomposition */
/*------------------------------------------------------------*/
mpi_Dims[0] = 1;
mpi_Dims[1] = mpi_Size; //Number of Y processes, since Y-processes divide X-axis, we have X-decomposition
mpi_Dims[2] = 1;
mpi_Is_periodic[0]=0;
mpi_Is_periodic[1]=0;
mpi_Is_periodic[2]=0;
mpi_Reorder = 0;
/*------------------------------------------------------------*/
/* Create MPI topology, find coords */
/*------------------------------------------------------------*/
MPI_Cart_create(mpi_Comm, 3, mpi_Dims, mpi_Is_periodic, mpi_Reorder, &mpi_Cart_comm);
MPI_Cart_coords(mpi_Cart_comm, mpi_Rank, 3, mpi_Coords);
/*----------------------*/
/* For time measurement */
/*----------------------*/
double t = 0.0;
double t_output_interval = 1.0;
double t_next_output_time = 0;
......@@ -118,24 +176,40 @@ int main( int argc, char* argv[] )
double dt = 0.01;
double t_max = 2.0;
double mesh_resolution=10.0;
/*--------------------------------------------------------------------*/
/* One input is expected i.e. the HALF length of the CUBIC domain */
/* Suppose we give 100 then domain in all directions is -100 to +100 */
/*--------------------------------------------------------------------*/
double domain_half_side= strtod(argv[1], NULL);
num_substrates=1;
/*------------------------------------------------------------------------------------*/
/* For the parallel version, we are going to update the central voxel of EACH process */
/*------------------------------------------------------------------------------------*/
std::vector<double> center(3);
center[0] = 0; center[1] = 0; center[2] = 0;
// openmp setup
omp_set_num_threads(omp_num_threads);
/*------------------------------------------------------------------*/
/* We control the threads in the SLURM submission script, hence no */
/* need for statement below UNLESS you want to reduce the number */
/* of threads from the maximum per process (i.e. tpp which is 24) */
/*------------------------------------------------------------------*/
//omp_set_num_threads(omp_num_threads);
// create a microenvironment;
Microenvironment microenvironment;
microenvironment.set_density(0, "substrate0" , "dimensionless" );
microenvironment.diffusion_coefficients[0] = 100000;
microenvironment.decay_rates[0] = 0.1;
double minX=-domain_half_side, minY=-domain_half_side, minZ=-domain_half_side, maxX=domain_half_side, maxY=domain_half_side, maxZ=domain_half_side;//, mesh_resolution=10;
microenvironment.resize_space_uniform( minX,maxX,minY,maxY,minZ,maxZ, mesh_resolution );
microenvironment.resize_space_uniform( minX,maxX,minY,maxY,minZ,maxZ, mesh_resolution, mpi_Dims, mpi_Coords );
// register the diffusion solver
microenvironment.diffusion_decay_solver = diffusion_decay_solver__constant_coefficients_LOD_3D;
......@@ -147,32 +221,50 @@ int main( int argc, char* argv[] )
{
microenvironment.density_vector(i)[0]= 100;
}
center_voxel_index= microenvironment.nearest_voxel_index(center);
// display information
microenvironment.display_information( std::cout );
center_voxel_index = (microenvironment.mesh.x_coordinates.size() * microenvironment.mesh.y_coordinates.size() * microenvironment.mesh.z_coordinates.size())/2;
// display information
if(mpi_Rank == 0)
microenvironment.display_information( std::cout );
//Timing measurements will start, hence a barrier
MPI_Barrier(MPI_COMM_WORLD);
BioFVM::RUNTIME_TIC();
BioFVM::TIC();
int output_index=0;
std::cout<<"number of agents: "<< all_basic_agents.size() <<std::endl;
while( t < t_max )
{
if( fabs( t - t_next_output_time ) < dt/10.0 )
{
process_output(t, dt, mesh_resolution);
if(mpi_Rank == 0)
process_output(t, dt, mesh_resolution);
t_next_output_time += t_output_interval;
next_output_index++;
}
// simulate microenvironment
microenvironment.simulate_bulk_sources_and_sinks( dt );
microenvironment.simulate_diffusion_decay( dt );
microenvironment.simulate_diffusion_decay( dt, mpi_Size, mpi_Rank, mpi_Coords, mpi_Dims, mpi_Cart_comm );
t += dt;
output_index++;
}
process_output(t_max, dt, mesh_resolution);
if(mpi_Rank == 0)
process_output(t_max, dt, mesh_resolution);
BioFVM::RUNTIME_TOC();
write_report(microenvironment.number_of_voxels(), BioFVM::runtime_stopwatch_value());
std::cout << "done!" << std::endl;
MPI_Barrier(MPI_COMM_WORLD);
if(mpi_Rank == 0)
write_report(microenvironment.number_of_voxels() * mpi_Size, BioFVM::runtime_stopwatch_value());
MPI_Finalize();
return 0;
}
\ No newline at end of file
......@@ -458,9 +458,9 @@ int main( int argc, char* argv[] )
MPI_Abort(mpi_Comm, -1);
}
/*------------------------------------------------------------*/
/* MPI Comm Size and Rank */
/* Right now it should be 1-D Z decomposition */
/*------------------------------------------------------------*/
/* MPI Comm Size and Rank */
/* Right now it should be 1-D X decomposition */
/*------------------------------------------------------------*/
......
#!/bin/bash
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=2
#SBATCH --cpus-per-task=24
#SBATCH -t 00:60:00
#SBATCH -o output-%j
#SBATCH -e error-%j
#SBATCH --exclusive
export OMP_DISPLAY_ENV=true
export OMP_NUM_THREADS=24
export OMP_PROC_BIND=spread
export OMP_PLACES=threads
mpiexec --map-by ppr:1:socket:pe=24 --report-bindings ./examples/perfVoxels 100
Markdown is supported
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