TABLE OF CONTENTS


ABINIT/m_outxml [ Modules ]

[ Top ] [ Modules ]

NAME

 m_outxml

FUNCTION

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt.
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt.

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_outxml
28 
29  use defs_basis
30  use defs_abitypes
31  use m_abicore
32  use m_errors
33 
34  use m_io_tools,    only : open_file
35  use m_geometry,    only : xcart2xred, xred2xcart
36  use m_results_gs , only : results_gs_type
37 
38  implicit none
39 
40  private
41 
42  public :: outxml_open
43  public :: outxml_finalise
44  public :: out_resultsgs_XML
45  public :: out_geometry_XML

m_outxml/out_geometry_xml [ Functions ]

[ Top ] [ m_outxml ] [ Functions ]

NAME

 out_geometry_xml

FUNCTION

 Output in the XML file, the box size and the atomic position.
 (see extras/post_processing/abinitRun.dtd)

INPUTS

  dtset <type(dataset_type)>=all input variables for this dataset
  level=use for indentation of the XML, two spaces are added by level.
  natom=number of atoms.
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  xred(3,natom)=reduced dimensionless atomic coordinates

PARENTS

      scfcv

CHILDREN

      xred2xcart

SOURCE

310 subroutine out_geometry_XML(dtset, level, natom, rprimd, xred)
311 
312 
313 !This section has been created automatically by the script Abilint (TD).
314 !Do not modify the following lines by hand.
315 #undef ABI_FUNC
316 #define ABI_FUNC 'out_geometry_XML'
317 !End of the abilint section
318 
319  implicit none
320 
321 !Arguments ------------------------------------
322 !scalars
323  integer,intent(in) :: level,natom
324  type(dataset_type),intent(in) :: dtset
325 !arrays
326  real(dp),intent(in) :: rprimd(3,3)
327  real(dp),intent(inout) :: xred(3,natom)
328 
329 !Local variables -------------------------
330   character(len = 1), parameter :: rprimd_names(3) = (/ "x", "y", "z" /)
331 !scalars
332  integer :: i,j
333  character(len=128) :: spacer,value
334 !arrays
335  real(dp),allocatable :: xcart(:,:)
336 
337 ! *************************************************************************
338 
339 !Compute the spacer to put before each markup
340  write(spacer, "(A,I0)") "A", 2 * level
341  write(ab_xml_out, "("//trim(spacer)//",A)") " ", '<geometry>'
342 !Compute the cartesian coordinates of atoms
343  ABI_ALLOCATE(xcart,(3, natom))
344  call xred2xcart(natom, rprimd, xcart, xred)
345 !Ouput the rprimd matrix
346  write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", '  <rprimd'
347  do i = 1, 3, 1
348    do j = 1, 3, 1
349      write(value, "(es20.8)") rprimd(i, j)
350      write(ab_xml_out, "(A,A,I0,A,A,A)", advance = "NO") ' ', rprimd_names(i), j, '="', trim(value) ,'"'
351    end do
352  end do
353  write(ab_xml_out, "(A)") ' />'
354 !Output the atom position
355  do i = 1, natom, 1
356    write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", '  <position'
357    write(ab_xml_out, "(A,I0,A,I0,A)", advance = "NO") ' atom="a_', dtset%jdtset ,'_', i ,'"'
358    do j = 1, 3, 1
359      write(value, "(es20.8)") xcart(j, i)
360      write(ab_xml_out, "(A,A,A,A,A)", advance = "NO") ' ', rprimd_names(j), '="', trim(value) ,'"'
361    end do
362    write(ab_xml_out, "(A)") ' />'
363  end do
364  ABI_DEALLOCATE(xcart)
365  write(ab_xml_out, "("//trim(spacer)//",A)") " ", '</geometry>'
366 
367 end subroutine out_geometry_XML

m_outxml/out_resultsgs_xml [ Functions ]

[ Top ] [ m_outxml ] [ Functions ]

NAME

 out_resultsgs_xml

FUNCTION

 Output in the XML file, the decomposition of the energy and
 the forces after a scfcv loop.
 (see extras/post_processing/abinitRun.dtd)
 Warning : this method is not thread safe and should be called only by one thread.

INPUTS

  dtset <type(dataset_type)>=all input variables for this dataset
  level=use for indentation of the XML, two spaces are added by level.
  results_gs <type(results_gs_type)>=results (energy and its components,
   forces and its components, the stress tensor) of a ground-state computation.
  usepaw= 0 for non paw calculation; =1 for paw calculation

PARENTS

      scfcv

CHILDREN

      xred2xcart

SOURCE

187 subroutine out_resultsgs_XML(dtset, level, results_gs, usepaw)
188 
189 
190 !This section has been created automatically by the script Abilint (TD).
191 !Do not modify the following lines by hand.
192 #undef ABI_FUNC
193 #define ABI_FUNC 'out_resultsgs_XML'
194 !End of the abilint section
195 
196  implicit none
197 
198 !Arguments ------------------------------------
199 !scalars
200  integer,intent(in) :: level,usepaw
201  type(dataset_type),intent(in) :: dtset
202  type(results_gs_type),intent(inout) :: results_gs
203 
204 !Local variables -------------------------
205   character(len = 1), parameter :: axes_names(3) = (/ "x", "y", "z" /)
206 !scalars
207  integer :: i,j
208  character(len=128) :: spacer,value
209 
210 ! *************************************************************************
211 
212 !Compute the spacer to put before each markup
213  write(spacer, "(A,I0)") "A", 2 * level
214 
215 !Begin with the energy part
216  write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", '<energy'
217  if (dtset%iscf < 10) then
218    write(ab_xml_out, "(A)", advance = "NO") ' type="direct"'
219    write(value, "(es20.8)") results_gs%energies%e_kinetic
220    write(ab_xml_out, "(A,A,A)", advance = "NO") ' kinetic="', trim(value) ,'"'
221    write(value, "(es20.8)") results_gs%energies%e_localpsp
222    write(ab_xml_out, "(A,A,A)", advance = "NO") ' local="', trim(value) ,'"'
223    write(value, "(es20.8)") results_gs%energies%e_nonlocalpsp
224    write(ab_xml_out, "(A,A,A)", advance = "NO") ' non-local="', trim(value) ,'"'
225    if (usepaw == 1) then
226      write(value, "(es20.8)") results_gs%energies%e_paw
227      write(ab_xml_out, "(A,A,A)", advance = "NO") ' paw="', trim(value) ,'"'
228    end if
229  else
230    write(ab_xml_out, "(A)", advance = "NO") ' type="double-counting"'
231    write(value, "(es20.8)") results_gs%energies%e_eigenvalues
232    write(ab_xml_out, "(A,A,A)", advance = "NO") ' eigen-values="', trim(value) ,'"'
233    write(value, "(es20.8)") results_gs%energies%e_xcdc
234    write(ab_xml_out, "(A,A,A)", advance = "NO") ' xcdc="', trim(value) ,'"'
235    if (usepaw == 1) then
236      write(value, "(es20.8)") results_gs%energies%e_pawdc
237      write(ab_xml_out, "(A,A,A)", advance = "NO") ' pawdc="', trim(value) ,'"'
238    end if
239  end if
240  if (dtset%berryopt == 4 .or. dtset%berryopt == 6 .or. dtset%berryopt == 7 .or. &
241 & dtset%berryopt ==14 .or. dtset%berryopt ==16 .or. dtset%berryopt ==17) then
242    write(value, "(es20.8)") results_gs%energies%e_elecfield
243    write(ab_xml_out, "(A,A,A)", advance = "NO") ' electric-field="', trim(value) ,'"'
244  end if
245  if(dtset%occopt >= 3 .and. dtset%occopt <= 8) then
246    write(value, "(es20.8)") results_gs%energies%e_entropy
247    write(ab_xml_out, "(A,A,A)", advance = "NO") ' entropy="', trim(value) ,'"'
248  end if
249  write(value, "(es20.8)") results_gs%energies%e_ewald
250  if (dtset%icoulomb == 0) then
251    write(ab_xml_out, "(A,A,A)", advance = "NO") ' ewald="', trim(value) ,'"'
252  else
253    write(ab_xml_out, "(A,A,A)", advance = "NO") ' ion-ion="', trim(value) ,'"'
254  end if
255  write(value, "(es20.8)") results_gs%energies%e_chempot
256  write(ab_xml_out, "(A,A,A)", advance = "NO") ' chempot="', trim(value) ,'"'
257  write(value, "(es20.8)") results_gs%energies%e_hartree
258  write(ab_xml_out, "(A,A,A)", advance = "NO") ' hartree="', trim(value) ,'"'
259  write(value, "(es20.8)") results_gs%energies%e_corepsp
260  write(ab_xml_out, "(A,A,A)", advance = "NO") ' core="', trim(value) ,'"'
261  write(value, "(es20.8)") results_gs%energies%e_xc
262  write(ab_xml_out, "(A,A,A)", advance = "NO") ' xc="', trim(value) ,'"'
263  write(value, "(es20.8)") results_gs%etotal
264  write(ab_xml_out, "(A,A,A)", advance = "NO") ' total="', trim(value) ,'"'
265  write(ab_xml_out, "(A)") ' />'
266 
267 
268 !finish with the forces part
269  if (dtset%optforces == 1) then
270    write(ab_xml_out, "("//trim(spacer)//",A)") " ", '<forces>'
271    do i = 1, dtset%natom, 1
272      write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", '  <force'
273      write(ab_xml_out, "(A,I0,A,I0,A)", advance = "NO") ' atom="a_', dtset%jdtset ,'_', i ,'"'
274      do j = 1, 3, 1
275        write(value, "(es20.8)") results_gs%fcart(j, i)
276        write(ab_xml_out, "(A,A,A,A,A)", advance = "NO") ' ', axes_names(j), '="', trim(value) ,'"'
277      end do
278      write(ab_xml_out, "(A)") ' />'
279    end do
280    write(ab_xml_out, "("//trim(spacer)//",A)") " ", '</forces>'
281  end if
282 
283 end subroutine out_resultsgs_XML

m_outxml/outxml_finalise [ Functions ]

[ Top ] [ m_outxml ] [ Functions ]

NAME

 outxml_finalise

FUNCTION

 Close the XML log file, and write the timing information.
 (see extras/post_processing/abinitRun.dtd)
 Warning : this method is not thread safe and should be called
 only by one thread.

INPUTS

  tsec=the cpu time and the wall time in seconds.
  values=the date values returned by date_and_time() intrinsic Fortran routine.

PARENTS

      abinit

CHILDREN

      xred2xcart

SOURCE

127 subroutine outxml_finalise(tsec, values)
128 
129 
130 !This section has been created automatically by the script Abilint (TD).
131 !Do not modify the following lines by hand.
132 #undef ABI_FUNC
133 #define ABI_FUNC 'outxml_finalise'
134 !End of the abilint section
135 
136   implicit none
137 
138 !Arguments -------------------------------
139   integer, intent(in) :: values(8)
140   real(dp), intent(in) :: tsec(2)
141 !Local variables -------------------------
142   character(len=500) :: message
143 
144 ! *********************************************************************
145 
146  write(ab_xml_out, "(A)") '  <!-- timeInfo markup : cpu and wall attributes are given in seconds. -->'
147  write(ab_xml_out, "(A,I0,A,I0,A,I0,A)", advance = "NO") &
148 & '  <timeInfo day="', values(3) , &
149 & '" month="', values(2) ,'" year="', values(1) ,'" '
150  write(message, *) tsec(1)
151  write(ab_xml_out, "(A,A,A)", advance = "NO") 'cpu="', trim(message) ,'"'
152  write(message, *) tsec(2)
153  write(ab_xml_out, "(A,A,A)") ' wall="', trim(message) ,'" />'
154  write(ab_xml_out, "(A)") "</abinitRun>"
155 !ABINIT has been compiled with XML output support, then we close the channel of the XML output file.
156  close(unit = ab_xml_out)
157 
158 end subroutine outxml_finalise

m_outxml/outxml_open [ Functions ]

[ Top ] [ m_outxml ] [ Functions ]

NAME

 outxml_open

FUNCTION

 Open the XML log file, and write the header inside.
 (see extras/post_processing/abinitRun.dtd)
 Warning: this method is not thread safe and should be called only by one thread.

INPUTS

  filename=the name of the file to write to.

PARENTS

      abinit

CHILDREN

      xred2xcart

SOURCE

 70 subroutine outxml_open(filename)
 71 
 72 
 73 !This section has been created automatically by the script Abilint (TD).
 74 !Do not modify the following lines by hand.
 75 #undef ABI_FUNC
 76 #define ABI_FUNC 'outxml_open'
 77 !End of the abilint section
 78 
 79   implicit none
 80 
 81 !Arguments -------------------------------
 82   character(len = *), intent(in) :: filename
 83 !Local variables -------------------------
 84  character(len=500) :: msg
 85 
 86 ! *********************************************************************
 87 
 88 !ABINIT has been compiled with XML output support, then we open the
 89 !channel of the XML output file.
 90  if (open_file(trim(filename)//"_LOG.xml", msg, unit=ab_xml_out, form="formatted", action="write") /= 0) then
 91    MSG_ERROR(msg)
 92  end if
 93 
 94  write(ab_xml_out, "(A)") '<?xml version="1.0" encoding="utf-8"?>'
 95 !DTD declaration : root element is "abinitRun", and DTD file is not public
 96 !but given in the distribution by the file abinitRun.dtd.
 97  write(ab_xml_out, "(A)") '<!DOCTYPE abinitRun SYSTEM "extras/post_processing/abinitRun.dtd">'
 98 !Creating the root markup.
 99  write(ab_xml_out, "(A)") '<abinitRun>'
100 
101 end subroutine outxml_open