File: lesson_ffield.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 (612 lines) | stat: -rw-r--r-- 28,401 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
<html>
<head>
<title>Tutorial on polarization and finite electric field.</title>
</head>
<body bgcolor="#ffffff">
<hr>
<h1>ABINIT, lesson on polarization and finite electric field</h1>
<h2>
AlAs : polarization, and responses to finite electric fields</h2>
<hr>

<p>This lesson aims at showing how to get the following physical properties, for an insulator :
<ul>
 <li> The polarization.</li>
 <li> The Born effective charge (by finite differences of polarization)</li>
 <li> The Born effective charge (by finite differences of forces)</li>
 <li> The dielectric constant </li>
 <li> The proper piezoelectric tensor (clamped and relaxed ions) </li>
</ul>

<p>
The case of the linear responses (Born effective charge, dielectric constant, piezoelectric tensor)
is treated independently in other tutorials
(<a href="lesson_rf1.html">lesson Response-Function 1</a>,
 <a href="lesson_elastic.html">lesson on Elastic properties</a>),
using Density-Functional Perturbation Theory.
You will learn here how to get these quantities using
the finite difference techniques within ABINIT.
To that end, we will describe how to compute the polarization, in the Berry phase formulation,
and how to perform finite electric field calculations.

<p>This lesson should take about 1 hour and 30 minutes. </p>

<p>The basic theory for Berry phase computation of the polarization
was proposed by R. D. King-Smith and D. Vanderbilt, <i>Phys. Rev. B</i> <b>47</b>, 1651 (1993).
The longer (excellent) paper D. Vanderbilt and R. D. King-Smith, <i>Phys. Rev. B</i> <b>48</b>, 
4442 (1993) clarifies many aspects of this theory (especially in view of application to AlAs, as in
this tutorial).  One might benefit also from a reading of 
R. Resta, <i>Rev. Mod. Physics</i> <b>66</b>, 899 (1994).

<p>
In order to gain the theoretical background needed to perform
a calculation with a finite electric field, you should consider reading the following papers:
<ul>
 <li>I. Souza, J. Iniguez and D. Vanderbilt, <i>Phys. Rev. Lett</i> <b>89</b>, 117602 (2002)</li>
 <li>R. Nunes and X. Gonze, <i>Phys. Rev. B</i> <b>63</b>, 155107 (2001)</li>
 <li>M. Veithen, <a href="http://ftp.abinit.org/docs/PhD-M.Veithen.pdf">PhD thesis</a></li>
</ul>

<p> Finally, the extension to the PAW formalism specifically in ABINIT is discussed in
Gonze <i>et al.,</i> <i>Comp. Phys. Commun.</i> <b> 180</b>, 2582 (2009) and
Zwanziger <i>et al.,</i> <i>Comp. Mater. Sci.</i> <b> 58</b>, 113 (2012).


<h5>Copyright (C) 2005-2014 ABINIT group (PhG, MVeithen, XG)
<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>
 
<h3><b>Contents of lesson on polarization and finite electric field :</b></h3>

<ul>
  <li><a href="lesson_ffield.html#0">Part-0. Ground-state properties of AlAs and general parameters</a>
  <li><a href="lesson_ffield.html#1">Part-1. Finite difference of polarization from Berry phase calculations in zero field</a>
  <li><a href="lesson_ffield.html#2">Part-2. Finite electric field calculations</a>
</ul>
<hr>

<a name="0">&nbsp;</a>
<h3><b>
Part-0. Ground-state properties of AlAs and general parameters
</b></h3>

<p><i>Before beginning, you might consider working in a different
subdirectory, as for the other lessons.
Why not create "Work-ffield" in ~abinit/tests/tutorespfn/Input ?
</i>

<p>
In this tutorial we will assume that the ground-state properties of AlAs have been previously obtained,
and that the corresponding convergence studies have been done.
We will adopt the following set of generic parameters :

<pre>
   acell               10.53
   ixc                 3
   ecut                2.8     (results with ecut = 5 are also reported
                               in the discussion)
   ecutsm              0.5
   dilatmx             1.05
   nband               4 (=number of occupied bands)
   ngkpt               6 6 6
   nshiftk             4
   shiftk   0.5 0.5 0.5
            0.5 0.0 0.0
            0.0 0.5 0.0
            0.0 0.0 0.5

   pseudopotentials    13al.pspnc
                       33as.pspnc
</pre>

