Re: How to calculate the temperature for a subset of my simulation

From: Josh Vermaas (vermaas2_at_illinois.edu)
Date: Tue Oct 20 2015 - 14:03:32 CDT

Hi Jhonatam,

Do the units make sense? I've never used a veldcd before, but mass *
velocity squared does not by default give you something in kcal/mol, and
I'm not sure if NAMD does any sort of conversion for you. That would be
the first place I'd look, since the rough idea behind the script
(temperature is related to kinetic energy) is sound.

-Josh Vermaas

On 10/20/2015 01:47 PM, Jhonatam Cordeiro Rodrigues wrote:
> Dear NAMD users
>
> After using the right files (.veldcd). The script generated values
> that are also incorrect. After running the script I got values like this:
>
> nDoF: 481968, esum: 638805356005.8174
> Frame: 201 ## Temp: 4.08764523
>
> Please, what am I missing?
>
> Sincerely
> Jhonatam Cordeiro,
> Ph.D. candidate
> Department of Industrial & Systems Engineering
> North Carolina A&T State University
> 1601, East Market st.
> Greensboro, NC 27411
> Phone: +1 336 706-9203
> Fax: +1 336 334-7729
> email: jcrodrig_at_aggies.nc
> <https://urldefense.proofpoint.com/v2/url?u=http-3A__goog-5F74760233&d=BQMFaQ&c=8hUWFZcy2Z-Za5rBPlktOQ&r=zFfoK61upjM5BwyoRAsX8dLq7rwWm8aw7r7dqtjgcCE&m=3C6-QVPKgRQ7M8Bs-Gvy9p2AW6cxmLMcotfGrLka3m0&s=Mc3QeXblLqH_NG4qH4RFbU6jdQKUSphWh_8UmG399Ec&e=>at.edu
> <https://urldefense.proofpoint.com/v2/url?u=http-3A__at.edu&d=BQMFaQ&c=8hUWFZcy2Z-Za5rBPlktOQ&r=zFfoK61upjM5BwyoRAsX8dLq7rwWm8aw7r7dqtjgcCE&m=3C6-QVPKgRQ7M8Bs-Gvy9p2AW6cxmLMcotfGrLka3m0&s=aYU4-hwn6__w1Nz0UXBEJ_Xsd6iN6muL6oXbnqhHreA&e=>
>
> On Tue, Oct 20, 2015 at 1:20 PM, Jhonatam Cordeiro Rodrigues
> <jcrodrig_at_aggies.ncat.edu <mailto:jcrodrig_at_aggies.ncat.edu>> wrote:
>
> Hello Norman
>
>
> Thank you for the suggestion. I was using the wrong file. I was
> loading the .dcd instead of the .veldcd.
>
>
> Thankfully
> Jhonatam Cordeiro,
> Ph.D. candidate
> Department of Industrial & Systems Engineering
> North Carolina A&T State University
> 1601, East Market st.
> Greensboro, NC 27411
> Phone: +1 336 706-9203 <tel:%2B1%20336%20706-9203>
> Fax: +1 336 334-7729
> email: jcrodrig_at_aggies.nc
> <https://urldefense.proofpoint.com/v2/url?u=http-3A__goog-5F74760233&d=BQMFaQ&c=8hUWFZcy2Z-Za5rBPlktOQ&r=zFfoK61upjM5BwyoRAsX8dLq7rwWm8aw7r7dqtjgcCE&m=3C6-QVPKgRQ7M8Bs-Gvy9p2AW6cxmLMcotfGrLka3m0&s=Mc3QeXblLqH_NG4qH4RFbU6jdQKUSphWh_8UmG399Ec&e=>at.edu
> <https://urldefense.proofpoint.com/v2/url?u=http-3A__at.edu&d=BQMFaQ&c=8hUWFZcy2Z-Za5rBPlktOQ&r=zFfoK61upjM5BwyoRAsX8dLq7rwWm8aw7r7dqtjgcCE&m=3C6-QVPKgRQ7M8Bs-Gvy9p2AW6cxmLMcotfGrLka3m0&s=aYU4-hwn6__w1Nz0UXBEJ_Xsd6iN6muL6oXbnqhHreA&e=>
>
> On Tue, Oct 20, 2015 at 2:23 AM, Norman Geist
> <norman.geist_at_uni-greifswald.de
> <mailto:norman.geist_at_uni-greifswald.de>> wrote:
>
> Did you actually used a veldcd ?
>
> Norman Geist
>
> *Von:*owner-namd-l_at_ks.uiuc.edu
> <mailto:owner-namd-l_at_ks.uiuc.edu>
> [mailto:owner-namd-l_at_ks.uiuc.edu
> <mailto:owner-namd-l_at_ks.uiuc.edu>] *Im Auftrag von *Jhonatam
> Cordeiro Rodrigues
> *Gesendet:* Montag, 19. Oktober 2015 22:08
> *An:* Namd Mailing List <namd-l_at_ks.uiuc.edu
> <mailto:namd-l_at_ks.uiuc.edu>>
> *Betreff:* namd-l: How to calculate the temperature for a
> subset of my simulation
>
> Dear NAMD users
>
> I am running simulations that have different temperature
> couplings on different subsets of the simulation. Please, how
> can I read the temperature for a given selection of atoms?
>
> I found this nice script written by Shailesh Pandey, but when
> I run it, I get temperature readings that are certainly
> incorrect (over 1Bi K).
>
> Is there something wrong with the script? or is there a script
> out there that would do the trick?
>
> ########### TCL script to calculate temperature from velocity
> trajectorySTART #########
>
> # rigidbonds none=0; water=1; all=2
>
> proc getTemp { frame rigidbonds } {
>
> set all [atomselect top all frame $frame]
> set eng {};
> set esum 0.0;
> set nDoF 0;
> set c1 [[ atomselect top "((type 'HW' or hydrogen) or
> water)" ] get index ]
> set c2 [[ atomselect top "water" ] get index ]
>
>
> foreach m [$all get mass] v [$all get {x y z}] id [$all
> get index] {
>
> set e [expr {0.5 * $m * [vecdot $v $v]}]
> # lappend eng $e
> set esum [expr {$esum + $e}]
> # set all [atomselect top all frame $frame]
> set DoF 3
>
> #SHAKE?
> if { ($rigidbonds == 2) && ($id in $c1) } {
> set DoF 2
> } elseif { ($rigidbonds == 1) && ($id in $c2) } {
> set DoF 2
> }
>
> set nDoF [expr {$nDoF + $DoF}]
> }
> $all delete
> puts "nDoF: $nDoF, esum: $esum"
> set kB 0.00198657 ;# kcal/(mol.K)
> set T [expr { $esum * 2.0 / ($nDoF * $kB) }]
> puts "Frame: $frame ## Temp: $T"
> }
>
> ########### TCL script to calculate temperature from velocity
> trajectory END #########
>
> After the running the script for water I get results like this:
> nDoF: 481968, esum: 638805356005.8174
> Frame: 201 ## Temp: 1334370608.0432146
>
> Sincerely
>
> Jhonatam Cordeiro,
> Ph.D. candidate
> Department of Industrial & Systems Engineering
> North Carolina A&T State University
>
> 1601, East Market st.
>
> Greensboro, NC 27411
> Phone: +1 336 706-9203 <tel:%2B1%20336%20706-9203>
> Fax: +1 336 334-7729 <tel:%2B1%20336%20334-7729>
>
> email: jcrodrig_at_aggies.nc
> <https://urldefense.proofpoint.com/v2/url?u=http-3A__goog-5F74760233&d=BQMFaQ&c=8hUWFZcy2Z-Za5rBPlktOQ&r=zFfoK61upjM5BwyoRAsX8dLq7rwWm8aw7r7dqtjgcCE&m=3C6-QVPKgRQ7M8Bs-Gvy9p2AW6cxmLMcotfGrLka3m0&s=Mc3QeXblLqH_NG4qH4RFbU6jdQKUSphWh_8UmG399Ec&e=>at.edu
> <https://urldefense.proofpoint.com/v2/url?u=http-3A__at.edu&d=BQMFaQ&c=8hUWFZcy2Z-Za5rBPlktOQ&r=zFfoK61upjM5BwyoRAsX8dLq7rwWm8aw7r7dqtjgcCE&m=3C6-QVPKgRQ7M8Bs-Gvy9p2AW6cxmLMcotfGrLka3m0&s=aYU4-hwn6__w1Nz0UXBEJ_Xsd6iN6muL6oXbnqhHreA&e=>
>
>
>

This archive was generated by hypermail 2.1.6 : Tue Dec 27 2016 - 23:21:26 CST