File: lesson_rf1.html

package info (click to toggle)
abinit 7.8.2-2
  • links: PTS, VCS
  • area: main
  • in suites: jessie, jessie-kfreebsd
  • size: 278,292 kB
  • ctags: 19,095
  • sloc: f90: 463,759; python: 50,419; xml: 32,095; perl: 6,968; sh: 6,209; ansic: 4,705; fortran: 951; objc: 323; makefile: 43; csh: 42; pascal: 31
file content (757 lines) | stat: -rw-r--r-- 32,857 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
<html>
<head>
<title>Tutorial Response-Function 1</title>
</head>
<body bgcolor="#ffffff">

<hr>
<h1>ABINIT, lesson response-function 1 (RF1) of the tutorial: </h1>
<h2>Dynamical and dielectric properties of AlAs.
</h2>
<hr>
<p>This lesson aims at showing how to get
the following physical properties, for an insulator :
<ul>
<li>the phonon frequencies and eigenvectors at Gamma
<li>the dielectric constant
<li>the Born effective charges
<li>the LO-TO splitting
<li>the phonon frequencies and eigenvectors at other q-points in the Brillouin Zone
</ul>
You will learn to use of response-function features of ABINIT.
In order to learn the use of the associated codes Mrgddb and Anaddb,
to produce efficiently phonon band structures and the associated thermodynamical properties,
please see the <a href="lesson_rf2.html"> tutorial response-function 2</a>.

<p>This lesson should take about 2 hours.

<h5>Copyright (C) 2000-2014 ABINIT group (XG,RC)
<br> This file is distributed under the terms of the GNU General Public License, see
~abinit/COPYING or <a href="http://www.gnu.org/copyleft/gpl.txt">
http://www.gnu.org/copyleft/gpl.txt </a>.
<br> For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
</h5>

<script type="text/javascript" src="list_internal_links.js"> </script>

<br>&nbsp;
<h3><b>Content of lesson RF1</b></h3>

<ul>
  <li><a href="lesson_rf1.html#1">1</a> The ground-state geometry of AlAs.</li>
  <li><a href="lesson_rf1.html#2">2</a> Frozen-phonon calculation of a second derivative of the total energy.</li>
  <li><a href="lesson_rf1.html#3">3</a> Response-function calculation of a second derivative of the total energy.</li>
  <li><a href="lesson_rf1.html#4">4</a> Response-function calculation of the dynamical matrix at Gamma.</li>
  <li><a href="lesson_rf1.html#5">5</a> Response-function calculation of the effect of an homogeneous electric field.</li>
  <li><a href="lesson_rf1.html#6">6</a> Response-function calculation of phonon frequencies at non-zero q.</li>
</ul>

<hr>

<p><a name="1"></a><br>

<h3><b>1. The ground-state geometry of AlAs.</b> </h3>
<h3>&nbsp;</h3>

<p>

<i>Before beginning, you might consider to work in a different subdirectory as
for the other lessons. Why not create "Work_rf1" in ~abinit/tests/tutorespfn/Input ? </i> <br>

<p>NOTE : the reference directory that contains
  the example files for the tutorial is no more ~abinit/tests/tutorial (as for the
  basic lessons and the specialized, non-response ones), but
  ~abinit/tests/tutorespfn . This will be the case for all the response-function
  based part of the tutorial.

<p>The file ~abinit/tests/tutorespfn/Input/trf1_x.files lists the file
  names and root names.
  You can copy it in the Work_rf1 directory (and change it, as usual).
  Note that two pseudopotentials are mentioned in this "files" file:
  one for the Aluminum atom, and one for the Arsenic atom. The first
  to be mentioned, for Al, will define the first type of atom. The second
  to be mentioned, for As, will define the second type of atom.
  It might the first time that you encounter this situation (more than one type of atoms) in the tutorials, at variance with the four "basic" lessons.
<p>
  You can also copy the file ~abinit/tests/tutorespfn/Input/trf1_1.in in Work_rf1.
  This is your input file. You
  should edit it, and read it carefully.
<p>
  Because of the use of
  two types of atoms, have a look at the following input variables :