<p>
In principle, the <a href="../input_variables/varbas.html#acell" target="kwimg">acell</a> to be used should be the one corresponding to the optimized structure
at the <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>, and
<a href="../input_variables/varbas.html#ngkpt" target="kwimg">ngkpt</a> combined with
<a href="../input_variables/varbas.html#nshiftk" target="kwimg">nshiftk</a> and
<a href="../input_variables/varbas.html#shiftk" target="kwimg">shiftk</a>, chosen for the calculations.
Unfortunately,  for the purpose of this tutorial, in order to limit the duration of the runs,
we have to work at an unusually low cutoff of 2.8 Ha for which the optimized lattice constant
is equal to 7.45*0.707*2 Bohr=10.53 Bohr (instead of the converged value of 10.64 Bohr).
For comparison, results with <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5 are also reported and, in that case,
were obtained at the optimized lattice constant of 10.64 bohr.
For those who would like to try later,
convergence tests and structural optimizations can be done using the file ~abinit/tests/tutorespfn/Input/tnlo_1.in.
Before going further, you might refresh your memory concerning the other variables :
<a href="../input_variables/varbas.html#ixc" target="kwimg">ixc</a>,
<a href="../input_variables/varrlx.html#ecutsm" target="kwimg">ecutsm</a>,
<a href="../input_variables/varrlx.html#dilatmx" target="kwimg">dilatmx</a>.
<p>

<hr>

<a name="1">&nbsp;</a>
<h3><b>
Part-1. Berry phase calculation in zero field
</b></h3>

<p> In this section, you will learn how to perform a Berry phase calculation
of the polarization.  As a practical problem we will try to compute the
Born effective charges from finite difference of the polarization (under
finite atomic displacements), for AlAs.

<p> You can now copy the file ~abinit/tests/tutorespfn/Input/tffield_1.in
and ~abinit/tests/tutorespfn/Input/tffield_x.files in Work-ffield, and
modify the latter as usual (for example, edit it so that the various file names it
contains refer to tffield_1 rather than tffield_x). 
<p>
  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>
  Because of the use of
  two types of atoms, have also a look at the following input variables present, in the "input" file:
<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>
You can start the calculation. It should take
90 seconds on a PC 3GHz.  Then, examine the tffield_1.in file. It is made
of three datasets corresponding to the reference optimized structure
(tau=0) and to structure with the Al atom displaced from 0.01 bohr
right and left (referred to as tau = +0.01 and tau = -0.01).  This is
typically the amplitude of atomic displacement to be considered in this
kind of computations.  Notice also that the displacements are given using
<a href="../input_variables/varbas.html#xcart" target="kwimg">xcart</a>,
that is, explicitly in Cartesian directions, rather than the primitive
cell axes (using <a href="../input_variables/varbas.html#xred"
target="kwimg">xred</a>).  This makes the correspondence with the
polarization output in Cartesian directions much simpler to understand.

<p> There are two implementations of the Berry phase
within ABINIT.  One corresponds to positive values of <a
href="../input_variables/varff.html#berryopt" target="kwimg">berryopt</a> and
was implemented by Na Sai.  The other one corresponds to negative values of
<a href="../input_variables/varff.html#berryopt" target="kwimg">berryopt</a>
was implemented by Marek Veithen.  Both are suitable to compute the
polarization.  Here we will focus on the implementation of Marek Veithen
for two reasons.  First, the results are directly provided in Cartesian
coordinates at the end of the run (while the implementation of Na Sai
reports them in reduced coordinates).  Second the implementation of Marek
Veithen is the one to be used for the finite electric field calculation
as described in the next Section.  Finally, note also that Veithen's
implementation works with <a href="../input_variables/varbas.html#kptopt"
target="kwimg">kptopt</a> = 1 or 2 while Na Sai implementation is restricted
to <a href="../input_variables/varbas.html#kptopt" target="kwimg">kptopt</a>
= 2, which is less convenient.

<p>
The file is the one of a usual self-consistent calculation.
On top of the usual variables, for the Berry phase calculation we simply need to define
<a href="../input_variables/varff.html#berryopt" target="kwimg">berryopt</a>
and <a href="../input_variables/varrf.html#rfdir" target="kwimg">rfdir</a>:
<pre>
        berryopt      -1
        rfdir          1 1 1
</pre>
<p>
Make the run, then open the output file and look for the occurrence "Berry".
The output reports values of the Berry phase for individual k-point strings.

