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

From: Norman Geist (norman.geist_at_uni-greifswald.de)
Date: Wed Oct 21 2015 - 02:24:09 CDT

The veldcd requires a conversion factor to bring the velocities to angstrom/picosecond.

 

[wrap]

set PDBVELFACTOR 20.45482706

foreach m [$all get mass] v [$all get {x y z}] id [$all get index] {

        set $v [vecscale $v $PDBVELFACTOR]
        set e [expr {0.5 * $m * [vecdot $v $v]}]
        # lappend eng $e
[wrap]

 

Maybe this works?

 

Norman Geist

 

Von: owner-namd-l_at_ks.uiuc.edu [mailto:owner-namd-l_at_ks.uiuc.edu] Im Auftrag von Josh Vermaas
Gesendet: Dienstag, 20. Oktober 2015 21:04
An: namd-l_at_ks.uiuc.edu; Jhonatam Cordeiro Rodrigues <jcrodrig_at_aggies.ncat.edu>; Norman Geist <norman.geist_at_uni-greifswald.de>
Betreff: Re: namd-l: How to calculate the temperature for a subset of my simulation

 

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 <mailto:jcrodrig_at_aggies.nc> 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 <mailto:jcrodrig_at_aggies.nc> 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