<ul>
  <li><a href="../input_variables/varbas.html#ntypat" target="kwimg">ntypat</a>
  <li><a href="../input_variables/varbas.html#typat" target="kwimg">typat</a>
</ul>
<p>
  Note that the value of <a href="../input_variables/varbas.html#tolvrs" target="kwimg">tolvrs</a>
  is rather stringent. This is because the wavefunctions determined by the
  present run will be used later as starting point of the response-function
  calculation. However, the number of steps, <a href="../input_variables/varbas.html#nstep" target="kwimg">nstep</a>,
  in this example file has been set to 15, and you will see that this is not enough to reach
  the target <a href="../input_variables/varbas.html#tolvrs" target="kwimg">tolvrs</a>. In production runs,
  you should choose a larger value of <a href="../input_variables/varbas.html#nstep" target="kwimg">nstep</a>,
  sufficient to reach your target <a href="../input_variables/varbas.html#tolvrs" target="kwimg">tolvrs</a>.
  In the present tutorial, due to portability concerns related to automatic testing,
  we could not allow a larger <a href="../input_variables/varbas.html#nstep" target="kwimg">nstep</a> value.
  This minor problem with some tutorial examples was mentioned briefly in
  <a href="lesson_base1.html#aq1">a side note to the answer to question 1 of lesson 1</a>.
  So, do not follow blindly all examples in the tutorials : check by yourself
  the convergence of your runs !
<p>
You will work at fixed <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>
(=3Ha) and k-point grid, defined
by <a href="../input_variables/vargs.html#kptrlatt" target="kwimg">kptrlatt</a>
(the 8x8x8 Monkhorst-Pack grid).
It is implicit that in "real life", you should do a convergence test with
respect to both convergence parameters.
We postpone the discussion of the accuracy of these choices and
the choice of pseudopotential to the end of the <a href="#5">fifth section of this tutorial</a>.
They give acceptable results, not very accurate,
but, more important, the speed is reasonable for a tutorial.
<br>

<p>
You should make the run (a few seconds),
and obtain the following value for the energy,
in the final echo section :
<pre>
 etotal   -9.7626837450E+00
</pre>
However, we will rely later on a more accurate (more digits) value of
this total energy, that can be found about a dozen of lines
before this final echo :
<pre>
    >>>>>>>>> Etotal= -9.76268374500289E+00
</pre>

<p>
The output file also mentions that the forces on both atoms vanish.

<p>The run that you just made will be considered as defining a ground-state
configuration, on top of which response functions will be computed.
The main output of this ground-state run is the wavefunction file
trf1_1o_WFK, that you can already rename as trf1_1i_WFK.
Indeed, it will be used in the next runs, as an input file.
<b>So, in the corresponding "files" file for all the following runs, at third line,
pay attention TO KEEP "trf1_1i", even if you change
the root name for output files (fourth line) to "trf1_2o" or "trf1_3o", as well as
the first, second and fifth lines of this file.</b>


<p>
<hr>
<p><a name="2"></a><br>

<h3><b>2. Frozen-phonon calculation of a second derivative of the total energy.</b> </h3>
<h3>&nbsp;</h3>

<p>
We will now aim at computing the second derivative of the total energy
with respect to an atomic displacement by different means.
For that purpose, you must first read
<a href="../users/respfn_help.html#0" target="helpsimg">sections 0 and
the first paragraph of section 1</a>
of the respfn_help file (the auxiliary help file,
that deals specifically with the response function features).

<p>
We will explain later, in more detail, the signification
of the different input parameters introduced in section 1 of the
respfn_help file.

<p>
For now, in order to be able to perform a direct comparison with the
result of a response-function calculation, we choose as a perturbation
the displacement of the
Al atom along the first axis of the reduced coordinates.

<p>
You can copy the file ~abinit/tests/tutorespfn/Input/trf1_2.in in Work_rf1.
This is your input file. You should edit it and briefly
look at the two changes with respect to the file
 ~abinit/tests/tutorespfn/Input/trf1_1.in : the change of