<pre>
 Computing the polarization (Berry phase) for reciprocal vector:
  0.16667  0.00000  0.00000 (in reduced coordinates)
 -0.01583  0.01583  0.01583 (in cartesian coordinates - atomic units)
 Number of strings:   144
 Number of k points in string:    6

 Summary of the results
 Electronic Berry phase     2.025853583E-03
            Ionic phase    -7.500000000E-01
            Total phase    -7.479741464E-01
    Remapping in [-1,1]    -7.479741464E-01

           Polarization    -1.557862805E-02 (a.u. of charge)/bohr^2
           Polarization    -8.913274881E-01 C/m^2
</pre>

<p> The "Remapping in [-1,1]"  is there to avoid the quantum
of polarization.  As discussed in H. Djani, E. Bousquet, A. Kellou,
Ph. Ghosez, <I>Phys. Rev. B</I> <B>86</B>, 054107 (2012), the indeterminacy
of the quantum phase, directly related to the quantum of polarization,
can lead to spurious effects (see Fig. 2 of the above-mentioned paper).
By remapping on the [-1,1] interval, any indeterminacy is removed. However,
removing such a quantum of polarization between two calculations might
give the false impression that one is on the same polarization branch in
the two calculations, while actually the branch is made different by this
remapping.  Cross-checking the polarization results by computing the Born
effective charge, further multiplied by the displacements between the two
geometries is an excellent way to estimate the amplitude of the polarization.

<p> Other subtleties of Berry phases, explained in D. Vanderbilt and
R. D. King-Smith, <i>Phys. Rev. B</i> <b>48</b>, 4442 (1993), also apply.
First, note that neither the electronic Berry phase nor the ionic phase
vanish in this highly symmetric case, contrary to intuition. Even though AlAs
does not have inversion symmetry, it does have tetrahedral symmetry, which
would be enough to make an ordinary vector vanish.  But a lattice-valued
vector does not have to vanish: the lattice just has to transform into itself
under the tetrahedral point group.  The ionic phase corresponds actually to a
lattice-valued vector (-3/4 -3/4 -3/4).  Concerning the electronic phase, it
does not exactly vanish, unless the sampling of k points becomes continuous.

If you go further in the file you will find
the final results in cartesian coordinates. You can collect them for the different values of tau.
<p>
tau = 0
<pre>
 Polarization in cartesian coordinates (a.u.):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.730821479E-04   0.730821479E-04   0.730821479E-04
     Ionic:                        -0.270560574E-01  -0.270560574E-01  -0.270560574E-01
     Total:                        -0.269829753E-01  -0.269829753E-01  -0.269829753E-01

 Polarization in cartesian coordinates (C/m^2):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.418137766E-02   0.418137766E-02   0.418137766E-02
     Ionic:                        -0.154800587E+01  -0.154800587E+01  -0.154800587E+01
     Total:                        -0.154382450E+01  -0.154382450E+01  -0.154382450E+01
</pre>
<p>
tau = +0.01
<pre>
 Polarization in cartesian coordinates (a.u.):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.410028779E-04   0.730922946E-04   0.730922946E-04
     Ionic:                        -0.269532804E-01  -0.270560574E-01  -0.270560574E-01
     Total:                        -0.269122775E-01  -0.269829652E-01  -0.269829652E-01

 Polarization in cartesian coordinates (C/m^2):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.234596988E-02   0.418195820E-02   0.418195820E-02
     Ionic:                        -0.154212551E+01  -0.154800587E+01  -0.154800587E+01
     Total:                        -0.153977954E+01  -0.154382392E+01  -0.154382392E+01
</pre>
<p>
tau = -0.01
<pre>
 Polarization in cartesian coordinates (a.u.):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.105181709E-03   0.730922946E-04   0.730922946E-04
     Ionic:                        -0.271588345E-01  -0.270560574E-01  -0.270560574E-01
     Total:                        -0.270536528E-01  -0.269829652E-01  -0.269829652E-01

 Polarization in cartesian coordinates (C/m^2):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.601794637E-02   0.418195820E-02   0.418195820E-02
     Ionic:                        -0.155388624E+01  -0.154800587E+01  -0.154800587E+01
     Total:                        -0.154786829E+01  -0.154382392E+01  -0.154382392E+01
