RE: Output of thermodynamic integration

From: Okazaki, Keiichi (NIH/NIDDK) [F] (keiichi.okazaki_at_nih.gov)
Date: Sun Jan 29 2012 - 14:22:26 CST

Dear Jerome and developers,

Thank you for your help!
I think there is a similar problem also for vdW_dU/dl.
Since vdW is turned on first for 0<=lambda<=0.5 in my case, I expect that vdW doesn't contribute for lambda>0.5, where LJ potential is independent on lambda.
But it's not in output file, for example, lambda=0.9 shown below.
If my guess is right, I think it is not a serious bug, but confusing.
Can I just ignore Elec_dU/dl for lambda<0.5 and vdW_dU/dl for lambda>0.5??

Regards,
Kei-ichi

## output file ###################################################################################
# STEP Elec_dU/dl Elec_avg vdW_dU/dl vdw_avg Elec_dU/dl Elec_avg vdW_dU/dl vdw_avg PME_dU/dl PME_avg
# <---------------------PARTITION 1------------------------> <---------------------PARTITION 2--------------------->
#NEW TI WINDOW: LAMBDA 0.9
#PARTITION 1 VDW LAMBDA 1
#PARTITION 1 ELEC LAMBDA 0.8
#PARTITION 2 VDW LAMBDA 0.2
#PARTITION 2 ELEC LAMBDA 0

TI 0 -452.9840 -452.9840 97.2694 97.2694 0.0000 0.0000 0.0000 0.0000
TI 1000 -390.6807 -390.7429 66.7038 66.7344 0.0000 0.0000 0.0000 0.0000

________________________________________
From: heninj_at_gmail.com [heninj_at_gmail.com] On Behalf Of Jérôme Hénin [jhenin_at_ifr88.cnrs-mrs.fr]
Sent: Saturday, January 28, 2012 1:52 PM
To: Okazaki, Keiichi (NIH/NIDDK) [F]
Cc: namd-l_at_ks.uiuc.edu
Subject: Re: namd-l: Output of thermodynamic integration

You are right, I should have seen it in the output you attached. This
is beginning to look like a bug. I don't know the TI code well, and I
don't have the time to go into the details right now. I hope that
another developer can chime in...

Best,
Jerome

On 28 January 2012 17:31, Okazaki, Keiichi (NIH/NIDDK) [F]
<keiichi.okazaki_at_nih.gov> wrote:
> Dear Jerome,
>
> Thank you for your suggestion.
> However, 1 is set in my alchFile for the flag, that is, it goes in the same direction with the global lambda.
> So, I am still not sure why Elec_dU/dl takes some values even though ele is not turned on yet.
>
> Regards,
> Kei-ichi
> ________________________________________
> From: heninj_at_gmail.com [heninj_at_gmail.com] On Behalf Of Jérôme Hénin [jhenin_at_ifr88.cnrs-mrs.fr]
> Sent: Saturday, January 28, 2012 5:25 AM
> To: Okazaki, Keiichi (NIH/NIDDK) [F]
> Cc: namd-l_at_ks.uiuc.edu
> Subject: Re: namd-l: Output of thermodynamic integration
>
> Dear Kei-ichi,
>
> If you are decoupling the solute from the solvent, it probably has -1
> flags in the alchFile. In this case, the lambdas for the solute are
> reversed with respect to the "global" lambda, and the first thing that
> happens is that electrostatics gets decoupled (which is the right
> thing to do).
>
> Best,
> Jerome
>
>
> On 28 January 2012 01:16, Kei-ichi Okazaki <keiichi.okazaki_at_nih.gov> wrote:
>> Dear NAMD users,
>>
>> I am trying to calculate solvation free energy of a small molecule by
>> thermodynamic integration (TI).
>> And I have a question about output file of TI.
>> I set up configuration file for TI with soft-core potential as below.
>> I set "alchElecLambdaStart 0.5" and "alchVdwLambdaEnd 0.5" to separate the
>> scaling of vdW and electrostatics, that is, vdW is turned on first for
>> 0<=lambda<=0.5 and then, ele is turned on for 0.5<=lambda<=1.0.
>> And in the output, Elec_dU/dl and vdW_dU/dl are reported as below.
>> When I saw the output of lambda=0.1, for example, I was confused because
>> Elec_dU/dl takes some values even though ele is not turned on yet.
>> I suspect that it just reports Elec_dU/dl for lambda=0.0, but it cannot be
>> used for TI calculation because there should not be any contribution from
>> ele for lambda<0.5.
>> Is my guess right? I appreciate any suggestions.
>>
>> Thanks in advance.
>> Kei-ichi Okazaki, PhD
>>
>> ## configuration file #########################################
>> # Free energy perturbation
>> alch On # Enable alchemical simulation module
>> alchType ti # Set method to thermodynamic integration
>> alchFile ../../Common/h2p.box20.alchem.pdb # PDB file with
>> perturbation flags
>> alchCol B # Perturbation flags in Beta column
>> alchOutfile $outputname.out
>> alchOutFreq 1000
>> alchEquilSteps 500000
>>
>> alchVdWShiftCoeff 6.0 # Enable soft-core vdW potential
>> alchElecLambdaStart 0.5 # Introduce electrostatics for lambda > 0.0
>> alchVdwLambdaEnd 0.5
>> alchLambda 0.1
>> alchLambda2 0.1
>>
>> ## output file #################################################
>> # STEP Elec_dU/dl Elec_avg vdW_dU/dl vdw_avg
>> Elec_dU/dl Elec_avg vdW_dU/dl vdw_avg PME_dU/dl
>> PME_avg
>> # <---------------------PARTITION 1------------------------>
>> <---------------------PARTITION 2--------------------->
>> #NEW TI WINDOW: LAMBDA 0.1
>> #PARTITION 1 VDW LAMBDA 0.2
>> #PARTITION 1 ELEC LAMBDA 0
>> #PARTITION 2 VDW LAMBDA 1
>> #PARTITION 2 ELEC LAMBDA 0.8
>>
>> TI 0 -452.9840 -452.9840 -9.6522 -9.6522
>> 0.0000 0.0000 0.0000 0.0000
>> TI 1000 -228.3854 -228.6098 1.7409 1.7296
>> 0.0000 0.0000 0.0000 0.0000
>>
>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 05:24:41 CST