<a href="../input_variables/varbas.html#xred" target="kwimg">xred</a>, and the
reading of the wavefunction file, using the
<a href="../input_variables/varfil.html#irdwfk" target="kwimg">irdwfk</a> input variable.
Then, you can make the run. The symmetry is lowered with respect
to the ground-state geometry, so that the number of k-points
increases a lot, and of course, the CPU time. 

<p>
From this run, it is possible to get the values of the total energy,
and the value of the gradient of the total energy (dE) with respect to
change of reduced coordinate (dt):
<pre>
 rms dE/dt=  3.5517E-03; max dE/dt=  5.0079E-03; dE/dt below (all hartree)
    1       0.005007930232      0.002526304574      0.002526304574
    2      -0.005007868293     -0.002526274885     -0.002526274885
 ...
    >>>>>>>>> Etotal= -9.76268124105780E+00
</pre>
<p>
The change of reduced coordinate of the Al atom along the first axis
was rather small (1/1000), and we can make an estimate of the second derivative
of the total energy with respect to the reduced coordinate thanks
to finite-difference formulas.

<p> We start first from the total-energy difference. The total energy is symmetric with respect to that
perturbation, so that it has no linear term. The difference
between the ground-state value (-9.76268374500102E+00 hartree)
of the previous run,
and the perturbed value (-9.76268124105590E+00 hartree) of the present
one, is thus
one half of the square of the coordinate change (1/1000) times
the second derivative of total energy (2DTE). From these number, the 2DTE is 5.00791024 hartree.

<p> Alternatively, we can start from the reduced gradients.
The value of the reduced gradient with respect to
a displacement of the Al atom along the first reduced axis
is 0.005007930232 (hartree). At first order, this
quantity is the product of the 2DTE by the reduced coordinate
difference. The estimate of the 2DTE is thus 5.007930232 hartree.
The agreement with the other estimate is rather good (2.10^-5 Hartree).

<p> However, it is possible to do much better, thanks to the
use of a higher-order finite-difference formula.
For this purpose, one can perform another calculation, in
which the change of reduced coordinate along the first
axis is 0.002, instead of 0.001. The doubling
of the perturbation allows for a rather simple higher-order
estimation, as we will see later.
The results of this calculation are as follows :
<pre>
 rms dE/dt=  7.1249E-03; max dE/dt=  1.0016E-02; dE/dt below (all hartree)
    1       0.010016299779      0.005097509981      0.005097509981
    2      -0.010016176675     -0.005097455174     -0.005097455174
 ...
</pre>
<p>
From these data, taking into account that the perturbation
was twice stronger, the same procedure as above leads
to the values 5.00800270 hartree (from finite difference of energy)
and 5.009149889 hartree (from finite difference of forces,
the value 0.010016299779 has to be multiplied by 1000/2).
The combination of these data with the previous estimate
can be done thanks to an higher-order finite-difference formula,
in which the difference of estimations (the largest perturbation
minus the smallest one) is divided by three,
and then subtracted from the smallest estimation.
As far as the total-energy estimation is concerned, the difference
is 0.0001104 Ha, which divided by three,
and subtracted from 5.00789230 hartree,
gives 5.0078555 hartree. The same higher-order procedure for
force estimates gives 5.0078552 hartree. So, the agreement
between total-energy estimate and force estimate of the 2DTE can be
observed up to the 7th digit, inclusive.

<p> Before comparing this result with the 2DTE directly
computed from the response-function capabilities of ABINIT,
a last comment is in order. One can observe that the
action-reaction law is fulfilled only approximately
by the system. Indeed, the force created on the second atom,
should be exactly equal in magnitude to the force on the first atom.
The values of dE/dt, mentioned above, for example for the
doubled displacement :
<pre>
 rms dE/dt=  7.1249E-03; max dE/dt=  1.0016E-02; dE/dt below (all hartree)
    1       0.010016299779      0.005097509981      0.005097509981
    2      -0.010016176675     -0.005097455174     -0.005097455174
