
C
C  Take measurements of energy and magnetization.
C
	SUBROUTINE init_measure
	INCLUDE "common.h"

C ----  Initialize sums and sum of squares.
	total_energy = 0.
	total_enersq = 0.
	total_mag   = 0.
	total_magsq = 0.
	total_accept = 0.
	RETURN 
	END


C----------------------------------------------------------------------------


C
C  Take measurements of energy and magnetization.
C
	SUBROUTINE measure
	INCLUDE "common.h"
	INTEGER i,j, ener,magn 
	REAL energy,magnetization
        INTEGER ctemp(3),type,bytes
	INTEGER intadd, d1
	EXTERNAL intadd

	ener = 0
	magn = 0
	DO i=1,isiz
	DO j=1,jsiz
	  ener = ener + spin(i,j) * ( spin(i+1,j) + spin(i,j+1) )
	  magn = magn + spin(i,j)
	ENDDO
	ENDDO
 
C ----  Sum the energy and magnetization over processors.
	type = COMB
	ctemp(1) = ener
	ctemp(2) = magn
	ctemp(3) = accept
	bytes = 4
        d1 = KXCOMB(ctemp,intadd,bytes,3,IALNOD,0,type)
	ener = ctemp(1)
	magn = ctemp(2)
	accept = ctemp(3)

C ----  Normalize.
	IF (magn.LT.0) THEN
	  magn = - magn
	ENDIF
	energy = - FLOAT(ener) / volume
	magnetization = FLOAT(magn) / volume
 
C ----  Sums and sum of squares.
	total_energy = total_energy + energy
	total_enersq = total_enersq + (energy * energy)
	total_mag    = total_mag    + magnetization
	total_magsq  = total_magsq  + (magnetization * magnetization)
	total_accept = total_accept + FLOAT(accept) / volume

C       WRITE(6,*) 'Energy: ', energy
C       WRITE(6,*) 'Magnet: ', magnetization
C       WRITE(6,*) 'Accept: ', FLOAT(accept) / volume
C	CALL KFLUSH(6)

	RETURN 
	END

	
C----------------------------------------------------------------------------

	
C
C  Print out averages and errors for energy and magnetization, and
C  the average acceptance for the Metropolis update.
C
	SUBROUTINE output
	INCLUDE "common.h"
	REAL err_ener, err_mag
	REAL*8 specific_heat, susceptibility
	
C ----  Normalize.
	total_energy = total_energy / FLOAT(measurements)
	total_enersq = total_enersq / FLOAT(measurements)
	total_mag    = total_mag    / FLOAT(measurements)
	total_magsq  = total_magsq  / FLOAT(measurements)
	total_accept = total_accept / FLOAT(measurements)

C ----  Variances and errors.
	specific_heat = total_enersq - (total_energy * total_energy)
	susceptibility = total_magsq - (total_mag * total_mag)
	err_ener = SQRT( specific_heat  / measurements )
	err_mag  = SQRT( susceptibility / measurements )
	specific_heat  = specific_heat  * volume
	susceptibility = susceptibility * volume

C ----  Print out results.
	WRITE(6,*) ' '
	WRITE(6,10) 'Average energy  = ',total_energy,'  +/- ',err_ener
	WRITE(6,10) 'Average magnet  = ',total_mag,'  +/- ',err_mag
	WRITE(6,*) ' '
	WRITE(6,20) 'Specific Heat   = ',specific_heat
	WRITE(6,20) 'Susceptibility  = ',susceptibility
	WRITE(6,*) ' '
	WRITE(6,20) 'Average accept  = ',total_accept
  10    FORMAT(A,F10.6,A,F10.6)
  20    FORMAT(A,F10.6)

	RETURN 
	END

	
C----------------------------------------------------------------------------


C
C  Integer add over processors.
C
	INTEGER FUNCTION intadd(i1,i2,size)
	INTEGER i1, i2, size
	i1 = i1 + i2
	intadd = 0
	RETURN 
	END

