00001 00007 // written by David Hurwitz, March to May 1998. 00008 00009 #include <memory.h> 00010 //#include <iomanip.h> 00011 #include <stdio.h> 00012 #include "InfoStream.h" 00013 #include "FreeEnergyEnums.h" 00014 #include "FreeEnergyAssert.h" 00015 #include "FreeEnergyLambda.h" 00016 #include "FreeEnergyLambdMgr.h" 00017 #include "FreeEnergyGroup.h" 00018 00019 ALambdaManager::ALambdaManager() { 00020 //------------------------------------------------------------------------ 00021 // make room for some LambdaControl objects. 00022 //------------------------------------------------------------------------ 00023 m_ActiveIndex = 0; 00024 m_NumObjects = 0; 00025 m_pPmfBlocks = new ALambdaControl[kLambdaNumToStart]; 00026 m_MaxNum = kLambdaNumToStart; 00027 } 00028 00029 00030 ALambdaManager::~ALambdaManager() { 00031 //------------------------------------------------------------------------ 00032 // return borrowed memory to the free store. 00033 //------------------------------------------------------------------------ 00034 delete []m_pPmfBlocks; 00035 } 00036 00037 00038 int ALambdaManager::GetNum_dU_dLambda() { 00039 //------------------------------------------------------------------------ 00040 // get number times dU_dLambda accumulated from active lambda object. 00041 //------------------------------------------------------------------------ 00042 ASSERT((*this)[m_ActiveIndex].IsActive()); 00043 return((*this)[m_ActiveIndex].GetNum_dU_dLambda()); 00044 } 00045 00046 00047 void ALambdaManager::Integrate_MCTI() { 00048 //------------------------------------------------------------------------ 00049 // integrate <dU/dLambda> for MCTI 00050 //------------------------------------------------------------------------ 00051 ASSERT((*this)[m_ActiveIndex].IsActive()); 00052 (*this)[m_ActiveIndex].Integrate_MCTI(); 00053 } 00054 00055 00056 void ALambdaManager::Accumulate(double dU_dLambda) { 00057 //------------------------------------------------------------------------ 00058 // add dU_dLambda to the active accumulator (in one of the lambda objects) 00059 //------------------------------------------------------------------------ 00060 ASSERT((*this)[m_ActiveIndex].IsActive()); 00061 (*this)[m_ActiveIndex].Accumulate(dU_dLambda); 00062 } 00063 00064 00065 double ALambdaManager::GetAccumulation() { 00066 //------------------------------------------------------------------------ 00067 // get the accumulation of dU_dLambda from active lambda object. 00068 //------------------------------------------------------------------------ 00069 ASSERT((*this)[m_ActiveIndex].IsActive()); 00070 return((*this)[m_ActiveIndex].GetAccumulation()); 00071 } 00072 00073 00074 double ALambdaManager::GetIntegration() { 00075 //------------------------------------------------------------------------ 00076 // get accumulation of <dU_dLambda> * dLambda from active lambda object. 00077 //------------------------------------------------------------------------ 00078 ASSERT((*this)[m_ActiveIndex].IsActive()); 00079 return((*this)[m_ActiveIndex].GetIntegration()); 00080 } 00081 00082 00083 void ALambdaManager::ZeroAccumulator() { 00084 //------------------------------------------------------------------------ 00085 // zero accumulation of dU_dLambda in the active lambda object. 00086 //------------------------------------------------------------------------ 00087 ASSERT((*this)[m_ActiveIndex].IsActive()); 00088 (*this)[m_ActiveIndex].ZeroAccumulator(); 00089 } 00090 00091 00092 Bool_t ALambdaManager::IsFirstStep() { 00093 //------------------------------------------------------------------------ 00094 // ASSUMING that the m_ActiveIndex'th LambdaControl is active, 00095 // decide if it's time to print restraint information 00096 //------------------------------------------------------------------------ 00097 ASSERT((*this)[m_ActiveIndex].IsActive()); 00098 return((*this)[m_ActiveIndex].IsFirstStep()); 00099 } 00100 00101 00102 Bool_t ALambdaManager::IsTimeToPrint() { 00103 //------------------------------------------------------------------------ 00104 // ASSUMING that the m_ActiveIndex'th LambdaControl is active, 00105 // decide if it's time to print restraint information 00106 //------------------------------------------------------------------------ 00107 ASSERT((*this)[m_ActiveIndex].IsActive()); 00108 return((*this)[m_ActiveIndex].IsTimeToPrint()); 00109 } 00110 00111 00112 Bool_t ALambdaManager::IsTimeToPrint_dU_dLambda() { 00113 //------------------------------------------------------------------------ 00114 // ASSUMING that the m_ActiveIndex'th LambdaControl is active, 00115 // decide if it's time to print du/dLambda information 00116 //------------------------------------------------------------------------ 00117 ASSERT((*this)[m_ActiveIndex].IsActive()); 00118 return((*this)[m_ActiveIndex].IsTimeToPrint_dU_dLambda()); 00119 } 00120 00121 00122 Bool_t ALambdaManager::IsTimeToClearAccumulator() { 00123 //------------------------------------------------------------------------ 00124 // ASSUMING that the m_ActiveIndex'th LambdaControl is active, 00125 // decide if it's time to start accumulating dU/dLambda from zero 00126 //------------------------------------------------------------------------ 00127 ASSERT((*this)[m_ActiveIndex].IsActive()); 00128 return((*this)[m_ActiveIndex].IsTimeToClearAccumulator()); 00129 } 00130 00131 00132 Bool_t ALambdaManager::IsEndOf_MCTI_Step() { 00133 //------------------------------------------------------------------------ 00134 // ASSUMING that the m_ActiveIndex'th LambdaControl is active, 00135 // decide if this is the last time step of an MCTI step 00136 //------------------------------------------------------------------------ 00137 ASSERT((*this)[m_ActiveIndex].IsActive()); 00138 return((*this)[m_ActiveIndex].IsEndOf_MCTI_Step()); 00139 } 00140 00141 00142 Bool_t ALambdaManager::IsEndOf_MCTI() { 00143 //------------------------------------------------------------------------ 00144 // ASSUMING that the m_ActiveIndex'th LambdaControl is active, 00145 // decide if this is the last time step of an MCTI block 00146 //------------------------------------------------------------------------ 00147 ASSERT((*this)[m_ActiveIndex].IsActive()); 00148 return((*this)[m_ActiveIndex].IsEndOf_MCTI()); 00149 } 00150 00151 00152 void ALambdaManager::PrintLambdaHeader(double dT) { 00153 //------------------------------------------------------------------------ 00154 // print header for a new lambda control object 00155 //------------------------------------------------------------------------ 00156 ASSERT((*this)[m_ActiveIndex].IsActive()); 00157 (*this)[m_ActiveIndex].PrintLambdaHeader(dT); 00158 } 00159 00160 00161 void ALambdaManager::PrintHeader(double dT) { 00162 //------------------------------------------------------------------------ 00163 // print information about current time step 00164 //------------------------------------------------------------------------ 00165 ASSERT((*this)[m_ActiveIndex].IsActive()); 00166 (*this)[m_ActiveIndex].PrintHeader(dT); 00167 } 00168 00169 00170 int ALambdaManager::GetNumStepsSoFar() { 00171 //------------------------------------------------------------------------ 00172 // return the number of steps taken in the active LambdaControl block 00173 //------------------------------------------------------------------------ 00174 ASSERT((*this)[m_ActiveIndex].IsActive()); 00175 return((*this)[m_ActiveIndex].GetNumStepsSoFar()); 00176 } 00177 00178 00179 int ALambdaManager::GetNumAccumStepsSoFar() { 00180 //------------------------------------------------------------------------ 00181 // return the total number of steps dU/dLambda has been accumulated 00182 //------------------------------------------------------------------------ 00183 ASSERT((*this)[m_ActiveIndex].IsActive()); 00184 return((*this)[m_ActiveIndex].GetNumAccumStepsSoFar()); 00185 } 00186 00187 00188 void ALambdaManager::PrintSomeSpaces() { 00189 //------------------------------------------------------------------------ 00190 // some stuff to make the output look nice 00191 //------------------------------------------------------------------------ 00192 #if !defined(_VERBOSE_PMF) 00193 iout << " "; 00194 #endif 00195 } 00196 00197 00198 void ALambdaManager::Print_dU_dLambda_Summary(double Sum_dU_dLambdas) { 00199 //------------------------------------------------------------------------ 00200 // print sum of dU/dLambda's for current time-step and 00201 // the accumulation of the above. 00202 //------------------------------------------------------------------------ 00203 char Str[100]; 00204 00205 #if defined(_VERBOSE_PMF) 00206 iout << "FreeEnergy: "; 00207 iout << "For all forcing restraints, dU/dLambda = "; 00208 iout << Sum_dU_dLambdas << std::endl << endi; 00209 iout << "FreeEnergy: "; 00210 iout << "For all forcing restraints, Free Energy = "; 00211 iout << GetAccumulation(); 00212 iout << " for " << GetNum_dU_dLambda() << " steps" << std::endl << endi; 00213 #else 00214 sprintf(Str, "%10.2e", GetAccumulation()); 00215 iout << Str << " "; 00216 sprintf(Str, "%6d", GetNum_dU_dLambda()); 00217 iout << Str << " "; 00218 #endif 00219 } 00220 00221 00222 void ALambdaManager::Print_MCTI_Integration() { 00223 //------------------------------------------------------------------------ 00224 // print the integral of: <dU/dLambda> * dLambda 00225 //------------------------------------------------------------------------ 00226 iout << "FreeEnergy: "; 00227 iout << "For MCTI, Free Energy Integral = "; 00228 iout << GetIntegration(); 00229 iout << " for " << GetNumAccumStepsSoFar() << " steps" << std::endl << endi; 00230 } 00231 00232 00233 Bool_t ALambdaManager::GetLambdas(double& LambdaKf, double& LambdaRef) { 00234 //------------------------------------------------------------------------ 00235 // get LambdaKf and LambdaRef from the active LambdaControl 00236 // 00237 // return(kTrue) if an active LambdaControl is found 00238 // return(kFalse) if all LambdaControls have expired 00239 //------------------------------------------------------------------------ 00240 // don't continue if all LamdaControl's have expired 00241 if (m_ActiveIndex == m_NumObjects) { 00242 return(kFalse); 00243 } 00244 00245 // if the m_ActiveIndex'th LambdaControl is no longer active 00246 if ( !(*this)[m_ActiveIndex].IsActive()) { 00247 // move on to the next one 00248 m_ActiveIndex++; 00249 // if there is no next object, return KFalse 00250 if (m_ActiveIndex == m_NumObjects) { 00251 return(kFalse); 00252 } 00253 // otherwise, make sure the next one's active 00254 else { 00255 ASSERT( (*this)[m_ActiveIndex].IsActive() ); 00256 } 00257 } 00258 // return LambdaKf and LambdaRef from the active LambdaControl 00259 LambdaKf = (*this)[m_ActiveIndex].GetLambdaKf(); 00260 LambdaRef = (*this)[m_ActiveIndex].GetLambdaRef(); 00261 return(kTrue); 00262 } 00263 00264 00265 void ALambdaManager::Clear() { 00266 //------------------------------------------------------------------------ 00267 // leave memory allocation alone. 00268 //------------------------------------------------------------------------ 00269 m_NumObjects = 0; 00270 } 00271 00272 00273 int ALambdaManager::Add(ALambdaControl& PmfBlock) { 00274 //------------------------------------------------------------------------ 00275 // add an object to the list. if there's not enough room, make room. 00276 // return an index to the added oject. 00277 //------------------------------------------------------------------------ 00278 ALambdaControl* pPmfBlocks; 00279 00280 // if there's no room for another object 00281 if (m_NumObjects == m_MaxNum) { 00282 // create an array with more space 00283 m_MaxNum *= kLambdaMultiplier; 00284 pPmfBlocks = new ALambdaControl[m_MaxNum]; 00285 // copy from the full array to the new one 00286 for (int i=0; i<m_NumObjects; i++) { 00287 pPmfBlocks[i] = m_pPmfBlocks[i]; 00288 } 00289 // return the space used for the full array 00290 delete []m_pPmfBlocks; 00291 // point to the bigger array 00292 m_pPmfBlocks = pPmfBlocks; 00293 } 00294 // add the object to the array 00295 m_pPmfBlocks[m_NumObjects] = PmfBlock; 00296 m_NumObjects++; 00297 return(m_NumObjects-1); 00298 } 00299 00300 00301 ALambdaControl& ALambdaManager::operator[] (int Index) { 00302 //------------------------------------------------------------------------ 00303 // return an object from this group of objects. 00304 //------------------------------------------------------------------------ 00305 ASSERT((Index>=0) && (Index<m_NumObjects)); 00306 return(m_pPmfBlocks[Index]); 00307 } 00308 00309 00310 int ALambdaManager::GetTotalNumSteps() { 00311 //------------------------------------------------------------------------ 00312 // calculate and return the total number of steps needed for all 00313 // pmf and mcti blocks 00314 //------------------------------------------------------------------------ 00315 int Total, i; 00316 00317 Total = 0; 00318 for (i=0; i<m_NumObjects; i++) { 00319 Total += m_pPmfBlocks[i].GetNumSteps(); 00320 } 00321 return(Total); 00322 } 00323
1.3.9.1