[Wien] Problems when trying to plot E vs c/a

[Fecher, Gerhard]

Dear c/a fitters,
This concerns the latest Wien2k version
I receive only the content of 
       test_opt.analysis 
when I try with w2web to plot E vs c/a 
but neither the result of the fit nor the plot are shown,
this seems to be a problem with the present version of the 
    eplot
script

when I use the eplot script of version 14.2 it is nearly ok, however,
there are still two issues: instead of the result of the fit, the "correlation matrix of the fit parameters" is shown
and the figure is missing.
Reason is that eplot and optimize.pl do not work well together:
     optimize.pl
prints the last 5 lines of fit.log (but the result is before these lines) and expects the graph as case.c_over_a.png (but has a different name .coa.)

this can be solved by changing the two lines
line 169        change "CASE.c_over_a.png"
	    $umps = qx(cp $DIR/$CASE.coa.png $tempdir/$SID-$$.png);    
line 173        change "tail -5"
	    $OUT .= qx(cd $DIR;echo '  ';echo "Fit of:  E = a1 + a2*x + a3*x^2 + a4*x^3 + a5*x^4";tail -15 fit.log);

or indeed, by changing eplot (I just did not find fast how to supress the output of the correlation matrix)

It would also nice to have the minimum (should be prepared by findmincoa called in eplot) in the fit.log and w2web output.



Ciao
Gerhard

DEEP THOUGHT in D. Adams; Hitchhikers Guide to the Galaxy:
"I think the problem, to be quite honest with you,
is that you have never actually known what the question is."

====================================
Dr. Gerhard H. Fecher
Institut of Inorganic and Analytical Chemistry
Johannes Gutenberg - University
55099 Mainz
and
Max Planck Institute for Chemical Physics of Solids
01187 Dresden
________________________________________
Von: Wien [wien-bounces at zeus.theochem.tuwien.ac.at] im Auftrag von Gavin Abo [gsabo at crimson.ua.edu]
Gesendet: Dienstag, 15. Mai 2018 06:40
An: wien at zeus.theochem.tuwien.ac.at
Betreff: Re: [Wien] Finite temperature DFT: electronic free energy

I think you found a bug in init_lapw.  The fix seems like it should be simple though.  On line 498 in the init_lapw script in the sed command, it looks like $fermit just needs changed to $fermits.

I made a patch file called init_lapw.patch [ https://github.com/gsabo/WIEN2k-Patches/tree/master/17.1 ], which you should be able to apply and try simply using:

username at computername:~/Desktop$ cd $WIENROOT
username at computername:~/WIEN2k$ wget https://raw.githubusercontent.com/gsabo/WIEN2k-Patches/master/17.1/init_lapw.patch
username at computername:~/WIEN2k$ patch -b init_lapw init_lapw.patch
patching file init_lapw

Maybe F is the :ENE (TOTAL ENERGY) in case.scf.

I tried following what a previous user did [ https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg10319.html ] but with a quick TiC calculation using run_lapw in WIEN2k 17.1:

T = 1000 K = 0.00633 Ry
init_lapw -b -fermits 0.00633
TiC.scf::ENE  : ********** TOTAL ENERGY IN Ry =        -1783.95041939
TiC.output2:          Kb * T            =   0.00633000
TiC.output2:          -T*Entr           =  -0.00054181

T = 0 K = 0.0 Ry
init_lapw -b -fermits 0.0
TiC.scf::ENE  : ********** TOTAL ENERGY IN Ry =        -1783.94928338
TiC.output2:          Kb * T            =   0.00180000
TiC.output2:          -T*Entr           =  -0.00005329

F_1000K - F_0K = (E - T*S) - (E - 0) = -T*S = -1783.95041939 - (-1783.94928338) = -0.001136 => -T*S = -0.001136

If -T*S divided by 2:

-T*S = -0.001136/2 = - 0.000568 <- This seems reasonably close to the above -0.00054181 for T = 1000 K.

As I recall, the -T*S term was doubled in WIEN2k versions after 2014 [ https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg10319.html ].

Though, there may be a flaw in my above calculation.

On 5/14/2018 5:35 PM, Sabry Moustafa wrote:

Hi;

I am interested in the electronic free energy (F=E-TS) using the finite temperature DFT using Fermi-Dirac statistics (i.e., Mermin's extension to 0 K DFT) -- i.e., I am not just looking for Fermi-Dirac as a broadening/smearing technique. So,how this can be done in WIEN2k?

It looks like I have to define "-fermits X" in the init_lapw command (where X= kT, in Ry). But when I did (T=1000K in my case):

init_lapw -b -fermits 0.00633

I got "fermit: Undefined variable." at the end. This is fixed though when defined fermit as well:

init_lapw -b -fermit 0.00633 -fermits 0.00633

So, do I have to define both?


Also, where can I find this "free energy F". It does not seem to be the "TOTAL ENERGY" in case.scf. Seems like I need to add this "TOTAL ENERGY" (E) to "-T*Entr" (-TS) in case.output2 file to get F?


Thanks ;

Sabry

:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:
Sabry G. Moustafa, Ph.D.
Department of Chemical and Biological Engineering,
University at Buffalo, The State University of New York.
511 Furnas Hall
Buffalo, NY 14260-4200
716-645-1186 (office)
716-239-8543 (cell)
E-mail: sabrygad at buffalo.edu<mailto:sabrygad at buffalo.edu>