better printE

energy.c

@@ -181,6 +181,13 @@ double targetenergy(Chain *chain)
 	//return chain->Erg(0, 0);
 }
 
+double cyclicenergy(Chain *chain)
+{
+
+	return chain->Erg(1, chain->NAA - 1);
+	//return chain->Erg(0, 0);
+}
+
 
 /* Print the energy matrix of a chain */
 void energy_matrix_print(Chain *chain, Biasmap *biasmap, model_params *mod_params) {

energy.h

@@ -20,6 +20,7 @@ double totenergy(Chain *chain);
 double locenergy(Chain *chain);
 double extenergy(Chain *chain);
 double targetenergy(Chain *chain);
+double cyclicenergy(Chain *chain);
 int bestRot(Chain *chain);
 void energy_matrix_print(Chain *,Biasmap *, model_params *mod_params);
 void biasmap_initialise(Chain *,Biasmap *, model_params *mod_params);

main.c

@@ -335,7 +335,6 @@ void simulate(Chain * chain, Chaint *chaint, Biasmap* biasmap, simulation_params
 				fprintf(stderr, "swap out no improv curr %g swap %g best %g \n", currTargetEnergy, swapEnergy[swapInd], targetBest);
 			}
 			else if (currTargetEnergy - targetBest <= 5.0) {
-
 				for (ind = 0; ind < swapLength; ++ind) {
 					if (swapEnergy[ind] == 9999.) swapInd = ind;
 					if (calculateRMSD(swapChains[ind], chain) < 2.0) {
@@ -378,7 +377,7 @@ void simulate(Chain * chain, Chaint *chaint, Biasmap* biasmap, simulation_params
 					//sim_params->protein_model.external_k[0] = external_k;
 					//bestIndex = i*sim_params->pace + j;
 				}
-				else if (currIndex - lastIndex > 200000 && rand() % 1000000 < 1) {
+				else if (currIndex - lastIndex > 200000 && rand() % 100 < 5) {
 					if (1 || external_k == sim_params->protein_model.external_k[0]) {
 
 						swapInd = rand() % (swapLength + 1);
@@ -390,6 +389,9 @@ void simulate(Chain * chain, Chaint *chaint, Biasmap* biasmap, simulation_params
 						//sim_params->protein_model.external_k[0] = external_k;
 					}
 				}
+				lastIndex = currIndex;
+
+
 			}
 			else if (currIndex - lastIndex > 200000 && rand() % 100 < 5) {
 				if (1 || external_k == sim_params->protein_model.external_k[0]) {
@@ -713,14 +715,14 @@ void simulate(Chain * chain, Chaint *chaint, Biasmap* biasmap, simulation_params
 			if (sim_params->thermobeta != sim_params->beta1)
 				continue;
 #endif
-			fprintf(stderr, "curr cyclic energy  %g \n", extenergy(chain));
+			if (sim_params->protein_model.external_potential_type2 == 4) fprintf(stderr, "curr cyclic energy  %g \n", extenergy(chain));
 			fprintf(sim_params->outfile, "-+- TEST BLOCK %5d -+-\n", i);
 			tests(chain, biasmap, sim_params->tmask, sim_params, 0x11, NULL);
 		}
 	}
 
 	freemem_chain(chain2); free(chain2);
-
+	energy_matrix_print(chain, biasmap, sim_params);
 }
 
 char *read_options(int argc, char *argv[], simulation_params *sim_params)

metropolis.c

@@ -79,16 +79,12 @@ static int allowed(Chain *chain, Chaint *chaint, Biasmap* biasmap, int start, in
 	double internalloss = loss;
 	double external_k = 1.0;
 	if (sim_params->protein_model.external_potential_type == 5)	external_k = sim_params->protein_model.external_k[0];
-	if (q > 0) external_k = 0.02;
+	if (q > 5) external_k = 0.02;
 	//loss is negative!! if loss is negative, it's worse, bad
 	/* Metropolis criteria */
 	//loss += q - chain->Erg(0, 0);
 	//loss = loss/sqrt(chain->NAA) + externalloss;
 	loss = loss + externalloss;
-	if (external_k == 0.5) {
-		fprintf(stderr, "heat starts \n");
-		sim_params->protein_model.external_k[0] = 0.51;
-	}
 
 
 
@@ -472,8 +468,12 @@ static int crankshaft(Chain * chain, Chaint *chaint, Biasmap *biasmap, double am
 	}
 
 	//for cyclic peptide Gary Hack
-	while (sim_params->protein_model.external_potential_type2 == 4 && (start <= 1 || end >= sim_params->NAA-1) && ((toss >>2 )%1000) > chain->Erg(0, 0)) {
+	while (sim_params->protein_model.external_potential_type2 == 4 && (start <= 1 || end >= sim_params->NAA-1) && ((toss >>2 )%10) > chain->Erg(0, 0)) {
 		toss = rand();
+		len = toss & 0x3;	/* segment length minus one */
+		if (len > chain->NAA - 2)
+			len = chain->NAA - 2;
+		N_len = sim_params->MC_lookup_table_n[len];
 		start = sim_params->MC_lookup_table[len*N + (toss >> 2) % N_len];
 		if ((sim_params->protein_model).fix_CA_atoms) {
 			end = start + 1;

probe.c

@@ -1917,8 +1917,8 @@ void cm_ideal_4(Chain *chain, Biasmap *biasmap,simulation_params *sim_params, vo
 void ergtot(Chain *chain,Biasmap *biasmap, simulation_params *sim_params, void *mpi_comm)
 {	/*changed, can be removed */ 
 	double tote = totenergy(chain);
-	fprintf(sim_params->outfile,"Energy = %.6f ( %g %g ) Rotmers (%i)\n", tote, locenergy(chain),extenergy(chain),bestRot(chain));
-	fprintf(stderr,"%.6f\n ", tote);
+	fprintf(sim_params->outfile, "Energy = totalE %.6f ( diagnolE %.6f extE %.6f cycE %.6f) Rotmers (%i)\n", tote, locenergy(chain), extenergy(chain), cyclicenergy(chain), bestRot(chain));
+	fprintf(stderr, "totalE %.6f extE %.6f \n", tote, extenergy(chain));
 }
 
 void vdw_max_gamma(Chain *chain, Biasmap *biasmap,simulation_params *sim_params, void *mpi_comm) {