</pre>
<p>
show a small, but non-negligible difference between the
two atoms. Actually, the forces should cancel each other exactly
if the translation symmetry were perfect. This is not the case,
but the breaking of this symmetry can be shown to arise <b>only</b>
from the presence of the exchange-correlation grid of points.
This grid does not move when atoms are displaced, and so there
is a very small variation of the total energy when the
system is moved as a whole. It is easy to restore the
action-reaction law, by subtracting from every force component
the mean of the forces on all atoms. This is actually done when
the gradient with respect to reduced coordinates are transformed
into forces, and specified in cartesian coordinates,
as can be seen in the output file for the small displacement :
<pre>
 cartesian forces (hartree/bohr) at end:
    1     -0.00000421123276    -0.00047199792687    -0.00047199792687
    2      0.00000421123276     0.00047199792687     0.00047199792687
</pre>
This effect will be seen also at the level of 2DTE. The so-called
"acoustic sum rule", imposing that the frequency
of three modes (called acoustic modes) tend
to zero with vanishing wavevector, will also be slightly broken.
In this case also, it will be rather easy to reimpose the
acoustic sum rule. In any case, taking a finer XC grid will
allow one to reduce this effect.

<p>
<hr>

<p><a name="3"></a><br>

<h3><b>3. Response-function calculation of a second derivative of the total energy.</b> </h3>
<h3>&nbsp;</h3>

<p>
We now compute the second derivative of the total energy
with respect to the same atomic displacement, using
the response-function capabilities of ABINIT.

<p>
You can copy the file ~abinit/tests/tutorespfn/Input/trf1_3.in in Work_rf1.
This is your input file. You should examine it. The changes
with respect to the file ~abinit/tests/tutorespfn/Input/trf1_1.in
are all gathered in the first part of this file,
before
<pre>
#######################################################################
#Common input variables
</pre>

<p>
Accordingly, you should get familiarized with the new input variables :
<a href="../input_variables/varrf.html#rfphon" target="kwimg">rfphon</a>,
<a href="../input_variables/varrf.html#rfatpol" target="kwimg">rfatpol</a>,
<a href="../input_variables/varrf.html#rfdir" target="kwimg">rfdir</a>.
Then, pay attention to the special use of the
<a href="../input_variables/varbas.html#kptopt" target="kwimg">kptopt</a>
input variable. It will be explained in more detail later.

<p>
When you have understood the purpose of the input variable values
specified before the "Common input variables" section, you can make
the code run. 

<p>
Then, we need to analyze the different output files. For that purpose,
you should read the content of the
<a href="../users/respfn_help.html#6" target="helpsimg"> section 6</a> of the respfn_help
file. Read it quickly, as we will come back to the most important
points hereafter.

<p> ABINIT has created four different files :
<ul>
 <li>trf1_3.log  (the log file)
 <li>trf1_3.out  (the output file)
 <li>trf1_3o_1WF1 (the 1st-order wavefunction file)
 <li>trf1_3o_DDB (the derivative database)
</ul>

<p> Let us have a look at the output file. You can follow the description
provided in the <a href="../users/respfn_help.html#6.2" target="helpsimg"> section 6.2</a>
of the respfn_help
file. You should be able to find the place where the iterations for the
minimisation (with respect to the unique perturbation) take place :
<pre>
      iter   2DEtotal(Ha)       deltaE(Ha) residm    vres2
 ETOT  1   6.5156051312863    -1.464E+01 1.146E-02 1.947E+02
 ETOT  2   5.0216331638978    -1.494E+00 9.267E-04 2.027E+00
 ETOT  3   5.0082678671217    -1.337E-02 4.772E-06 7.929E-02
 ETOT  4   5.0078677958105    -4.001E-04 1.980E-07 2.712E-03
 ETOT  5   5.0078558860285    -1.191E-05 6.103E-09 5.074E-05
 ETOT  6   5.0078557520180    -1.340E-07 1.258E-10 9.613E-06
 ETOT  7   5.0078557017091    -5.031E-08 2.768E-11 3.841E-07
 ETOT  8   5.0078557001188    -1.590E-09 2.983E-12 6.089E-09
</pre>
<p>
From these data, you can see that the 2DTE determined by the
response-function technique is in excellent agreement with the
higher-order finite-difference values for the 2DTE, determined
in the previous section :
5.0078555 hartree from the energy differences,
and 5.0078552 hartree from the force differences.

