ctest for remd

checks average ACC, plots histograms of energy for each bath
temperature and energy v.s bath temperature

diff --git a/ctest/1L2Y_remd.inp b/ctest/1L2Y_remd.inp
new file mode 100644 (file)
index 0000000..35d6493
--- /dev/null
+++ b/ctest/1L2Y_remd.inp
@@ -0,0 +1,17 @@
+1L2Y small remd test
+SEED=-6086648 RE pdbref refstr ONE_LETTER  rescale_mode=2 out1file extconf
+nstep=200000 ntwe=100  ntwx=1000 lang=0 tbf t_bath=300 tau_bath=1.0            &
+dt=0.10 reset_moment=1000 reset_vel=1000 damax=4.0 restart
+nrep=8 nstex=1000 tlist mlist traj1file rest1file sync nsyn=1000 
+ 240 260 280 300 320 340 360 390 
+   1   1   1   1   1   1   1   1 
+WSC=1.00000 WSCP=1.23315 WELEC=0.84476 WBOND=1.00000 WANG=0.62954              &
+WSCLOC=0.10554 WTOR=1.84316 WTORD=1.26571 WCORRH=0.19212 WCORR5=0.00000        &
+WCORR6=0.00000 WEL_LOC=0.37357 WTURN3=1.40323 WTURN4=0.64673 WTURN6=0.00000    &
+WVDWPP=0.23173 WHPB=1.00000 WSCCOR=0.25                                        &
+CUTOFF=7.00000 WCORR4=0.00000
+1L2Y.pdb
+22
+XNLYIQWLKDGGPSSGRPPPSX
+0
+0

diff --git a/ctest/matplotlib_hist.py b/ctest/matplotlib_hist.py
new file mode 100755 (executable)
index 0000000..86881a6
--- /dev/null
+++ b/ctest/matplotlib_hist.py
@@ -0,0 +1,43 @@
+#! /usr/bin/env python
+
+import matplotlib
+#matplotlib.use('GTK')
+matplotlib.use('Agg')
+import matplotlib.pyplot as plt
+import matplotlib.cm as cm
+import numpy as np
+
+x,y= np.loadtxt('remd_all.stat',usecols=(11,3),unpack=True)
+
+h,bin=np.histogram(y,bins=40,density=False)
+
+plt.xlim(min(bin), max(bin))
+plt.ylim(0,max(h)/4)
+plt.ylabel('number of samples')
+plt.xlabel('potential energy [kcal/mol]')
+
+Tremd=[240, 260, 280, 300, 320, 340, 360, 390]
+colors = cm.rainbow(np.linspace(0, 1, len(Tremd)))
+for T,c in zip(Tremd,colors):
+ yt=y[x==T]
+ h,bin=np.histogram(yt,bins=40,density=False)
+ center = (bin[:-1] + bin[1:]) / 2
+ plt.plot(center,h,'-',color=c)
+# plt.bar(bin[:-1], h, width = bin[2]-bin[1],color=c)
+
+plt.savefig('1L2Y_remd_ene_hist.png')
+#plt.show()   
+
+plt.clf()
+plt.xlabel('bath temperature [K]')
+plt.ylabel('potential energy [kcal/mol]')
+
+plt.xlim(230, 400)
+Tremd=[240, 260, 280, 300, 320, 340, 360, 390]
+colors = cm.rainbow(np.linspace(0, 1, len(Tremd)))
+for T,c in zip(Tremd,colors):
+ yt=y[x==T]
+ xt=x[x==T]
+ plt.plot(xt,yt,'.',color=c)
+    
+plt.savefig('1L2Y_remd_Tene.png')

diff --git a/ctest/prota_unres_energy_check.sh b/ctest/prota_unres_energy_check.sh
index 3aa859d..868163f 100755 (executable)
--- a/ctest/prota_unres_energy_check.sh
+++ b/ctest/prota_unres_energy_check.sh
@@ -93,6 +93,23 @@ elif [ "$1" == "1L2Y_B" ]; then
     exit 0
  fi
 
+
+elif [ "$1" == "1L2Y_remd" ]; then
+ tail -q -n +5 1L2Y_remd*.stat >remd_all.stat
+ ./matplotlib_hist.py
+
+ echo  "<DartMeasurementFile name=\"Histograms $1\" type=\"image/png\">`pwd`/1L2Y_remd_ene_hist.png</DartMeasurementFile>"
+ echo  "<DartMeasurementFile name=\"Energy $1\" type=\"image/png\">`pwd`/1L2Y_remd_Tene.png</DartMeasurementFile>" 
+ grep ACC 1L2Y_remd.out_GB000 |tail -7
+ exchange=`grep ACC 1L2Y_remd.out_GB000 |tail -7|awk '{a=a+$4}END{print a/NR}'`
+ echo "average exchange = ${exchange}"
+ if [ `echo "a=${exchange};a<0.1"|bc -l` != "0" ]; then
+    echo 'ERROR average exchange smaller than 0.1'
+    exit 1
+ else
+    exit 0
+ fi
+               
 else
  exit 1
 fi