</pre>
<p>
From the previous data, we can extract the Born effective charge of Al.
Values to be used are those in a.u., in order to find the charge in electron unit.
It corresponds to (the volume of the primitive unit cell must be specified in Bohr too) :
<pre>
        Z* = \Omega_0  (P[tau=+0.01] - P[tau=-0.01]) / (2*tau)
           =  291.89   (-0.269122775E-01 - -0.270536528E-01) / 0.02
           = 2.06
</pre>
<br>
For comparison, the calculation using Density-Functional Perturbation Theory (DFPT)
 can be done by using the file 
~abinit/tests/tutorespfn/Input/tffield_2.in.
Actually, the file tffield_2.in not only leads to the computation of the
Born effective charges, but also the computation of the piezoelectric constants (see later).

<p>
You can review how to use DFPT thanks to the
<a href="lesson_rf1.html">lesson Response-function 1</a> and
<a href="lesson_rf2.html">lesson Response-function 2</a> tutorials.
For now, go ahead and run the input file, and have a look at the output file,
to identify the place where the Born effective charge is written (search for
the phrase "Effective charges" in the output file).
The value we get from DFPT is 2.06, in agreement with the
above-mentioned value of 2.06. This level of agreement is fortuitous for 
unconverged calculations, though
both methods (finite-difference and DFPT) will tend to the same value for better converged calculations.

<p>
The DDB generated by ~abinit/tests/tutorespfn/Input/tffield_2.in 
can be fed in anaddb, thanks to the ~abinit/tests/tutorespfn/Input/tffield_3.in
input file and the tffield_3.files file. Note that tffield_3.files is expecting the DDB to be
named tffield_2o_DS3_DDB, so either adjust tffield_x.files before running tffield_2.in to match this,
or else change the name of the DDB file after generation. In any case, the DFPT
calculation gives for the piezoelectric constants (as well as the elastic constants), as found in
tffield_3.out, the following:
<pre>
  Proper piezoelectric constants (clamped ion) (unit:c/m^2)
 ...
      -0.64879283      0.00000000      0.00000000
 ...
  Proper piezoelectric constants (relaxed ion) (unit:c/m^2)
 ...
      0.04574987      0.00000000      0.00000000
</pre>
Restarting case ffield_2 and ffield_3 with
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5 (three minutes on a PC 3 GHz),
one gets
the following values
<pre>
  Proper piezoelectric constants (clamped ion) (unit:c/m^2)
 ...
     -0.68377667      0.00000000      0.00000000
 ...
  Proper piezoelectric constants (relaxed ion) (unit:c/m^2)
 ...
      0.03384772      0.00000000      0.00000000
</pre>
The latter value, where the ionic relaxation largely suppresses the electronic piezoelectricity,
will be much more difficult to converge than the clamped ion result.

<p> Using the Berry phase approach it is also possible to compute
the piezoelectric constant from finite difference of polarization
with respect to strains.  This can be done considering clamped ions
or relaxed ions configurations.  For this purpose, have a look at the
files ~abinit/tests/tutorespfn/Input/tffield_4.in (clamped ions) and
~abinit/tests/tutorespfn/Input/tffield_5.in (relaxed ions).  Notice how
in the relaxed ion case, the input file includes
<a href="../input_variables/varrlx.html#ionmov" target="kwimg">ionmov</a>=2  
and 
<a href="../input_variables/varrlx.html#optcell" target="kwimg">optcell</a>=0,
in order to relax the ion positions at fixed cell geometry.  
These calculations should
give the following final results (obtained by taking finite difference
expressions of the strains for different electric fields)
<pre>
Proper piezoelectric constants(clamped ion)(Unit:c/m^2) = -0.6491
Proper piezoelectric constants(relaxed ion)(Unit:c/m^2) =  0.0430
</pre>
For example, the clamped ion piezoelectric constant was obtained from tffield_4.out :
<pre>
 Polarization in cartesian coordinates (C/m^2):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:       -0.230212263E-02   0.421671951E-02   0.421671986E-02
     Ionic:                        -0.154692152E+01  -0.155466099E+01  -0.155466099E+01
     Total:                        -0.154922364E+01  -0.155044427E+01  -0.155044427E+01
</pre>
and
<pre>
 Polarization in cartesian coordinates (C/m^2):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic berry phase:        0.106812992E-01   0.417465885E-02   0.417465885E-02
     Ionic:                        -0.154692152E+01  -0.153919172E+01  -0.153919172E+01
     Total:                        -0.153624022E+01  -0.153501707E+01  -0.153501707E+01