<p> Now, you can read the remaining of
the <a href="../users/respfn_help.html#6.2" target="helpsimg"> section 6.2</a>
of the respfn_help file. Then, you should also edit the
trf1_3o_DDB file, and read the corresponding
<a href="../users/respfn_help.html#6.5" target="helpsimg"> section 6.5</a>
of the respfn_help file.

<p> Finally, the excellent agreement between the finite-difference
formula and the response-function approach calls for
some accuracy considerations. These can be found
in <a href="../users/respfn_help.html#7" target="helpsimg"> section 7</a>
of the respfn_help file.

<p>
<hr>

<p><a name="4"></a><br>

<h3><b>4. Response-function calculation of the dynamical matrix at Gamma.</b> </h3>
<h3>&nbsp;</h3>

<p> We are now in the position to compute the full dynamical matrix at Gamma (q=0).
You can copy the file ~abinit/tests/tutorespfn/Input/trf1_4.in in Work_rf1.
This is your input file. You should edit it.
As for test rf1_3, the changes
with respect to the file ~abinit/tests/tutorespfn/Input/trf1_1.in
are all gathered in the first part of this file.
Moreover, the changes with respect to trf1_3.in concern
only the input variables
<a href="../input_variables/varrf.html#rfatpol" target="kwimg">rfatpol</a>, and
<a href="../input_variables/varrf.html#rfdir" target="kwimg">rfdir</a>.
Namely, all the atoms will be displaced, in all the directions.

<p>
There are six perturbations to consider.
So, one might think that the CPU time will raise accordingly.
This is not true, as ABINIT is able to determine
which perturbations are the symmetric of another perturbation,
see <a href="../users/respfn_help.html#3" target="helpsimg"> section 3</a>
of the respfn_help file.

<p>
Now, you can make the run. 
You edit the file trf1_4.out, and notice that
the response to two perturbations were computed explicitly,
while the response to the other four could be deduced by using
the symmetries. 
<pre>
 The list of irreducible perturbations for this q vector is:
    1)    idir= 1    ipert=   1
    2)    idir= 1    ipert=   2
</pre>

<p>Nothing mysterious : one of
the two irreducible perturbations is for the Al atom,
placed in a rather symmetric local site, and the other
perturbation is for the As atom.

<p> The phonon frequencies, obtained by diagonalizing
the dynamical matrix (where the atomic masses have been taken
into account, see <a href="../input_variables/varrlx.html#amu" target="kwimg">amu</a>),
are given as follows :
<pre>
  Phonon wavevector (reduced coordinates) :  0.00000  0.00000  0.00000
 Phonon energies in Hartree :
   2.586632E-06  2.590723E-06  2.614440E-06  1.568560E-03  1.568560E-03
   1.568560E-03
 Phonon frequencies in cm-1    :
-  5.677000E-01  5.685980E-01  5.738033E-01  3.442590E+02  3.442590E+02
-  3.442590E+02
</pre>
<p><i>You might wonder about the dash sign present in the first column of the
two lines giving the frequencies in cm-1. The first column of the main ABINIT
output files is always dedicated to signs needed to automatically treat the
comparison with respect to reference files. Except if you become a
ABINIT developer, you should ignore these signs. In the present case,
they should not be interpreted as a minus sign for the floating numbers
that follow them.</i></p>

<p>There is a good news about this result, and a bad news. The good
news is that there are indeed three acoustic modes, with frequency
rather close to zero (less than 1 cm^-1, which is rather good !).
The bad news comes when the three other frequencies are compared with
experimental results, or other theoretical results. Indeed,
in the present run, one obtains three degenerate modes, while there should
be a (2+1) splitting. This can be seen in the paper
<cite> Ab initio calculation of phonon dispersions in semiconductors,
by P. Giannozzi, S. de Gironcoli, P. Pavone and S. Baroni,
Phys. Rev. B 43, 7231 (1991) </cite>, especially Fig. 2.</p>

<p>Actually, we have forgotten to take into account the
coupling between atomic displacements and the homogeneous electric
field, that exists in the case of polar insulators,
for so-called "Longitudinal Optic (LO) modes".
A splitting appears between these modes and the "Transverse Optic (TO) modes".
This splitting (Lyddane-Sachs-Teller LO-TO splitting)
is presented in simple terms in standard textbooks,
and should not be forgotten in doing Ab initio calculations
of phonon frequencies.