diff --git a/source/unres/src_MD/CMakeLists.txt b/source/unres/src_MD/CMakeLists.txt
index c90d228..1c07c18 100644 (file)
--- a/source/unres/src_MD/CMakeLists.txt
+++ b/source/unres/src_MD/CMakeLists.txt
@@ -413,6 +413,9 @@ FILE(COPY ${CMAKE_SOURCE_DIR}/ctest/1L2Y_NH.inp
 FILE(COPY ${CMAKE_SOURCE_DIR}/ctest/1L2Y.pdb
         DESTINATION ${CMAKE_CURRENT_BINARY_DIR} )
 
+FILE(COPY ${CMAKE_SOURCE_DIR}/ctest/1L2Y_remd.inp
+        DESTINATION ${CMAKE_CURRENT_BINARY_DIR} )
+
 
 FILE(COPY ${CMAKE_SOURCE_DIR}/ctest/prota_unres_energy_check.sh
         DESTINATION ${CMAKE_CURRENT_BINARY_DIR} 
@@ -424,6 +427,12 @@ FILE(COPY ${CMAKE_SOURCE_DIR}/ctest/matplotlib_fit_hist.py
         FILE_PERMISSIONS OWNER_READ OWNER_WRITE OWNER_EXECUTE GROUP_READ GROUP_EXECUTE WORLD_READ WORLD_EXECUTE
 )
 
+FILE(COPY ${CMAKE_SOURCE_DIR}/ctest/matplotlib_hist.py
+        DESTINATION ${CMAKE_CURRENT_BINARY_DIR} 
+        FILE_PERMISSIONS OWNER_READ OWNER_WRITE OWNER_EXECUTE GROUP_READ GROUP_EXECUTE WORLD_READ WORLD_EXECUTE
+)
+
+
 
 #=========================================
 #  test_prota_E0LL2Y.sh
@@ -464,6 +473,44 @@ FILE(      COPY ${CMAKE_CURRENT_BINARY_DIR}/scripts/test_prota_E0LL2Y.sh
        DESTINATION ${CMAKE_CURRENT_BINARY_DIR}
        FILE_PERMISSIONS OWNER_READ OWNER_WRITE OWNER_EXECUTE GROUP_READ GROUP_EXECUTE WORLD_READ WORLD_EXECUTE
 )
+
+FILE(WRITE ${CMAKE_CURRENT_BINARY_DIR}/scripts/test_mpi_E0LL2Y.sh
+"#!/bin/sh
+export POT=GB
+export FGPROCS=$2
+export PREFIX=$1
+#-----------------------------------------------------------------------------
+UNRES_BIN=${CMAKE_BINARY_DIR}/bin/${UNRES_BIN}
+#-----------------------------------------------------------------------------
+DD=${CMAKE_SOURCE_DIR}/PARAM
+export BONDPAR=$DD/bond_AM1.parm
+export THETPAR=$DD/theta_abinitio.parm
+export ROTPAR=$DD/rotamers_AM1_aura.10022007.parm
+export TORPAR=$DD/torsion_631Gdp.parm
+export TORDPAR=$DD/torsion_double_631Gdp.parm
+export ELEPAR=$DD/electr_631Gdp.parm
+export SIDEPAR=$DD/scinter_$POT.parm
+export FOURIER=$DD/fourier_opt.parm.1igd_hc_iter3_3
+export SCPPAR=$DD/scp.parm
+export SCCORPAR=$DD/sccor_am1_pawel.dat
+export THETPARPDB=$DD/thetaml.5parm
+export ROTPARPDB=$DD/scgauss.parm
+export PATTERN=$DD/patterns.cart
+#-----------------------------------------------------------------------------
+echo CTEST_FULL_OUTPUT
+mpiexec -np $3 $UNRES_BIN | grep -v traj1file
+./prota_unres_energy_check.sh $1 ${UNRES_BIN}
+")
+
+#
+# File permissions workaround
+#
+FILE(  COPY ${CMAKE_CURRENT_BINARY_DIR}/scripts/test_mpi_E0LL2Y.sh 
+       DESTINATION ${CMAKE_CURRENT_BINARY_DIR}
+       FILE_PERMISSIONS OWNER_READ OWNER_WRITE OWNER_EXECUTE GROUP_READ GROUP_EXECUTE WORLD_READ WORLD_EXECUTE
+)
+
+
 endif(UNRES_MD_FF STREQUAL "E0LL2Y")
 
 
@@ -510,6 +557,7 @@ else(NOT UNRES_WITH_MPI)
     add_test(NAME UNRES_Langevin COMMAND mpiexec ${boot_lam} -np 2 ${CMAKE_CURRENT_BINARY_DIR}/test_prota_E0LL2Y.sh 1L2Y_L 2 )
     add_test(NAME UNRES_NoseHoover COMMAND mpiexec ${boot_lam} -np 2 ${CMAKE_CURRENT_BINARY_DIR}/test_prota_E0LL2Y.sh 1L2Y_NH 2 )
     add_test(NAME UNRES_Berendsen COMMAND mpiexec ${boot_lam} -np 2 ${CMAKE_CURRENT_BINARY_DIR}/test_prota_E0LL2Y.sh 1L2Y_B 2 )
+    add_test(NAME UNRES_remd COMMAND sh ${CMAKE_CURRENT_BINARY_DIR}/test_mpi_E0LL2Y.sh 1L2Y_remd 1 8 )
   endif(UNRES_MD_FF STREQUAL "E0LL2Y")
 
 endif(NOT UNRES_WITH_MPI)