</pre>
The difference between -0.154922364E+01 (x-strain is +0.01) and -0.153624022E+01 (x-strain is -0.01)
gives the finite
difference -0.012983, that, divided by 0.02 (the change of strain) gives -0.6491, announced above.

<br>
Using to <a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5
gives
<pre>
Proper piezoelectric constants(clamped ion)(Unit:c/m^2) = -0.6841215
Proper piezoelectric constants(relaxed ion)(Unit:c/m^2) =  0.0313915
</pre>
in excellent agreement with the above-mentioned DFPT values.

<hr>

<a name="2">&nbsp;</a>
<h3><b>
Part-2. Finite electric field calculations
</b></h3>

<br>
In this section, you will learn how to perform a calculation with a finite electric field.
<p>
You can now copy the file ~abinit/tests/tutorespfn/Input/tffield_6.in in Work-ffield,
and modify the tffield_x.files accordingly.
You can start the run immediately, it will take one minute on a PC 3GHz.
It performs a finite field calculation at clamped atomic positions.
You can look at this input file to identify the features specific to the finite field calculation.
<p>
As general parameters, one has to specify
<a href="../input_variables/varbas.html#nband" target="kwimg">nband</a>,
<a href="../input_variables/vargs.html#nbdbuf" target="kwimg">nbdbuf</a>
and <a href="../input_variables/varbas.html#kptopt" target="kwimg">kptopt</a>:
<pre>
        nband          4
        nbdbuf         0
        kptopt         1
</pre>
<p>
As a first step (dataset 11), the code must perform a Berry phase calculation in zero field.
For that purpose, it is necessary to set the values of
<a href="../input_variables/varff.html#berryopt" target="kwimg">berryopt</a>
and <a href="../input_variables/varrf.html#rfdir" target="kwimg">rfdir</a>:
<pre>
        berryopt11     -1
        rfdir11        1 1 1
</pre>
<p>
<i>Important warning:</i> you cannot use berryopt to +1 to initiate a finite field calculation. You must
begin with berryopt -1.
<p>
After that, there are different steps corresponding to various value of the electric field,
thanks to <a href="../input_variables/varff.html#efield" target="kwimg">efield</a>.
For those steps, it is important to take care of the following parameters, shown here for dataset 21:
<pre>
        berryopt21     4
        efield21       0.0001  0.0001  0.0001
        getwfk21       11
</pre>
<p>
The electric field is applied here along the [111] direction.
It must be incremented step by step (it is not possible to go to high field directly).
At each step, the wavefunctions of the previous step must be used.
When the field gets larger, it is important to check that it does not significantly exceed
the critical field for Zener breakthrough
(the drop of potential over the Born-von Karmann supercell must be smaller than the gap).
In practice, the number of k-point must be enlarged to reach convergence.
However, at the same time, the critical field becomes smaller.
In practice, reasonable fields can still be reached for k-point grids providing
a reasonable degree of convergence. A compromise must however be found.
<p>
As these calculations are quite long, the input file has been limited to very small fields.
Three cases have been selected : E = 0, E = +0.0001 and E = -0.0001.
If you have time later, you can relax this constraint
and perform a more exhaustive calculations for a larger set of fields.
<p>
You can now start the run.
Various quantities can be extracted from the finite field calculation
at clamped ions using finite difference techniques:
the Born Effective Charge Z* can be extracted from the forces,
the optical dielectric constant can be deduced from the polarization,
and the clamped ion piezoelectric tensor can be deduced from the stress tensor.
As an illustration we will focus here on the computation of Z*.
<p>
You can look at your output file.
For each field, the file contains various quantities that can be collected.
For the present purpose, we can look at the evolution of the forces
with the field and extract the following results from the output file :
<p>
E=0
<pre>
cartesian forces (hartree/bohr) at end:
    1      0.00000000000000     0.00000000000000     0.00000000000000
    2      0.00000000000000     0.00000000000000     0.00000000000000
</pre>
<p>
E = +0.0001
<pre>
cartesian forces (hartree/bohr) at end:
    1      0.00020614538503     0.00020614538503     0.00020614538503
    2     -0.00020614538503    -0.00020614538503    -0.00020614538503
</pre>
<p>
E = -0.0001
<pre>
cartesian forces (hartree/bohr) at end:
    1     -0.00020650900243    -0.00020650900243    -0.00020650900243
    2      0.00020650900243     0.00020650900243     0.00020650900243