<p> Thus we have now to treat correctly the homogeneous electric
field type perturbation.


<p>
<hr>

<p><a name="5"></a><br>

<h3><b>5. Response-function calculation of the effect of an homogeneous electric field.</b> </h3>
<h3>&nbsp;</h3>

<p> The treatment of the homogeneous electric field perturbation is formally
much more complex than the treatment of atomic displacements. This is primarily
because the change of potential associated with an homogeneous electric field
is not periodic, and thus does not satisfy the Born-von Karman periodic
boundary conditions.

<p> For the purpose of the present tutorial, one should read the section II.C of the
above-mentioned paper <cite>
P. Giannozzi, S. de Gironcoli, P. Pavone and S. Baroni,
Phys. Rev. B 43, 7231 (1991) </cite>. The reader will find in
<cite> X. Gonze, Phys. Rev. B 55, 10337 (1997)</cite> and
<cite> X. Gonze and C. Lee, Phys. Rev. B 55, 10355 (1997)</cite> more detailed
information of this perturbation, closely related to the ABINIT implementation.
There is also an extensive discussion of the Born effective charges
by <cite> Ph. Ghosez, J.-P. Michenaud and X. Gonze, Phys. Rev. B 58, 6224 (1998)</cite>.

<p> In order to compute the response of solids to an homogeneous electric field,
as implemented in ABINIT, the remaining sections of the respfn_help file should
be read. These sections also present the information needed to compute phonons with
non-zero q wavevector, which will be the subject of the next section
of the present tutorial. The sections to be read are :
<ul>
 <li>the part of <a href="../users/respfn_help.html#1" target="helpsimg">section 1</a> that
     had not yet been read
 <li><a href="../users/respfn_help.html#2" target="helpsimg">section 2</a>
 <li><a href="../users/respfn_help.html#4" target="helpsimg">section 4</a>
 <li>and, for completeness,<a href="../users/respfn_help.html#5" target="helpsimg">section 5</a>
</ul>


<p> You are now in the position to compute the full dynamical matrix at Gamma (q=0),
including the coupling with an homogeneous electric field.
You can copy the file ~abinit/tests/tutorespfn/Input/trf1_5.in in Work_rf1.
This is your input file. You should edit it.
As for the other RF tests, the changes
with respect to the file ~abinit/tests/tutorespfn/Input/trf1_1.in
are all gathered in the first part of this file.
Unlike the other tests, however, the multi-dataset mode
was used, computing from scratch the ground-state properties,
then computing the effect of the ddk perturbation, then
the effect of all other perturbations (electric field as well as
atomic displacements). 

<p>The analysis of the output file is even more cumbersome than the
previous ones. Let us skip the first dataset. In the dataset 2 section,
one perturbation is correctly selected :
<pre>
 ==>  initialize data related to q vector <==

 The list of irreducible perturbations for this q vector is:
    1)    idir= 1    ipert=   3

================================================================================

--------------------------------------------------------------------------------
 Perturbation wavevector (in red.coord.)   0.000000  0.000000  0.000000
 Perturbation : derivative vs k along direction   1
</pre>
The analysis of the output for this particular perturbation is
not particularly interesting, except for the f-sum rule
ratio
<pre>
 loper3 : ek2=    1.6833336546E+01
          f-sum rule ratio=    1.0028582985E+00
</pre>
that should be close to 1, and becomes closer to it when
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a> is increased,
and the sampling of k points is improved.
(In the present status of ABINIT, the f-rule ratio
is not computed correctly when
<a href="../input_variables/varrlx.html#ecutsm" target="kwimg">ecutsm</a>/=0)

<p> In the third dataset section, three irreducible perturbations
are considered :
<pre>
 ==>  initialize data related to q vector <==

 The list of irreducible perturbations for this q vector is:
    1)    idir= 1    ipert=   1
    2)    idir= 1    ipert=   2
    3)    idir= 1    ipert=   4
