void wham_FEP()

in Documents/NamdSample/namd210/lib/replica/FEP_wca/wham/wham-fep.cpp [74:165]


void wham_FEP(void)
{
	int i, Time_1, Time_2, Time_3;
	double Temp, lambx, eppx;
	int Ntot1, Ntot2, ReadItem;
	int ioffset=0;	// the window used as the reference
	FILE *fIn, *fOut;
	
	Ntot1 = Ntot2 = 0;
	Temp = 300.0;	
	kbt=KBOLTZ*Temp;
	
	nb_wind = 0;
	NPoints[0] = 0;
	lambda[0] = 0.0;
	
	fIn = fopen(szData, "r");
	if(fIn == NULL)	{
		printf("Fail to open file: %s\nQuit\n", szData);
		Free_Memory();
		exit(1);
	}

	while(1)	{
		ReadItem = fscanf(fIn, "%lf%lf", &lambx, &eppx);
		if(ReadItem != 2)	{
			break;
		}
		
		if( lambx != lambda[nb_wind] )	{
			nb_wind++;
			if(nb_wind > max_wind)	{
				printf("Number of windows exceeds limit.\nQuit\n");
				exit(1);
			}
			NPoints[nb_wind] = 0;
			if(Ntot1 > Ntot2) Ntot2=Ntot1;
			Ntot1 = 0;
			lambda[nb_wind] = lambx;
		}
		
		
		Ntot1 ++;
		
		if(NPoints[nb_wind] < max_data)	{
			epprtd[nb_wind][NPoints[nb_wind]] = eppx;
			NPoints[nb_wind] = NPoints[nb_wind] + 1;
		}
	}
	fclose(fIn);
	
	if(Ntot2 > max_data)	{
		printf("Some data points have been left out.\n");
	}
	else	{
		printf("All data points have been read.\n");
	}
	
	printf("Window       #_of_Data       Lambda\n");
	
	for(i=0; i<=nb_wind; i++)	{
		printf("%6d       %10d      %.5lf\n", i+1, NPoints[i], lambda[i]);
        F[i] = 0.0;
        F2[i] = 0.0;
	}
	
//	Time_1 = time(NULL);
//	wham();		// the original wham by Benoit Roux
//	Time_2 = time(NULL);
//	printf("Ellapsed time %d seconds in wham().\n", Time_2-Time_1);

//	printf("\n\nWindow       #_of_Data     Lambda    F\n");
//	for(i=0; i<=nb_wind; i++)	{
//		printf("    %-d   %10d        %.5lf  %.5lf\n", i+1, NPoints[i], lambda[i], F2[i]-F2[ioffset]);
//	}

	Time_2 = time(NULL);
	wham_E_Histogram();	// loop over the histogram
	Time_3 = time(NULL);
	printf("\nEllapsed time %d seconds in wham().\n", Time_3-Time_2);

	printf("\n\nWindow       #_of_Data     Lambda    F\n");
	for(i=0; i<=nb_wind; i++)	{
		printf(" %4d   %10d        %.5lf  %.5lf\n", i+1, NPoints[i], lambda[i], F2[i]-F2[ioffset]);
	}

	fOut = fopen(szFile_FE, "w");
	fprintf(fOut, " 100  %.5lf\n", F2[nb_wind]-F2[0]);
	fclose(fOut);
	

}