13.5 The Metropolis Algorithm and the Two-dimensional Ising Model 441
#include
#include
#include
#include
#include"lib.h"
using namespacestd;
// output file
ofstream ofile;
// inline function for periodic boundary conditions
inline intperiodic(inti,intlimit,intadd){
return(i+limit+add) % (limit);
}
// Function to initialise energy and magnetization
voidinitialize(int,int,double&,double&);
// The metropolis algorithm
voidMetropolis(int,long&,int,double&,double&,double);
// prints to file the results of the calculations
voidoutput(int,int,double,double);
// Main program begins here
intmain(intargc,charargv[])
{
charoutfilename;
longidum;
int*spin_matrix, n_spins, mcs, my_rank, numprocs;
doublew[17], average[5], total_average[5],
initial_temp, final_temp, E, M, temp_step;
// MPI initializations
MPI_Init (&argc, &argv);
MPI_Comm_size (MPI_COMM_WORLD, &numprocs);
MPI_Comm_rank (MPI_COMM_WORLD, &my_rank);
if(my_rank == 0 && argc <= 1){
cout <<"Bad Usage: "<< argv[0] <<
" read output file"<< endl;
exit(1);
}
if(my_rank == 0 && argc > 1){
outfilename=argv[1];
ofile.open(outfilename);
}
n_spins = 40; mcs = 1000000; initial_temp = 2.4; final_temp = 2.7; temp_step =0.1;
/
Determine number of intervall which are used by all processes
myloop_begin gives the starting point on process my_rank
myloop_end gives the end point for summation on process my_rank
/
intno_intervalls = mcs/numprocs;
intmyloop_begin = my_rankno_intervalls + 1;
intmyloop_end = (my_rank+1)*no_intervalls;
if( (my_rank == numprocs-1) &&( myloop_end < mcs) ) myloop_end = mcs;
// broadcast to all nodes common variables
MPI_Bcast (&n_spins, 1, MPI_INT, 0, MPI_COMM_WORLD);
MPI_Bcast (&initial_temp, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast (&final_temp, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
MPI_Bcast (&temp_step, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);