</pre>
Much later, the dielectric tensor is given :
<pre>
  Dielectric tensor, in cartesian coordinates,
     j1       j2             matrix element
  dir pert dir pert     real part    imaginary part

   1    4   1    4    9.7606048428    0.0000000000
   1    4   2    4    0.0000000000    0.0000000000
   1    4   3    4    0.0000000000    0.0000000000

   2    4   1    4    0.0000000000    0.0000000000
   2    4   2    4    9.7606048428    0.0000000000
   2    4   3    4    0.0000000000    0.0000000000

   3    4   1    4    0.0000000000    0.0000000000
   3    4   2    4    0.0000000000    0.0000000000
   3    4   3    4    9.7606048428    0.0000000000
</pre>
It is diagonal and isotropic, and corresponds to a dielectric
constant of 9.7606048428 .

<p>
Then, the Born effective charges are given, either
computed from the derivative of the wavefunctions
with respect to the electric field, or
computed from the derivative of the wavefunctions
with respect to an atomic displacement, as explained
in section II of
<cite>X. Gonze, Phys. Rev. B55, 10355 (1997)</cite> :
<pre>
  Effective charges, in cartesian coordinates,
  (from electric field response)
  ...
</pre>
and
<pre>
  Effective charges, in cartesian coordinates,
  (from phonon response)
  ...
</pre>
Namely, the Born effective charge of the Al atom is
2.104, and the one of the As atom is
-2.127 . The charge neutrality sum rule
is not fulfilled exactly. When
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a> is increased,
and the sampling of k points is improved, the
sum of the two charges goes closer to zero.

<p>Finally, the phonon frequencies are computed:
<pre>
  Phonon wavevector (reduced coordinates) :  0.00000  0.00000  0.00000
 Phonon energies in Hartree :
   2.586632E-06  2.590723E-06  2.614440E-06  1.568560E-03  1.568560E-03
   1.568560E-03
 Phonon frequencies in cm-1    :
-  5.677000E-01  5.685980E-01  5.738033E-01  3.442590E+02  3.442590E+02
-  3.442590E+02

  Phonon at Gamma, with non-analyticity in the
  direction (cartesian coordinates)  1.00000  0.00000  0.00000
 Phonon energies in Hartree :
   2.590670E-06  2.590723E-06  4.101029E-06  1.568560E-03  1.568560E-03
   1.729575E-03
 Phonon frequencies in cm-1    :
-  5.685864E-01  5.685980E-01  9.000719E-01  3.442590E+02  3.442590E+02
-  3.795979E+02

  Phonon at Gamma, with non-analyticity in the
  direction (cartesian coordinates)  0.00000  1.00000  0.00000
 Phonon energies in Hartree :
   2.586632E-06  2.614440E-06  4.088526E-06  1.568560E-03  1.568560E-03
   1.729575E-03
 Phonon frequencies in cm-1    :
-  5.677000E-01  5.738033E-01  8.973277E-01  3.442590E+02  3.442590E+02
-  3.795979E+02

  Phonon at Gamma, with non-analyticity in the
  direction (cartesian coordinates)  0.00000  0.00000  1.00000
 Phonon energies in Hartree :
   2.590723E-06  2.610339E-06  4.088540E-06  1.568560E-03  1.568560E-03
   1.729575E-03
 Phonon frequencies in cm-1    :
-  5.685980E-01  5.729031E-01  8.973308E-01  3.442590E+02  3.442590E+02
-  3.795979E+02
</pre>
<p>
The first few lines discard any effect of the homogeneous electric field,
while the next sections consider it along the three
cartesian coordinates.

<p>In the present material, the directionality of the electric field
has no influence. We note that there are still three acoustic mode, below
1cm^-1, while the optic modes have the correct degeneracies : two TO
modes at 344.3 cm^-1, and one LO mode at 379.6 cm^-1 .

<p>These values can be compared to experimental (361 cm^-1 , 402 cm^-1)
as well as theoretical (363 cm^-1 , 400 cm^-1) values (again the
<cite>Gianozzi et al</cite> paper mentioned above).
Most of the discrepancy comes from the
too low value of <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>.
Using ABINIT with <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=6 hartree
gives (358.8 cm^-1 , 389.8 cm^-1). The remaining of the discrepancy
may come partly from the pseudopotentials, that are particularly soft.

