TABLE OF CONTENTS
- ABINIT/m_outxml
- m_outxml/out_geometry_xml
- m_outxml/out_resultsgs_xml
- m_outxml/outxml_finalise
- m_outxml/outxml_open
ABINIT/m_outxml [ Modules ]
NAME
m_outxml
FUNCTION
COPYRIGHT
Copyright (C) 1998-2022 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.
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_outxml 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 use m_dtset 28 29 use m_io_tools, only : open_file 30 use m_geometry, only : xcart2xred, xred2xcart 31 use m_results_gs , only : results_gs_type 32 33 implicit none 34 35 private 36 37 public :: outxml_open 38 public :: outxml_finalise 39 public :: out_resultsgs_XML 40 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
SOURCE
254 subroutine out_geometry_XML(dtset, level, natom, rprimd, xred) 255 256 !Arguments ------------------------------------ 257 !scalars 258 integer,intent(in) :: level,natom 259 type(dataset_type),intent(in) :: dtset 260 !arrays 261 real(dp),intent(in) :: rprimd(3,3) 262 real(dp),intent(inout) :: xred(3,natom) 263 264 !Local variables ------------------------- 265 character(len = 1), parameter :: rprimd_names(3) = (/ "x", "y", "z" /) 266 !scalars 267 integer :: i,j 268 character(len=128) :: spacer,value 269 !arrays 270 real(dp),allocatable :: xcart(:,:) 271 272 ! ************************************************************************* 273 274 !Compute the spacer to put before each markup 275 write(spacer, "(A,I0)") "A", 2 * level 276 write(ab_xml_out, "("//trim(spacer)//",A)") " ", '<geometry>' 277 !Compute the cartesian coordinates of atoms 278 ABI_MALLOC(xcart,(3, natom)) 279 call xred2xcart(natom, rprimd, xcart, xred) 280 !Ouput the rprimd matrix 281 write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", ' <rprimd' 282 do i = 1, 3, 1 283 do j = 1, 3, 1 284 write(value, "(es20.8)") rprimd(i, j) 285 write(ab_xml_out, "(A,A,I0,A,A,A)", advance = "NO") ' ', rprimd_names(i), j, '="', trim(value) ,'"' 286 end do 287 end do 288 write(ab_xml_out, "(A)") ' />' 289 !Output the atom position 290 do i = 1, natom, 1 291 write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", ' <position' 292 write(ab_xml_out, "(A,I0,A,I0,A)", advance = "NO") ' atom="a_', dtset%jdtset ,'_', i ,'"' 293 do j = 1, 3, 1 294 write(value, "(es20.8)") xcart(j, i) 295 write(ab_xml_out, "(A,A,A,A,A)", advance = "NO") ' ', rprimd_names(j), '="', trim(value) ,'"' 296 end do 297 write(ab_xml_out, "(A)") ' />' 298 end do 299 ABI_FREE(xcart) 300 write(ab_xml_out, "("//trim(spacer)//",A)") " ", '</geometry>' 301 302 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
SOURCE
146 subroutine out_resultsgs_XML(dtset, level, results_gs, usepaw) 147 148 !Arguments ------------------------------------ 149 !scalars 150 integer,intent(in) :: level,usepaw 151 type(dataset_type),intent(in) :: dtset 152 type(results_gs_type),intent(inout) :: results_gs 153 154 !Local variables ------------------------- 155 character(len = 1), parameter :: axes_names(3) = (/ "x", "y", "z" /) 156 !scalars 157 integer :: i,j 158 character(len=128) :: spacer,value 159 160 ! ************************************************************************* 161 162 !Compute the spacer to put before each markup 163 write(spacer, "(A,I0)") "A", 2 * level 164 165 !Begin with the energy part 166 write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", '<energy' 167 if (dtset%iscf < 10) then 168 write(ab_xml_out, "(A)", advance = "NO") ' type="direct"' 169 write(value, "(es20.8)") results_gs%energies%e_kinetic 170 write(ab_xml_out, "(A,A,A)", advance = "NO") ' kinetic="', trim(value) ,'"' 171 write(value, "(es20.8)") results_gs%energies%e_localpsp 172 write(ab_xml_out, "(A,A,A)", advance = "NO") ' local="', trim(value) ,'"' 173 write(value, "(es20.8)") results_gs%energies%e_nlpsp_vfock 174 write(ab_xml_out, "(A,A,A)", advance = "NO") ' non-local="', trim(value) ,'"' 175 if (usepaw == 1) then 176 write(value, "(es20.8)") results_gs%energies%e_paw 177 write(ab_xml_out, "(A,A,A)", advance = "NO") ' paw="', trim(value) ,'"' 178 end if 179 else 180 write(ab_xml_out, "(A)", advance = "NO") ' type="double-counting"' 181 write(value, "(es20.8)") results_gs%energies%e_eigenvalues 182 write(ab_xml_out, "(A,A,A)", advance = "NO") ' eigen-values="', trim(value) ,'"' 183 write(value, "(es20.8)") results_gs%energies%e_xcdc 184 write(ab_xml_out, "(A,A,A)", advance = "NO") ' xcdc="', trim(value) ,'"' 185 if (usepaw == 1) then 186 write(value, "(es20.8)") results_gs%energies%e_pawdc 187 write(ab_xml_out, "(A,A,A)", advance = "NO") ' pawdc="', trim(value) ,'"' 188 end if 189 end if 190 if (dtset%berryopt == 4 .or. dtset%berryopt == 6 .or. dtset%berryopt == 7 .or. & 191 & dtset%berryopt ==14 .or. dtset%berryopt ==16 .or. dtset%berryopt ==17) then 192 write(value, "(es20.8)") results_gs%energies%e_elecfield 193 write(ab_xml_out, "(A,A,A)", advance = "NO") ' electric-field="', trim(value) ,'"' 194 end if 195 if(dtset%occopt >= 3 .and. dtset%occopt <= 8) then 196 write(value, "(es20.8)") results_gs%energies%e_entropy 197 write(ab_xml_out, "(A,A,A)", advance = "NO") ' entropy="', trim(value) ,'"' 198 end if 199 write(value, "(es20.8)") results_gs%energies%e_ewald 200 if (dtset%icoulomb == 0) then 201 write(ab_xml_out, "(A,A,A)", advance = "NO") ' ewald="', trim(value) ,'"' 202 else 203 write(ab_xml_out, "(A,A,A)", advance = "NO") ' ion-ion="', trim(value) ,'"' 204 end if 205 write(value, "(es20.8)") results_gs%energies%e_chempot 206 write(ab_xml_out, "(A,A,A)", advance = "NO") ' chempot="', trim(value) ,'"' 207 write(value, "(es20.8)") results_gs%energies%e_hartree 208 write(ab_xml_out, "(A,A,A)", advance = "NO") ' hartree="', trim(value) ,'"' 209 write(value, "(es20.8)") results_gs%energies%e_corepsp 210 write(ab_xml_out, "(A,A,A)", advance = "NO") ' core="', trim(value) ,'"' 211 write(value, "(es20.8)") results_gs%energies%e_xc 212 write(ab_xml_out, "(A,A,A)", advance = "NO") ' xc="', trim(value) ,'"' 213 write(value, "(es20.8)") results_gs%etotal 214 write(ab_xml_out, "(A,A,A)", advance = "NO") ' total="', trim(value) ,'"' 215 write(ab_xml_out, "(A)") ' />' 216 217 218 !finish with the forces part 219 if (dtset%optforces == 1) then 220 write(ab_xml_out, "("//trim(spacer)//",A)") " ", '<forces>' 221 do i = 1, dtset%natom, 1 222 write(ab_xml_out, "("//trim(spacer)//",A)", advance = "NO") " ", ' <force' 223 write(ab_xml_out, "(A,I0,A,I0,A)", advance = "NO") ' atom="a_', dtset%jdtset ,'_', i ,'"' 224 do j = 1, 3, 1 225 write(value, "(es20.8)") results_gs%fcart(j, i) 226 write(ab_xml_out, "(A,A,A,A,A)", advance = "NO") ' ', axes_names(j), '="', trim(value) ,'"' 227 end do 228 write(ab_xml_out, "(A)") ' />' 229 end do 230 write(ab_xml_out, "("//trim(spacer)//",A)") " ", '</forces>' 231 end if 232 233 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.
SOURCE
101 subroutine outxml_finalise(tsec, values) 102 103 !Arguments ------------------------------- 104 integer, intent(in) :: values(8) 105 real(dp), intent(in) :: tsec(2) 106 !Local variables ------------------------- 107 character(len=500) :: message 108 109 ! ********************************************************************* 110 111 write(ab_xml_out, "(A)") ' <!-- timeInfo markup : cpu and wall attributes are given in seconds. -->' 112 write(ab_xml_out, "(A,I0,A,I0,A,I0,A)", advance = "NO") & 113 & ' <timeInfo day="', values(3) , & 114 & '" month="', values(2) ,'" year="', values(1) ,'" ' 115 write(message, *) tsec(1) 116 write(ab_xml_out, "(A,A,A)", advance = "NO") 'cpu="', trim(message) ,'"' 117 write(message, *) tsec(2) 118 write(ab_xml_out, "(A,A,A)") ' wall="', trim(message) ,'" />' 119 write(ab_xml_out, "(A)") "</abinitRun>" 120 !ABINIT has been compiled with XML output support, then we close the channel of the XML output file. 121 close(unit = ab_xml_out) 122 123 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.
SOURCE
59 subroutine outxml_open(filename) 60 61 !Arguments ------------------------------- 62 character(len = *), intent(in) :: filename 63 !Local variables ------------------------- 64 character(len=500) :: msg 65 66 ! ********************************************************************* 67 68 !ABINIT has been compiled with XML output support, then we open the 69 !channel of the XML output file. 70 if (open_file(trim(filename)//"_LOG.xml", msg, unit=ab_xml_out, form="formatted", action="write") /= 0) then 71 ABI_ERROR(msg) 72 end if 73 74 write(ab_xml_out, "(A)") '<?xml version="1.0" encoding="utf-8"?>' 75 !DTD declaration : root element is "abinitRun", and DTD file is not public 76 !but given in the distribution by the file abinitRun.dtd. 77 write(ab_xml_out, "(A)") '<!DOCTYPE abinitRun SYSTEM "extras/post_processing/abinitRun.dtd">' 78 !Creating the root markup. 79 write(ab_xml_out, "(A)") '<abinitRun>' 80 81 end subroutine outxml_open