1! ---
2! Copyright (C) 1996-2016	The SIESTA group
3!  This file is distributed under the terms of the
4!  GNU General Public License: see COPYING in the top directory
5!  or http://www.gnu.org/copyleft/gpl.txt .
6! See Docs/Contributors.txt for a list of contributors.
7! ---
8      subroutine get_target_stress(tp,tstres)
9      use precision, only: dp
10      use parallel, only: ionode
11      use fdf
12      use units, only: kBar
13      use sys,   only: die
14
15      implicit none
16
17      real(dp), intent(in) :: tp
18      real(dp), intent(out) :: tstres(3,3)
19
20      integer  :: i, j
21      real(dp) :: sxx, syy, szz, sxy, sxz, syz
22      logical  :: tarstr
23
24      type(block_fdf)            :: bfdf
25      type(parsed_line), pointer :: pline
26
27C Look for target stress and read it if found, otherwise generate it --------
28          tarstr = fdf_block('MD.TargetStress',bfdf)
29
30          if (tarstr) then
31            if (ionode) then
32              write(6,'(/a,a)')
33     $             'Broyden_optim: Reading %block MD.TargetStress',
34     .             ' (units of MD.TargetPressure).'
35            endif
36            if (.not. fdf_bline(bfdf,pline))
37     $        call die('get_target_stress: ERROR in ' //
38     .                 'MD.TargetStress block')
39            sxx = fdf_bvalues(pline,1)
40            syy = fdf_bvalues(pline,2)
41            szz = fdf_bvalues(pline,3)
42            sxy = fdf_bvalues(pline,4)
43            sxz = fdf_bvalues(pline,5)
44            syz = fdf_bvalues(pline,6)
45            call fdf_bclose(bfdf)
46
47            tstres(1,1) = - sxx * tp
48            tstres(2,2) = - syy * tp
49            tstres(3,3) = - szz * tp
50            tstres(1,2) = - sxy * tp
51            tstres(2,1) = - sxy * tp
52            tstres(1,3) = - sxz * tp
53            tstres(3,1) = - sxz * tp
54            tstres(2,3) = - syz * tp
55            tstres(3,2) = - syz * tp
56          else
57            if (ionode) then
58              write(6,'(/a,a)')
59     $              'Broyden_optim: No target stress found, ',
60     .              'assuming hydrostatic MD.TargetPressure.'
61            endif
62            do i= 1, 3
63              do j= 1, 3
64                tstres(i,j) = 0._dp
65              enddo
66              tstres(i,i) = - tp
67            enddo
68          endif
69
70C Write target stress down --------------------------------------------------
71          if (ionode) then
72            write(6,"(/a)") 'Broyden_optim: Target stress (kBar)'
73            do i = 1, 3
74              write(6,"(a,2x,3f12.3)")
75     .            'get_target_stress:',
76     .              tstres(i,1)/kBar, tstres(i,2)/kBar,
77     .            tstres(i,3)/kBar
78            enddo
79          endif
80
81          end subroutine get_target_stress
82