</pre>
<p>
In a finite electric field, the force on  atom A in direction i can be written as :
<pre>
F_Ai = Z*_A,ii * E + \Omega_0 dchi/dtau E^2
</pre>
<p>
The value for positive and negative fields above are nearly the same,
showing that the quadratic term is almost negligible.
This clearly appears in the figure below where the field dependence
of the force for a larger range of electric fields is plotted.
<p>
<img src="lesson_nlo/image001.gif">
<p>
We can therefore extract with a good accuracy the Born effective charge as :
<pre>
Z*_Al    = (F_Al[E=+0.0001] - F_Al[E=-0.0001])/(2*0.0001)
         = (0.00020614538503 + 0.00020650900243)/0.0002
         = 2.06
</pre>
<p>
This value is similar to the value reported from DFPT.
If you do calculations for more electric fields,
fitting them with the general expression of the force above (including the E^2 term), you can find
the dchi/dtau term.  At
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5 you should get dchi/dtau for Al = -0.06169.
This will be useful later for the
<a href="lesson_nlo.html">lesson on static Non-linear properties</a>.

<p>
Going back to the output file, you can also look at the evolution of the polarization with the field.
<p>
E=0
<pre>
Polarization in cartesian coordinates (a.u.):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic:  0.730821479E-04   0.730821479E-04   0.730821479E-04
     Ionic:      -0.270560574E-01  -0.270560574E-01  -0.270560574E-01
     Total:      -0.269829753E-01  -0.269829753E-01  -0.269829753E-01
</pre>
<p>
E = +0.0001
<pre>
Polarization in cartesian coordinates (a.u.):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic:  0.129869542E-03   0.129869542E-03   0.129869542E-03
     Ionic:      -0.270560574E-01  -0.270560574E-01  -0.270560574E-01
     Total:      -0.269261879E-01  -0.269261879E-01  -0.269261879E-01
</pre>
<p>
E = -0.0001
<pre>
Polarization in cartesian coordinates (a.u.):
 (the sum of the electronic and ionic Berry phase has been folded into [-1, 1])
     Electronic:  0.163575373E-04   0.163575372E-04   0.163575371E-04
     Ionic:      -0.270560574E-01  -0.270560574E-01  -0.270560574E-01
     Total:      -0.270396999E-01  -0.270396999E-01  -0.270396999E-01
</pre>
<p>
In a finite electric field, the polarization in terms of the linear and quadratic
susceptibilities as
<pre>
P = \chi * E + \chi(2)  E^2
</pre>
<p>
The change of polarization for positive and negative fields above are nearly
the same, showing again that the quadratic term is almost negligible.
This clearly appears in the figure below where
the field dependence of the polarization for a larger range of electric fields is plotted.
<p>
<img src="lesson_ffield/image004.gif">
<p>
We can therefore extract with a good accuracy the linear optical dielectric susceptibility :
<pre>
\chi(l)_11 = (P[E=+0.0001] - P[E=-0.0001])/(2*0.0001)
           = (-0.269261879E-01+ 0.270396999E-01)/0.0002
           = 0.56756
</pre>
<p>
The optical dielectric constant is :
<pre>
\epsilon_11 = 1+ 4 \pi \chi(1)_11
           = 8.13
</pre>
<p>
This value underestimates the value of 9.20 obtained by DFPT.
The agreement is not better at
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5. Typically, finite field calculations
converge with the density of the k point grid more slowly than DFPT calcuations.

<p>If you do calculations for more electric fields, fitting them with
the general expression of the force above (including the E^2 term)
you can find the non-linear optical susceptibilities \chi(2).
At
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5 you should get \chi(2) = 29.769 pm/V.
This result will be also be useful for the
<a href="lesson_nlo.html">lesson on static Non-linear properties</a>.
<p>
Looking at the evolution of the stress (see graphic below),
you should also be able to extract the piezoelectric constants.
You can try to do it as an exercise.
As the calculation here was at clamped ions,
you will get the clamped ion proper piezoelectric constants.
At
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5,
you should obtain -0.69088 C/m^2.
<p>
<img src="lesson_ffield/image007.gif">
<p>
Redoing the same kind of finite field calculation
when allowing the ions to relax
<!-- (using the input file ~abinit/tests/tutorespfn/Input/tffield_7.in)  -->
,
one can access the relaxed ion proper piezoelectric constant.
At
<a href="../input_variables/varbas.html#ecut" target="kwimg">ecut</a>=5 the result is -0.03259.

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

</body>
</html>