<p>The comparison of Born effective charges is also interesting.
After imposition of the neutrality sum rule, the Al Born
effective charge is 2.116 . The value from
<cite>Gianozzi et al</cite> is 2.17, the experimental value
is 2.18. Increasing <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>
to 6 hartree in ABINIT gives 2.168.

<p>For the dielectric tensor, it is more delicate.
The value from <cite>Gianozzi et al</cite> is 9.2, while the
experimental value is 8.2 . The agreement is not very good,
a fact that can be attributed to the LDA lack
of polarization-dependence
(<cite> X. Gonze, Ph. Ghosez and R. Godby, Phys. Rev. Lett. (1995)</cite>).
Still, the agreement of our calculation
with the theoretical result is not very good. With
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=3 hartree, we have
9.76. Changing it to 6 hartree gives 10.40 . A better k point
sampling (8x8x8), with
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=6 hartree,
reduces the value to 9.89 . Changing pseudopotentials
finally improves the agreement : with the much harder
13al.pspgth and 33as.psphgh pseudopotentials with
adequate <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=16 hartree
and 8x8x8 Monkhorst-Pack sampling, we reach a value of 9.37 .
This illustrates that the dielectric tensor is a much more sensitive
quantity than the others.

<p>
<hr>

<p><a name="6"></a><br>

<h3><b>6. Response-function calculation of phonon frequencies at non-zero q.</b> </h3>
<h3>&nbsp;</h3>

<p> The computation of phonon frequencies at non-zero q is actually
simpler than the one at Gamma. One must distinguish two cases.
Either the q wavevector connects k points that belong to the
same grid, or the wavevector q is general.
In any case, the computation within the response-function
formalism is more efficient than using the frozen-phonon
technique: the use of supercell is completely avoided.
For an explanation of this fact, see for example
section IV of <cite>X. Gonze, Phys. Rev. B55, 10337 (1997)</cite>.

<p>
You can copy the file ~abinit/tests/tutorespfn/Input/trf1_6.in in Work_rf1.
This is your input file. You should edit it.
As for the other RF tests, the changes
with respect to the file ~abinit/tests/tutorespfn/Input/trf1_1.in
are all gathered in the first part of this file.
The multi-dataset mode
is used, computing from scratch the ground-state properties,
then computing different dynamical matrices.
The run is about 2 minutes on a 2.8 GHz machine.
So, you would better leave your computer running, and either
read more of the ABINIT documentation
(why not the
<a href="../users/mrgddb_help.html" target="helpsimg">mrgddb_help</a>
 and the <a href="../users/anaddb_help.html" target="helpsimg">anaddb_help</a>),
or make a walk.

<p>The results of this simulation
can be compared to those provided in the
<cite>Gianozzi et al</cite> paper. The agreement is rather good,
despite the low cut-off energy, and different pseudopotentials.

<p>At X, they get 95cm^-1, 216cm^-1, 337cm^-1 and 393cm^-1,
while we get 92.5cm^-1, 204.6cm^-1, 313.9cm^-1 and 375.9cm^-1.
With <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=6 hartree,
we get 89.7cm^-1, 212.3cm^-1, 328.5cm^-1 and 385.8cm^-1.

<p>
At L, they get 71cm^-1, 212cm^-1, 352cm^-1 and 372cm^-1,
while we get 69.0cm^-1, 202.5cm^-1, 332.6cm^-1 and 352.3cm^-1.
With <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=6 hartree,
we get 68.1cm^-1, 208.5cm^-1, 346.7cm^-1 and 362.6cm^-1.

<p> At q=(0.1 0 0), we get
31.6cm^-1, 63.6cm^-1, 342.0cm^-1 and 379.7cm^-1.
The acoustic modes tends (nearly-)linearly to zero,
while the optic modes are close to their values at Gamma :
344.3 cm^-1 and 379.6 cm^-1.

<p>
<hr>

<p>This ABINIT tutorial is now finished. You are advised to
read now the
<a href="lesson_rf2.html"> second tutorial on response functions</a>

<script type="text/javascript" src="list_internal_links.js"> </script>

</body>
</html>