TABLE OF CONTENTS


ABINIT/m_outvar_o_z [ Modules ]

[ Top ] [ Modules ]

NAME

  m_outvar_o_z

FUNCTION

COPYRIGHT

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

SOURCE

16 #if defined HAVE_CONFIG_H
17 #include "config.h"
18 #endif
19 
20 #include "abi_common.h"
21 
22 module m_outvar_o_z
23 
24  use defs_basis
25 
26  use m_errors
27  use m_results_out
28  use m_abicore
29  use m_xmpi
30  use m_dtset
31 
32  use m_geometry,     only : mkrdim, xred2xcart
33  use m_parser,       only : prttagm, prttagm_images, ab_dimensions
34 
35  implicit none
36 
37  private

ABINIT/outvar_o_z [ Functions ]

[ Top ] [ Functions ]

NAME

 outvar_o_z

FUNCTION

 Echo variables between acell and gw_ ... (by alphabetic order) for the ABINIT code.

INPUTS

  choice= 1 if echo of preprocessed variables, 2 if echo after call driver
  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables
  iout=unit number for echoed output
  jdtset_(0:ndtset_alloc)=actual index of the dataset (equal to dtsets(:)%jdtset)
  marr=maximum number of numbers in an array (might need to be increased ... !)
  multivals= <type ab_dimensions>  either 0 or 1 , depending whether the
     dimension has different values for different datasets
  mxvals= <type ab_dimensions>
     maximum size of some arrays along all datasets, including
         lpawu      =maximal value of input lpawu for all the datasets
         gw_nqlwl   =maximal value of input gw_nqlwl for all the datasets
         mband      =maximum number of bands
         natom      =maximal value of input natom for all the datasets
         natpawu    =maximal value of number of atoms on which +U is applied for all the datasets
         natsph     =maximal value of input natsph for all the datasets
         natvshift  =maximal value of input natvshift for all the datasets
         nconeq     =maximal value of input nconeq for all the datasets
         nimage     =maximal value of input nimage for all the datasets
         nkptgw     =maximal value of input nkptgw for all the datasets
         nkpt       =maximal value of input nkpt for all the datasets
         nnos       =maximal value of input nnos for all the datasets
         nqptdm     =maximal value of input nqptdm for all the datasets
         nspinor    =maximal value of input nspinor for all the datasets
         nsppol     =maximal value of input nsppol for all the datasets
         nsym       =maximum number of symmetries
         ntypat     =maximum number of type of atoms
         nzchempot  =maximal value of input nzchempot for all the datasets
  ncid= NetCDF handler
  ndtset=number of datasets
  ndtset_alloc=number of datasets, corrected for allocation of at least
      one data set. Use for most dimensioned arrays.
  npsp=number of pseudopotentials
  prtvol_glob= if 0, minimal output volume, if 1, no restriction.
  results_out(0:ndtset_alloc)=<type results_out_type>contains the results
   needed for outvars, including evolving variables
  timopt=input variable to modulate the timing

OUTPUT

NOTES

 Note that this routine is called only by the processor me==0 .
 In consequence, no use of message and wrtout routine.
 The lines of code needed to output the defaults are preserved
 (see last section of the routine, but are presently disabled)

  Note that acell, occ, rprim, xred and vel might have been modified by the
  computation, so that their values if choice=1 or choice=2 will differ.

SOURCE

 104  subroutine outvar_o_z(choice,dtsets,iout,&
 105 & jdtset_,marr,multivals,mxvals,ncid,ndtset,ndtset_alloc,npsp,prtvol_glob,&
 106 & results_out,strimg,timopt)
 107 
 108 !Arguments ------------------------------------
 109 !scalars
 110  integer,intent(in) :: choice,iout,marr,ndtset
 111  integer,intent(in) :: ndtset_alloc,prtvol_glob,ncid,npsp,timopt
 112 !arrays
 113  integer,intent(in) :: jdtset_(0:ndtset_alloc)
 114  type(ab_dimensions),intent(in) :: multivals,mxvals
 115  type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc)
 116  type(results_out_type),intent(in) :: results_out(0:ndtset_alloc)
 117  character(len=8),intent(in) :: strimg(mxvals%nimage)
 118 
 119 !Local variables-------------------------------
 120 !scalars
 121  integer,parameter :: nkpt_max=50
 122  integer :: iat,icount,idtset,ii,iimage,ndtset_alloc_tmp
 123  integer :: narr
 124  integer :: multi_kptopt
 125  integer :: natom
 126  integer :: nimage,nnos,nsym
 127  integer :: ntypalch,ntypat,size1,size2,test_write,tnkpt,timopt_default,tmpimg0
 128  logical :: compute_static_images
 129  character(len=1) :: firstchar_gpu
 130 !arrays
 131  integer,allocatable :: narrm(:)
 132  integer,allocatable :: nimagem(:),prtimg(:,:)
 133  integer,allocatable :: intarr(:,:)
 134  real(dp) :: rprimd(3,3)
 135  real(dp),allocatable :: dprarr(:,:),dprarr_images(:,:,:)
 136  real(dp),allocatable :: xangst(:,:),xcart(:,:),xred(:,:)
 137  real(dp),allocatable :: xangst_(:,:,:,:),xcart_(:,:,:,:)
 138 
 139 ! *************************************************************************
 140 
 141 !###########################################################
 142 !### 01. Initial allocations and initialisations.
 143 
 144  ABI_MALLOC(dprarr,(marr,0:ndtset_alloc))
 145  ABI_MALLOC(dprarr_images,(marr,mxvals%nimage,0:ndtset_alloc))
 146  ABI_MALLOC(intarr,(marr,0:ndtset_alloc))
 147  ABI_MALLOC(narrm,(0:ndtset_alloc))
 148  ABI_MALLOC(nimagem,(0:ndtset_alloc))
 149  ABI_MALLOC(prtimg,(mxvals%nimage,0:ndtset_alloc))
 150 
 151  do idtset=0,ndtset_alloc
 152    nimagem(idtset)=dtsets(idtset)%nimage
 153  end do
 154 
 155  firstchar_gpu=' ';if (maxval(dtsets(1:ndtset_alloc)%use_gpu_cuda)>0) firstchar_gpu='-'
 156 
 157  natom=dtsets(1)%natom
 158  nimage=dtsets(1)%nimage
 159  nnos=dtsets(1)%nnos
 160  nsym=-1;nsym=dtsets(1)%nsym
 161  ntypalch=dtsets(1)%ntypalch
 162  ntypat=dtsets(1)%ntypat
 163 
 164 !###########################################################
 165 !### 02. Specific treatment for occopt, xangst, xcart, xred
 166 
 167 !Must compute xangst and xcart
 168  ABI_MALLOC(xangst_,(3,mxvals%natom,mxvals%nimage,0:ndtset_alloc))
 169  ABI_MALLOC(xcart_,(3,mxvals%natom,mxvals%nimage,0:ndtset_alloc))
 170  xangst_(:,:,:,:)=0.0_dp ; xcart_(:,:,:,:)=0.0_dp
 171 
 172  do idtset=1,ndtset_alloc
 173    natom=dtsets(idtset)%natom
 174    ABI_MALLOC(xred,(3,natom))
 175    ABI_MALLOC(xangst,(3,natom))
 176    ABI_MALLOC(xcart,(3,natom))
 177    do iimage=1,dtsets(idtset)%nimage
 178      xred(:,1:natom)=results_out(idtset)%xred(:,1:natom,iimage)
 179      call mkrdim(results_out(idtset)%acell(:,iimage),results_out(idtset)%rprim(:,:,iimage),rprimd)
 180 !    Compute xcart from xred and rprimd
 181      call xred2xcart(natom,rprimd,xcart,xred)
 182 !    Compute xangst from xcart
 183      xangst(:,:)=xcart(:,:)*Bohr_Ang
 184 !    Save the data
 185      xangst_(1:3,1:natom,iimage,idtset)=xangst(:,:)
 186      xcart_(1:3,1:natom,iimage,idtset)=xcart(:,:)
 187    end do
 188    if(dtsets(idtset)%nimage/=mxvals%nimage)then
 189      xangst_(1:3,1:natom,dtsets(idtset)%nimage+1:mxvals%nimage,idtset)=zero
 190      xcart_(1:3,1:natom,dtsets(idtset)%nimage+1:mxvals%nimage,idtset)=zero
 191    end if
 192    ABI_FREE(xred)
 193    ABI_FREE(xangst)
 194    ABI_FREE(xcart)
 195  end do
 196 
 197 !###########################################################
 198 !### 03. Print all the input variables (O)
 199 !##
 200 
 201 !occ
 202 !The use of prttagm for occ if occopt>=2 is not possible because
 203 !the different k-point and spins must be separated on different lines ...
 204 !Also, if prtvol=-1 and NC file has been created, only print one dataset
 205  ndtset_alloc_tmp=ndtset_alloc
 206  if(ncid<0)ndtset_alloc_tmp=1
 207  call prtocc(dtsets,iout,jdtset_,mxvals,ndtset_alloc_tmp,nimagem,prtvol_glob,results_out,strimg)
 208 
 209  intarr(1,:)=dtsets(:)%occopt
 210  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'occopt','INT',0)
 211 
 212  dprarr(1,:)=dtsets(:)%omegasimax
 213  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'omegasimax','ENE',0)
 214 
 215  dprarr(1,:)=dtsets(:)%omegasrdmax
 216  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'omegasrdmax','ENE',0)
 217 
 218  intarr(1,:)=dtsets(:)%optcell
 219  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'optcell','INT',0)
 220 
 221  intarr(1,:)=dtsets(:)%optdcmagpawu
 222  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'optdcmagpawu','INT',0)
 223 
 224  intarr(1,:)=dtsets(:)%optdriver
 225  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'optdriver','INT',0)
 226 
 227  intarr(1,:)=dtsets(:)%optforces
 228  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'optforces','INT',0)
 229 
 230  intarr(1,:)=dtsets(:)%optnlxccc
 231  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'optnlxccc','INT',0)
 232 
 233  intarr(1,:)=dtsets(:)%optstress
 234  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'optstress','INT',0)
 235 
 236  intarr(1,:)=dtsets(:)%orbmag
 237  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'orbmag','INT',0)
 238 
 239  intarr(1,:)=dtsets(:)%ortalg
 240  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ortalg','INT',0,firstchar=firstchar_gpu)
 241 
 242 !###########################################################
 243 !### 03. Print all the input variables (P)
 244 !##
 245 
 246  intarr(1,:)=dtsets(:)%paral_atom
 247  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'paral_atom','INT',0, firstchar="-")
 248  !call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'paral_atom','INT',0)
 249 
 250  intarr(1,:)=dtsets(:)%paral_kgb
 251  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'paral_kgb','INT',0)
 252 
 253  intarr(1,:)=dtsets(:)%paral_rf
 254  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'paral_rf','INT',0)
 255 
 256  intarr(1,:)=dtsets(:)%pawcpxocc
 257  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawcpxocc','INT',0)
 258 
 259  intarr(1,:)=dtsets(:)%pawcross
 260  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawcross','INT',0)
 261 
 262  dprarr(1,:)=dtsets(:)%pawecutdg
 263  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'pawecutdg','ENE',0)
 264 
 265  intarr(1,:)=dtsets(:)%pawfatbnd
 266  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawfatbnd','INT',0)
 267 
 268  intarr(1,:)=dtsets(:)%pawlcutd
 269  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawlcutd','INT',0)
 270 
 271  intarr(1,:)=dtsets(:)%pawlmix
 272  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawlmix','INT',0)
 273 
 274  intarr(1,:)=dtsets(:)%pawmixdg
 275  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawmixdg','INT',0)
 276 
 277  intarr(1,:)=dtsets(:)%pawnhatxc
 278  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawnhatxc','INT',0)
 279 
 280  intarr(1,:)=dtsets(:)%pawnphi
 281  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawnphi','INT',0)
 282 
 283  intarr(1,:)=dtsets(:)%pawntheta
 284  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawntheta','INT',0)
 285 
 286  intarr(1,:)=dtsets(:)%pawnzlm
 287  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawnzlm','INT',0)
 288 
 289  intarr(1,:)=dtsets(:)%pawoptmix
 290  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawoptmix','INT',0)
 291 
 292  intarr(1,:)=dtsets(:)%pawoptosc
 293  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawoptosc','INT',0)
 294 
 295  dprarr(1,:)=dtsets(:)%pawovlp
 296  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawovlp','DPR',0)
 297 
 298  intarr(1,:)=dtsets(:)%pawprtdos
 299  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawprtdos','INT',0)
 300 
 301  intarr(1,:)=dtsets(:)%pawprtvol
 302  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawprtvol','INT',0)
 303 
 304  intarr(1,:)=dtsets(:)%pawprtwf
 305  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawprtwf','INT',0)
 306 
 307  intarr(1,:)=dtsets(:)%pawprt_b
 308  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawprt_b','INT',0)
 309 
 310  intarr(1,:)=dtsets(:)%pawprt_k
 311  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawprt_k','INT',0)
 312 
 313  intarr(1,:)=dtsets(:)%pawspnorb
 314  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawspnorb','INT',0)
 315 
 316  intarr(1,:)=dtsets(:)%pawstgylm
 317  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawstgylm','INT',0)
 318 
 319  intarr(1,:)=dtsets(:)%pawsushat
 320  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawsushat','INT',0)
 321 
 322  intarr(1,:)=dtsets(:)%pawujat
 323  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawujat','INT',0)
 324 
 325  dprarr(1,:)=dtsets(:)%pawujrad
 326  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawujrad','LEN',0)
 327 
 328  dprarr(1,:)=dtsets(:)%pawujv
 329  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawujv','ENE',0)
 330 
 331  intarr(1,:)=dtsets(:)%pawusecp
 332  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawusecp','INT',0)
 333 
 334  intarr(1,:)=dtsets(:)%pawxcdev
 335  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pawxcdev','INT',0)
 336 
 337  intarr(1,:)=dtsets(:)%ph_intmeth
 338  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_intmeth','INT',0)
 339 
 340  intarr(1,:)=dtsets(:)%ph_ndivsm
 341  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ph_ndivsm','INT',0)
 342 
 343  do idtset=0,ndtset_alloc
 344    intarr(1:3,idtset)=dtsets(idtset)%ph_ngqpt
 345  end do
 346  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'ph_ngqpt','INT',0)
 347 
 348  intarr(1,:)=dtsets(:)%ph_nqpath
 349  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ph_nqpath','INT',0)
 350 
 351  intarr(1,:)=dtsets(:)%ph_nqshift
 352  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_nqshift','INT',0)
 353 
 354  dprarr(1,:)=dtsets(:)%ph_smear
 355  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_smear','ENE',0)
 356 
 357  dprarr(1,:)=dtsets(:)%ph_wstep
 358  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_wstep','ENE',0)
 359 
 360 !pimass
 361  icount=0
 362  do idtset=0, ndtset_alloc
 363    do ii = 1, ntypat
 364      dprarr(ii,idtset) = dtsets(idtset)%pimass(ii)
 365      if (dtsets(idtset)%pimass(ii)/=dtsets(idtset)%amu_orig(ii,1)) icount=1
 366    end do ! end loop over ntypat
 367  end do ! end loop over datasets
 368  if (icount/=0) then
 369    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,ntypat,narrm,ncid,ndtset_alloc,'pimass','DPR',0)
 370  end if
 371 
 372  intarr(1,:)=dtsets(:)%pimd_constraint
 373  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pimd_constraint','INT',0)
 374 
 375  intarr(1,:)=dtsets(:)%pitransform
 376  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'pitransform','INT',0)
 377 
 378  intarr(1,:)=dtsets(:)%plowan_bandi
 379  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'plowan_bandi','INT',0)
 380 
 381  intarr(1,:)=dtsets(:)%plowan_bandf
 382  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'plowan_bandf','INT',0)
 383 
 384  intarr(1,:)=dtsets(:)%plowan_compute
 385  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'plowan_compute','INT',0)
 386 
 387  intarr(1,:)=dtsets(:)%plowan_natom
 388  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'plowan_natom','INT',0)
 389 
 390  intarr(1,:)=dtsets(:)%plowan_nt
 391  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'plowan_nt','INT',0)
 392 
 393  intarr(1,:)=dtsets(:)%plowan_realspace
 394  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'plowan_realspace','INT',0)
 395 
 396 
 397 !plowan_it
 398  narr=100
 399  do idtset=0,ndtset_alloc       ! specific size for each dataset
 400    narrm(idtset)=3*dtsets(idtset)%plowan_nt
 401    !if(idtset==0)narrm(idtset)=100
 402    if (narrm(idtset)>0.and.dtsets(idtset)%plowan_compute>=0) then
 403      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%plowan_it(1:narrm(idtset))
 404    end if
 405  end do
 406  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'plowan_it','INT',1)
 407 
 408 
 409 !plowan_iatom
 410  narr=mxvals%natom
 411  do idtset=0,ndtset_alloc       ! specific size for each dataset
 412    narrm(idtset)=dtsets(idtset)%plowan_natom
 413    !if(idtset==0)narrm(idtset)=mxvals%natom
 414    if (narrm(idtset)>0.and.dtsets(idtset)%plowan_compute>=0) then
 415      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%plowan_iatom(1:narrm(idtset))
 416    end if
 417  end do
 418  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'plowan_iatom','INT',1)
 419 
 420 !plowan_nbl
 421  narr=mxvals%natom
 422  do idtset=0,ndtset_alloc       ! specific size for each dataset
 423    narrm(idtset)=dtsets(idtset)%plowan_natom
 424    !if(idtset==0)narrm(idtset)=mxvals%natom
 425    if (narrm(idtset)>0.and.dtsets(idtset)%plowan_compute>=0) then
 426      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%plowan_nbl(1:narrm(idtset))
 427    end if
 428  end do
 429  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'plowan_nbl','INT',1)
 430 
 431 !plowan_lcalc
 432  narr=12*mxvals%natom
 433  do idtset=0,ndtset_alloc       ! specific size for each dataset
 434    narrm(idtset)=sum(dtsets(idtset)%plowan_nbl(1:dtsets(idtset)%plowan_natom))
 435    !if(idtset==0)narrm(idtset)=12*mxvals%natom
 436    if (narrm(idtset)>0.and.dtsets(idtset)%plowan_compute>=0) then
 437      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%plowan_lcalc(1:narrm(idtset))
 438    end if
 439  end do
 440  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'plowan_lcalc','INT',1)
 441 
 442 !plowan_projcalc
 443  narr=12*mxvals%natom
 444  do idtset=0,ndtset_alloc       ! specific size for each dataset
 445    narrm(idtset)=sum(dtsets(idtset)%plowan_nbl(1:dtsets(idtset)%plowan_natom))
 446    !if(idtset==0)narrm(idtset)=12*mxvals%natom
 447    if (narrm(idtset)>0.and.dtsets(idtset)%plowan_compute>=0) then
 448      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%plowan_projcalc(1:narrm(idtset))
 449    end if
 450  end do
 451  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'plowan_projcalc','INT',1)
 452 
 453  dprarr(1,:)=dtsets(:)%polcen(1)
 454  dprarr(2,:)=dtsets(:)%polcen(2)
 455  dprarr(3,:)=dtsets(:)%polcen(3)
 456  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'polcen','DPR',0)
 457 
 458  intarr(1,:)=dtsets(:)%posdoppler
 459  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'posdoppler','INT',0)
 460 
 461  intarr(1,:)=dtsets(:)%positron
 462  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'positron','INT',0)
 463 
 464  intarr(1,:)=dtsets(:)%posnstep
 465  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'posnstep','INT',0)
 466 
 467  dprarr(1,:)=dtsets(:)%posocc
 468  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'posocc','DPR',0)
 469 
 470  dprarr(1,:)=dtsets(:)%postoldfe
 471  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'postoldfe','ENE',0)
 472 
 473  dprarr(1,:)=dtsets(:)%postoldff
 474  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'postoldff','DPR',0)
 475 
 476  dprarr(1,:)=dtsets(:)%ppmfrq
 477  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ppmfrq','ENE',0)
 478 
 479  intarr(1,:)=dtsets(:)%ppmodel
 480  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ppmodel','INT',0)
 481 
 482  intarr(1,:)=dtsets(:)%prepalw
 483  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prepalw','INT',0)
 484 
 485  intarr(1,:)=dtsets(:)%prepanl
 486  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prepanl','INT',0)
 487 
 488  intarr(1,:)=dtsets(:)%prepgkk
 489  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prepgkk','INT',0)
 490 
 491 !prtatlist
 492  if(multivals%natom==0)then
 493    do idtset=0,ndtset_alloc
 494      intarr(1:natom,idtset)=dtsets(idtset)%prtatlist(1:natom)
 495    end do
 496    intarr(1:mxvals%natom,0)=(/ (ii,ii=1,mxvals%natom) /)
 497    call prttagm(dprarr,intarr,iout,jdtset_,4,marr,natom,narrm,ncid,ndtset_alloc,'prtatlist','INT',0)
 498  else
 499 !  This thing will disapear with new generalized prttagm
 500  end if
 501 
 502  intarr(1,:)=dtsets(:)%prtbbb
 503  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtbbb','INT',0)
 504 
 505  intarr(1,:)=dtsets(:)%prtbltztrp
 506  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtbltztrp','INT',0)
 507 
 508  intarr(1,:)=dtsets(:)%prtchkprdm
 509  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtchkprdm','INT',0)
 510 
 511  intarr(1,:)=dtsets(:)%prtcif
 512  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtcif','INT',0)
 513 
 514  intarr(1,:)=dtsets(:)%prtden
 515  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtden','INT',0)
 516 
 517  intarr(1,:)=dtsets(:)%prtdensph
 518  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtdensph','INT',0)
 519 
 520  intarr(1,:)=dtsets(:)%prtdos
 521  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtdos','INT',0)
 522 
 523  intarr(1,:)=dtsets(:)%prtdosm
 524  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtdosm','INT',0)
 525 
 526  intarr(1,:)=dtsets(:)%prtebands
 527  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtebands','INT',0)
 528 
 529  intarr(1,:)=dtsets(:)%prtefmas
 530  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtefmas','INT',0)
 531 
 532  intarr(1,:)=dtsets(:)%prteig
 533  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prteig','INT',0)
 534 
 535  intarr(1,:)=dtsets(:)%prtelf
 536  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtelf','INT',0)
 537 
 538  intarr(1,:)=dtsets(:)%prteliash
 539  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prteliash','INT',0)
 540 
 541  intarr(1,:)=dtsets(:)%prtfull1wf
 542  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtfull1wf','INT',0)
 543 
 544  intarr(1,:)=dtsets(:)%prtfsurf
 545  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtfsurf','INT',0)
 546 
 547  intarr(1,:)=dtsets(:)%prtgden
 548  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtgden','INT',0)
 549 
 550  intarr(1,:)=dtsets(:)%prtgeo
 551  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtgeo','INT',0)
 552 
 553  intarr(1,:)=dtsets(:)%prtgkk
 554  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtgkk','INT',0)
 555 
 556  intarr(1,:)=dtsets(:)%prtgsr
 557  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtgsr','INT',0)
 558 
 559  intarr(1,:)=dtsets(:)%prtkden
 560  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtkden','INT',0)
 561 
 562  intarr(1,:)=dtsets(:)%prtlden
 563  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtlden','INT',0)
 564 
 565  intarr(1,:)=dtsets(:)%prtnabla
 566  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtnabla','INT',0)
 567 
 568  intarr(1,:)=dtsets(:)%prtnest
 569  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtnest','INT',0)
 570 
 571  intarr(1,:)=dtsets(:)%prtocc
 572  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtocc','INT',0)
 573 
 574  intarr(1,:)=dtsets(:)%prtphbands
 575  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtphbands','INT',0)
 576 
 577  intarr(1,:)=dtsets(:)%prtphdos
 578  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtphdos','INT',0)
 579 
 580  intarr(1,:)=dtsets(:)%prtphsurf
 581  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtphsurf','INT',0)
 582 
 583  intarr(1,:)=dtsets(:)%prtposcar
 584  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtposcar','INT',0)
 585 
 586  intarr(1,:)=dtsets(:)%prtprocar
 587  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtprocar','INT',0)
 588 
 589  intarr(1,:)=dtsets(:)%prtpot
 590  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtpot','INT',0)
 591 
 592  intarr(1,:)=dtsets(:)%prtpsps
 593  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtpsps','INT',0)
 594 
 595  intarr(1,:)=dtsets(:)%prtspcur
 596  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtspcur','INT',0)
 597 
 598  intarr(1,:)=dtsets(:)%prtstm
 599  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtstm','INT',0)
 600 
 601  intarr(1,:)=dtsets(:)%prtsuscep
 602  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtsuscep','INT',0)
 603 
 604  intarr(1,:)=dtsets(:)%prtvclmb
 605  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvclmb','INT',0)
 606 
 607  intarr(1,:)=dtsets(:)%prtvha
 608  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvha','INT',0)
 609 
 610  intarr(1,:)=dtsets(:)%prtvhxc
 611  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvhxc','INT',0)
 612 
 613  intarr(1,:)=dtsets(:)%prtkbff
 614  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtkbff','INT',0)
 615 
 616  intarr(1,:)=dtsets(:)%prtvol
 617  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvol','INT',0)
 618 
 619  intarr(1,:)=dtsets(:)%prtvolimg
 620  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvolimg','INT',0)
 621 
 622  intarr(1,:)=dtsets(:)%prtvpsp
 623  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvpsp','INT',0)
 624 
 625  intarr(1,:)=dtsets(:)%prtvxc
 626  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtvxc','INT',0)
 627 
 628  intarr(1,:)=dtsets(:)%prtwant
 629  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtwant','INT',0)
 630 
 631  intarr(1,:)=dtsets(:)%prtwf
 632  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtwf','INT',0)
 633 
 634  intarr(1,:)=dtsets(:)%prtwf_full
 635  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtwf_full','INT',0)
 636 
 637  intarr(1,:)=dtsets(:)%prtxml
 638  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prtxml','INT',0)
 639 
 640 !prt1dm
 641  intarr(1,:)=dtsets(:)%prt1dm
 642  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prt1dm','INT',0)
 643 
 644  !ptcharge
 645  do idtset=0, ndtset_alloc
 646    do ii = 1, ntypat
 647      dprarr(ii,idtset) = dtsets(idtset)%ptcharge(ii)
 648    end do ! end loop over ntypat
 649  end do ! end loop over datasets
 650  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,ntypat,narrm,ncid,ndtset_alloc,'ptcharge','DPR',0)
 651 
 652  intarr(1,:)=dtsets(:)%prt_lorbmag
 653  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'prt_lorbmag','INT',0)
 654 
 655  intarr(1,:)=dtsets(:)%ptgroupma
 656  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ptgroupma','INT',0)
 657 
 658  dprarr(1,:)=dtsets(:)%pvelmax(1)
 659  dprarr(2,:)=dtsets(:)%pvelmax(2)
 660  dprarr(3,:)=dtsets(:)%pvelmax(3)
 661  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'pvelmax','DPR',0)
 662 
 663  dprarr(1,:)=dtsets(:)%pw_unbal_thresh
 664  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'pw_unbal_thresh','DPR',0)
 665 
 666 !###########################################################
 667 !### 03. Print all the input variables (Q)
 668 !##
 669 
 670 !qmass
 671  narr=nnos ! default size for all datasets
 672  do idtset=0,ndtset_alloc       ! specific size for each dataset
 673    narrm(idtset)=dtsets(idtset)%nnos
 674    if(idtset==0)narrm(idtset)=mxvals%nnos
 675    if (narrm(idtset)>0) then
 676      dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%qmass(1:narrm(idtset))
 677    end if
 678  end do
 679  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'qmass','DPR',multivals%nnos)
 680 
 681  intarr(1,:)=dtsets(:)%qprtrb(1)
 682  intarr(2,:)=dtsets(:)%qprtrb(2)
 683  intarr(3,:)=dtsets(:)%qprtrb(3)
 684  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'qprtrb','INT',0)
 685 
 686  dprarr(1,:)=dtsets(:)%qptn(1)
 687  dprarr(2,:)=dtsets(:)%qptn(2)
 688  dprarr(3,:)=dtsets(:)%qptn(3)
 689  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'qpt','DPR',0)
 690 
 691 !qptdm
 692  narr=3*dtsets(1)%nqptdm ! default size for all datasets
 693  do idtset=0,ndtset_alloc       ! specific size for each dataset
 694    if(idtset/=0)then
 695      narrm(idtset)=3*dtsets(idtset)%nqptdm
 696      if (narrm(idtset)>0)&
 697 &     dprarr(1:narrm(idtset),idtset)=&
 698 &     reshape(dtsets(idtset)%qptdm(1:3,&
 699 &     1:dtsets(idtset)%nqptdm),&
 700 &     (/ narrm(idtset) /) )
 701    else
 702      narrm(idtset)=3*mxvals%nqptdm
 703      if (narrm(idtset)>0)&
 704 &     dprarr(1:narrm(idtset),idtset)=&
 705 &     reshape(dtsets(idtset)%qptdm(1:3,&
 706 &     1:mxvals%nqptdm),&
 707 &     (/ narrm(idtset) /) )
 708    end if
 709  end do
 710  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'qptdm','DPR',multivals%nqptdm)
 711 
 712  do idtset=0, ndtset_alloc
 713    do ii = 1, ntypat
 714      dprarr(ii,idtset) = dtsets(idtset)%quadmom(ii)
 715    end do ! end loop over ntypat
 716  end do ! end loop over datasets
 717  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,ntypat,narrm,ncid,ndtset_alloc,'quadmom','DPR',0)
 718 
 719  intarr(1,:)=dtsets(:)%quadquad
 720  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'quadquad','INT',0)
 721 
 722 
 723 !###########################################################
 724 !### 03. Print all the input variables (R)
 725 !##
 726 
 727 !variables used for the random positions in unit cell
 728  intarr(1,:)=dtsets(:)%random_atpos
 729  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'random_atpos','INT',0)
 730 
 731  dprarr(1,:)=dtsets(:)%ratsm
 732  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ratsm','LEN',0)
 733 
 734  do idtset=0, ndtset_alloc
 735    do ii = 1, ntypat
 736      dprarr(ii,idtset) = dtsets(idtset)%ratsph(ii)
 737    end do ! end loop over ntypat
 738  end do ! end loop over datasets
 739  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,ntypat,narrm,ncid,ndtset_alloc,'ratsph','LEN',0)
 740 
 741  dprarr = zero
 742  dprarr(1,:) = dtsets(:)%ratsph_extra
 743  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ratsph_extra','LEN',0)
 744 
 745  dprarr(1,:)=dtsets(:)%rcut
 746  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'rcut','LEN',0)
 747 
 748 !Variables used for recursion method
 749  dprarr(1,:)=dtsets(:)%recefermi
 750  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'recefermi','ENE',0)
 751 
 752  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%recgratio
 753  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'recgratio','INT',0)
 754 
 755  intarr(1,:)=dtsets(:)%recnpath
 756  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'recnpath','INT',0)
 757 
 758  intarr(1,:)=dtsets(:)%recnrec
 759  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'recnrec','INT',0)
 760 
 761  intarr(1,:)=dtsets(:)%recptrott
 762  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'recptrott','INT',0)
 763 
 764  dprarr(1,:)=dtsets(:)%recrcut
 765  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'recrcut','LEN',0)
 766 
 767  intarr(1,:)=dtsets(:)%rectesteg
 768  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rectesteg','INT',0)
 769 
 770  dprarr(1,:)=dtsets(:)%rectolden
 771  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'rectolden','DPR',0)
 772 
 773  dprarr(1,:)=dtsets(:)%red_dfield(1)    !!HONG
 774  dprarr(2,:)=dtsets(:)%red_dfield(2)
 775  dprarr(3,:)=dtsets(:)%red_dfield(3)
 776  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'red_dfield','DPR',0)
 777 
 778  dprarr(1,:)=dtsets(:)%red_efield(1)    !!HONG
 779  dprarr(2,:)=dtsets(:)%red_efield(2)
 780  dprarr(3,:)=dtsets(:)%red_efield(3)
 781  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'red_efield','DPR',0)
 782 
 783  dprarr(1,:)=dtsets(:)%red_efieldbar(1)   !!HONG
 784  dprarr(2,:)=dtsets(:)%red_efieldbar(2)
 785  dprarr(3,:)=dtsets(:)%red_efieldbar(3)
 786  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'red_efieldbar','DPR',0)
 787 
 788  intarr(1,:)=dtsets(:)%restartxf
 789  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'restartxf','INT',0)
 790 
 791  intarr(1,:)=dtsets(:)%rfasr
 792  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfasr','INT',0)
 793 
 794  test_write=0
 795  do idtset=1,ndtset_alloc
 796    if(dtsets(idtset)%rfatpol(1)/=1 .or. dtsets(idtset)%rfatpol(2)/=dtsets(idtset)%natom)test_write=1
 797  enddo
 798  if(test_write==1)then
 799    intarr(1,:)=dtsets(:)%rfatpol(1)
 800    intarr(2,:)=dtsets(:)%rfatpol(2)
 801    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'rfatpol','INT',0)
 802  endif
 803 
 804  intarr(1,:)=dtsets(:)%rfddk
 805  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfddk','INT',0)
 806 
 807  intarr(1,:)=dtsets(:)%rfdir(1)
 808  intarr(2,:)=dtsets(:)%rfdir(2)
 809  intarr(3,:)=dtsets(:)%rfdir(3)
 810  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'rfdir','INT',0)
 811 
 812  intarr(1,:)=dtsets(:)%rfelfd
 813  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfelfd','INT',0)
 814 
 815  intarr(1,:)=dtsets(:)%rfmagn
 816  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfmagn','INT',0)
 817 
 818  intarr(1,:)=dtsets(:)%rfmeth
 819  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfmeth','INT',0)
 820 
 821  intarr(1,:)=dtsets(:)%rfphon
 822  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfphon','INT',0)
 823 
 824  intarr(1,:)=dtsets(:)%rfstrs
 825  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfstrs','INT',0)
 826 
 827  intarr(1,:)=dtsets(:)%rfstrs_ref
 828  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfstrs_ref','INT',0)
 829 
 830  intarr(1,:)=dtsets(:)%rfuser
 831  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rfuser','INT',0)
 832 
 833  intarr(1,:)=dtsets(:)%rf2_dkdk
 834  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rf2_dkdk','INT',0)
 835 
 836  intarr(1,:)=dtsets(:)%rf2_dkde
 837  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'rf2_dkde','INT',0)
 838 
 839  intarr(1,:)=dtsets(:)%rf2_pert1_dir(1)
 840  intarr(2,:)=dtsets(:)%rf2_pert1_dir(2)
 841  intarr(3,:)=dtsets(:)%rf2_pert1_dir(3)
 842  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'rf2_pert1_dir','INT',0)
 843 
 844  intarr(1,:)=dtsets(:)%rf2_pert2_dir(1)
 845  intarr(2,:)=dtsets(:)%rf2_pert2_dir(2)
 846  intarr(3,:)=dtsets(:)%rf2_pert2_dir(3)
 847  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'rf2_pert2_dir','INT',0)
 848 
 849  dprarr(1,:)=dtsets(:)%rhoqpmix
 850  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'rhoqpmix','DPR',0)
 851 
 852  dprarr(1,:)=dtsets(:)%rifcsph
 853  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'rifcsph','DPR',0)
 854 
 855  intarr(1,:)=dtsets(:)%rmm_diis
 856  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'rmm_diis','INT',0)
 857  intarr(1,:)=dtsets(:)%rmm_diis_savemem
 858  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'rmm_diis_savemem','INT',0)
 859 
 860 !rprim
 861  prtimg(:,:)=1
 862  do idtset=0,ndtset_alloc
 863    narrm(idtset)=9
 864    do iimage=1,nimagem(idtset)
 865      if (narrm(idtset)>0) then
 866        dprarr_images(1:narrm(idtset),iimage,idtset)=&
 867 &       reshape(results_out(idtset)%rprim(1:3,1:3,iimage), (/ narrm(idtset) /) )
 868      end if
 869    end do
 870  end do
 871  call prttagm_images(dprarr_images,iout,jdtset_,-2,marr,narrm,ncid,ndtset_alloc,'rprim','DPR',&
 872 & mxvals%nimage,nimagem,ndtset,prtimg,strimg,forceprint=2)
 873 
 874 
 875 !###########################################################
 876 !### 03. Print all the input variables (S)
 877 !##
 878 
 879 !shiftk (printed only when kptopt>0)
 880  if(sum((dtsets(1:ndtset_alloc)%kptopt)**2)/=0)then
 881    multi_kptopt=0
 882    dprarr(:,0)=0.0_dp
 883    narr=3*dtsets(1)%nshiftk ! default size for all datasets
 884    do idtset=1,ndtset_alloc       ! specific size for each dataset
 885      narrm(idtset)=3*dtsets(idtset)%nshiftk
 886      if (narrm(idtset)>0) then
 887        dprarr(1:narrm(idtset),idtset)=&
 888 &       reshape(dtsets(idtset)%shiftk(1:3,1:dtsets(idtset)%nshiftk),(/ narrm(idtset) /) )
 889      end if
 890      if(dtsets(idtset)%kptopt<=0)then
 891        narrm(idtset)=0
 892        multi_kptopt=1
 893      end if
 894    end do
 895    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'shiftk','DPR',multivals%nshiftk)
 896 !  End of test to see whether kptopt/=0 for some dataset
 897  end if
 898 
 899  intarr(1,:)=dtsets(:)%sigma_bsum_range(1)
 900  intarr(2,:)=dtsets(:)%sigma_bsum_range(2)
 901  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'sigma_bsum_range','INT',0)
 902 
 903  dprarr(1,:)=dtsets(:)%sigma_erange(1)
 904  dprarr(2,:)=dtsets(:)%sigma_erange(2)
 905  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'sigma_erange','ENE',0)
 906 
 907  intarr(1,:)=dtsets(:)%transport_ngkpt(1)
 908  intarr(2,:)=dtsets(:)%transport_ngkpt(2)
 909  intarr(3,:)=dtsets(:)%transport_ngkpt(3)
 910  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'transport_ngkpt','INT',0)
 911 
 912  intarr(1,:)=dtsets(:)%sigma_ngkpt(1)
 913  intarr(2,:)=dtsets(:)%sigma_ngkpt(2)
 914  intarr(3,:)=dtsets(:)%sigma_ngkpt(3)
 915  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'sigma_ngkpt','INT',0)
 916 
 917  intarr(1,:)=dtsets(:)%signperm
 918  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'signperm','INT',0)
 919 
 920  dprarr(1,:)=dtsets(:)%slabwsrad
 921  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'slabwsrad','DPR',0)
 922 
 923  dprarr(1,:)=dtsets(:)%slabzbeg
 924  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'slabzbeg','DPR',0)
 925 
 926  dprarr(1,:)=dtsets(:)%slabzend
 927  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'slabzend','DPR',0)
 928 
 929  intarr(1,:)=dtsets(:)%slk_rankpp
 930  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'slk_rankpp','INT',0)
 931 
 932  intarr(1,:)=dtsets(:)%smdelta
 933  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'smdelta','INT',0)
 934 
 935  do idtset=0,ndtset_alloc
 936    intarr(1:npsp,idtset)=dtsets(idtset)%so_psp(1:npsp)
 937  end do
 938  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,npsp,narrm,ncid,ndtset_alloc,'so_psp','INT',0)
 939 
 940  dprarr(1,:)=dtsets(:)%spbroad
 941  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'spbroad','ENE',0)
 942 
 943  intarr(1,:)=dtsets(:)%spgroup
 944  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'spgroup','INT',0)
 945 
 946 !spinat
 947  dprarr(:,0)=0.0_dp
 948  narr=3*natom ! default size for all datasets
 949  do idtset=1,ndtset_alloc       ! specific size for each dataset
 950    narrm(idtset)=3*dtsets(idtset)%natom
 951    if (narrm(idtset)>0) then
 952      dprarr(1:narrm(idtset),idtset)=reshape(dtsets(idtset)%spinat(1:3,1:dtsets(idtset)%natom), (/narrm(idtset)/))
 953    end if
 954    if(sum(abs( dtsets(idtset)%spinat(1:3,1:dtsets(idtset)%natom))) < tol12 ) narrm(idtset)=0
 955  end do
 956  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,narrm,ncid,ndtset_alloc,'spinat','DPR',multivals%natom)
 957 
 958  dprarr(1,:)=dtsets(:)%spinmagntarget
 959  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'spinmagntarget','DPR',0)
 960 
 961  intarr(1,:)=dtsets(:)%spmeth
 962  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'spmeth','INT',0)
 963 
 964  dprarr(1,:)=dtsets(:)%spnorbscl
 965  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'spnorbscl','DPR',0)
 966 
 967  dprarr(1,:)=dtsets(:)%stmbias
 968  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'stmbias','DPR',0)
 969 
 970  dprarr(1,:)=dtsets(:)%strfact
 971  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'strfact','DPR',0)
 972 
 973  intarr(1,:)=dtsets(:)%string_algo
 974  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'string_algo','INT',0)
 975 
 976  do ii=1,6
 977    dprarr(ii,:)=dtsets(:)%strtarget(ii)
 978  end do
 979  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,6,narrm,ncid,ndtset_alloc,'strtarget','DPR',0)
 980 
 981 !strten
 982  if(choice==2)then
 983    prtimg(:,:)=1
 984    do idtset=0,ndtset_alloc       ! specific size for each dataset
 985      compute_static_images=(dtsets(idtset)%istatimg>0)
 986      narrm(idtset)=6
 987      if(dtsets(idtset)%iscf>=0)then
 988        do iimage=1,dtsets(idtset)%nimage
 989          if (narrm(idtset)>0) then
 990            dprarr_images(1:narrm(idtset),iimage,idtset)=results_out(idtset)%strten(:,iimage)
 991          end if
 992          if(.not.(dtsets(idtset)%dynimage(iimage)==1.or.compute_static_images))then
 993            prtimg(iimage,idtset)=0
 994          end if
 995        end do
 996      else
 997        narrm(idtset)=0
 998      end if
 999    end do
1000 !  This is a trick to force printing of strten even if zero, still not destroying the value of nimagem(0).
1001    tmpimg0=nimagem(0)
1002    nimagem(0)=0
1003    call prttagm_images(dprarr_images,iout,jdtset_,2,marr,narrm,ncid,ndtset_alloc,'strten','DPR',&
1004 &   mxvals%nimage,nimagem,ndtset,prtimg,strimg)
1005    nimagem(0)=tmpimg0
1006  end if
1007 
1008 !symafm
1009  intarr(:,0)=1
1010  narr=nsym ! default size for all datasets
1011  do idtset=1,ndtset_alloc       ! specific size for each dataset
1012    narrm(idtset)=dtsets(idtset)%nsym
1013    if (narrm(idtset)>0) then
1014      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%symafm(1:narrm(idtset))
1015    end if
1016  end do
1017  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'symafm','INT', multivals%nsym)
1018 
1019  intarr(1,:)=dtsets(:)%symchi
1020  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'symchi','INT',0)
1021 
1022  intarr(1,:)=dtsets(:)%symdynmat
1023  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'symdynmat','INT',0)
1024 
1025  intarr(1,:)=dtsets(:)%symmorphi
1026  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'symmorphi','INT',0)
1027 
1028 !symrel
1029  intarr(1:9,0)=(/ 1,0,0, 0,1,0, 0,0,1 /)
1030  narr=9*nsym ! default size for all datasets
1031  do idtset=1,ndtset_alloc       ! specific size for each dataset
1032    narrm(idtset)=9*dtsets(idtset)%nsym
1033    if (narrm(idtset)>0) then
1034      intarr(1:narrm(idtset),idtset)=&
1035 &     reshape(dtsets(idtset)%symrel(1:3,1:3,1:dtsets(idtset)%nsym), [narrm(idtset)] )
1036    end if
1037  end do
1038  call prttagm(dprarr,intarr,iout,jdtset_,3,marr,narr,narrm,ncid,ndtset_alloc,'symrel','INT', multivals%nsym)
1039 
1040  intarr(1,:)=dtsets(:)%symsigma
1041  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'symsigma','INT',0)
1042 
1043  intarr(1,:)=dtsets(:)%symv1scf
1044  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'symv1scf','INT',0)
1045 
1046 !###########################################################
1047 !### 03. Print all the input variables (T)
1048 !##
1049 
1050  dprarr(1,:)=dtsets(:)%td_maxene
1051  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'td_maxene','DPR',0)
1052 
1053  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%td_mexcit
1054  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'td_mexcit','INT',0)
1055 
1056  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%tfkinfunc
1057  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'tfkinfunc','INT',0)
1058 
1059  dprarr(1,:)=dtsets(:)%tfw_toldfe
1060  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tfw_toldfe','ENE',0)
1061 
1062  intarr(1,:)=dtsets(:)%tim1rev
1063  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'tim1rev','INT',0)
1064 
1065 
1066 !timopt
1067  timopt_default=1; if(xmpi_paral==1) timopt_default=0
1068 
1069  if(timopt/=timopt_default)then
1070    intarr(1,:)=timopt
1071    intarr(1,0)=timopt_default
1072    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'timopt','INT',0)
1073  end if
1074 
1075 !WVL - tails related variables
1076  intarr(1,:)=dtsets(:)%tl_nprccg
1077  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'tl_nprccg','INT',0)
1078  dprarr(1,:)=dtsets(:)%tl_radius
1079  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tl_radius','DPR',0)
1080 
1081 !tnons
1082  dprarr(:,0)=0.0_dp
1083  narr=3*nsym ! default size for all datasets
1084  do idtset=1,ndtset_alloc       ! specific size for each dataset
1085    narrm(idtset)=3*dtsets(idtset)%nsym
1086    if (narrm(idtset)>0) then
1087      dprarr(1:narrm(idtset),idtset)=reshape(dtsets(idtset)%tnons(1:3,1:dtsets(idtset)%nsym), [narrm(idtset)])
1088    end if
1089  end do
1090  call prttagm(dprarr,intarr,iout,jdtset_,-3,marr,narr,narrm,ncid,ndtset_alloc,'tnons','DPR',multivals%nsym)
1091 
1092  dprarr(1,:)=dtsets(:)%tolcum
1093  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolcum','DPR',0)
1094 
1095  dprarr(1,:)=dtsets(:)%toldfe
1096  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'toldfe','ENE',0)
1097 
1098  dprarr(1,:)=dtsets(:)%tolmxde
1099  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolmxde','ENE',0)
1100 
1101  dprarr(1,:)=dtsets(:)%toldff
1102  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'toldff','DPR',0)
1103 
1104  dprarr(1,:)=dtsets(:)%tolimg
1105  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolimg','ENE',0)
1106 
1107  dprarr(1,:)=dtsets(:)%tolmxf
1108  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolmxf','DPR',0)
1109 
1110  dprarr(1,:)=dtsets(:)%tolrde
1111  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolrde','DPR',0)
1112 
1113  dprarr(1,:)=dtsets(:)%tolrff
1114  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolrff','DPR',0)
1115 
1116  dprarr(1,:)=dtsets(:)%tolsym
1117  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolsym','DPR',0)
1118 
1119  dprarr(1,:)=dtsets(:)%tolvrs
1120  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolvrs','DPR',0)
1121 
1122  dprarr(1,:)=dtsets(:)%tolwfr
1123  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tolwfr','DPR',0)
1124 
1125  dprarr(1,:)=dtsets(:)%tphysel
1126  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tphysel','ENE',0)
1127 
1128  dprarr(1,:) = dtsets(:)%tmesh(1); dprarr(2,:) = dtsets(:)%tmesh(2); dprarr(3,:) = dtsets(:)%tmesh(3)
1129  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'tmesh','DPR',0)
1130 
1131  dprarr(1,:)=dtsets(:)%tsmear
1132  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'tsmear','ENE',0)
1133 
1134 !typat
1135  narr=natom                      ! default size for all datasets
1136  do idtset=0,ndtset_alloc       ! specific size for each dataset
1137    narrm(idtset)=dtsets(idtset)%natom
1138    if(idtset==0)narrm(idtset)=mxvals%natom
1139    if (narrm(idtset)>0) then
1140      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%typat(1:narrm(idtset))
1141    end if
1142  end do
1143  call prttagm(dprarr,intarr,iout,jdtset_,4,marr,narr,narrm,ncid,ndtset_alloc,'typat','INT',multivals%natom,forceprint=2)
1144 
1145 !###########################################################
1146 !### 03. Print all the input variables (U)
1147 !##
1148  intarr(1,:)=dtsets(:)%ucrpa
1149  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ucrpa','INT',0)
1150 
1151  intarr(1,:)=dtsets(:)%ucrpa_bands(1)
1152  intarr(2,:)=dtsets(:)%ucrpa_bands(2)
1153  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'ucrpa_bands','INT',0)
1154 
1155  dprarr(1,:)=dtsets(:)%ucrpa_window(1)
1156  dprarr(2,:)=dtsets(:)%ucrpa_window(2)
1157  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'ucrpa_window','ENE',0)
1158 
1159 !upawu
1160  prtimg(:,:)=1
1161  do idtset=0,ndtset_alloc
1162    narrm(idtset)=dtsets(idtset)%ntypat
1163    if (idtset==0) narrm(idtset)=mxvals%ntypat
1164    do iimage=1,nimagem(idtset)
1165      if (narrm(idtset)>0) then
1166        dprarr_images(1:narrm(idtset),iimage,idtset)=dtsets(idtset)%upawu(1:narrm(idtset),iimage)
1167      end if
1168    end do
1169  end do
1170  call prttagm_images(dprarr_images,iout,jdtset_,1,marr,narrm,&
1171 & ncid,ndtset_alloc,'upawu','ENE',mxvals%nimage,nimagem,ndtset,prtimg,strimg)
1172 
1173  intarr(1,:)=dtsets(:)%usedmatpu
1174  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usedmatpu','INT',0)
1175 
1176  intarr(1,:)=dtsets(:)%usedmft
1177  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usedmft','INT',0)
1178 
1179  intarr(1,:)=dtsets(:)%useexexch
1180  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'useexexch','INT',0)
1181 
1182  intarr(1,:)=dtsets(:)%usefock
1183  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usefock','INT',0)
1184 
1185  intarr(1,:)=dtsets(:)%usepotzero
1186  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usepotzero','INT',0)
1187 
1188  intarr(1,:)=dtsets(:)%usekden
1189  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usekden','INT',0)
1190 
1191  intarr(1,:)=dtsets(:)%use_gemm_nonlop
1192  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'use_gemm_nonlop','INT',0)
1193 
1194  intarr(1,:)=dtsets(:)%useextfpmd
1195  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'useextfpmd','INT',0)
1196 
1197  intarr(1,:)=dtsets(:)%use_yaml
1198  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'use_yaml','INT',0)
1199 
1200  intarr(1,:)=dtsets(:)%use_nonscf_gkk
1201  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'use_nonscf_gkk','INT',0)
1202 
1203  intarr(1,:)=dtsets(:)%usepawu
1204  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usepawu','INT',0)
1205 
1206  intarr(1,:)=dtsets(:)%usepead
1207  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usepead','INT',0)
1208 
1209  intarr(1,:)=dtsets(:)%useria
1210  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'useria','INT',0)
1211 
1212  intarr(1,:)=dtsets(:)%userib
1213  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'userib','INT',0)
1214 
1215  intarr(1,:)=dtsets(:)%useric
1216  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'useric','INT',0)
1217 
1218  intarr(1,:)=dtsets(:)%userid
1219  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'userid','INT',0)
1220 
1221  intarr(1,:)=dtsets(:)%userie
1222  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'userie','INT',0)
1223 
1224  dprarr(1,:)=dtsets(:)%userra
1225  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'userra','DPR',0)
1226 
1227  dprarr(1,:)=dtsets(:)%userrb
1228  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'userrb','DPR',0)
1229 
1230  dprarr(1,:)=dtsets(:)%userrc
1231  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'userrc','DPR',0)
1232 
1233  dprarr(1,:)=dtsets(:)%userrd
1234  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'userrd','DPR',0)
1235 
1236  dprarr(1,:)=dtsets(:)%userre
1237  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'userre','DPR',0)
1238 
1239  intarr(1,:)=dtsets(:)%usewvl
1240  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usewvl','INT',0)
1241 
1242  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%usexcnhat_orig
1243  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'usexcnhat','INT',0)
1244 
1245  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%useylm
1246  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'useylm','INT',0,firstchar=firstchar_gpu)
1247 
1248  intarr(1,:)=dtsets(:)%use_gpu_cuda
1249  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'use_gpu_cuda','INT',0,firstchar=firstchar_gpu)
1250 
1251  intarr(1,:)=dtsets(:)%use_slk
1252  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'use_slk','INT',0, firstchar="-")
1253 
1254  intarr(1,:)=dtsets(:)%use_oldchi
1255  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'use_oldchi','INT',0)
1256 
1257 
1258 !###########################################################
1259 !### 03. Print all the input variables (V)
1260 !##
1261 
1262  dprarr(1,:)=dtsets(:)%vcutgeo(1)
1263  dprarr(2,:)=dtsets(:)%vcutgeo(2)
1264  dprarr(3,:)=dtsets(:)%vcutgeo(3)
1265  call prttagm(dprarr,intarr,iout,jdtset_,3,marr,3,narrm,ncid,ndtset_alloc,'vcutgeo','DPR',0)
1266 
1267  if(sum(dtsets(1:ndtset_alloc)%prtwant) >1)then
1268 !  van der Waals correction with MLWFs related variables
1269    if(any(dtsets(1:ndtset_alloc)%vdw_xc==10).or.any(dtsets(1:ndtset_alloc)%vdw_xc==11).or.&
1270 &   any(dtsets(1:ndtset_alloc)%vdw_xc==14))then
1271      intarr(1,:)=dtsets(:)%vdw_nfrag
1272      call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'vdw_nfrag','INT',0)
1273    end if !vdw_xc==10,11,14
1274    if(any(dtsets(1:ndtset_alloc)%vdw_xc==10).or.any(dtsets(1:ndtset_alloc)%vdw_xc==11).or.&
1275 &   any(dtsets(1:ndtset_alloc)%vdw_xc==14))then
1276      intarr(1,:)=dtsets(:)%vdw_supercell(1)
1277      intarr(2,:)=dtsets(:)%vdw_supercell(2)
1278      intarr(3,:)=dtsets(:)%vdw_supercell(3)
1279      call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'vdw_supercell','INT',0)
1280    end if !vdw_xc==10,11,14
1281  end if !prtwant>1
1282 
1283  dprarr(1,:)=dtsets(:)%vdw_tol
1284  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'vdw_tol','DPR',0)
1285  dprarr(1,:)=dtsets(:)%vdw_tol_3bt
1286  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'vdw_tol_3bt','DPR',0)
1287 
1288  if(sum(dtsets(1:ndtset_alloc)%prtwant) >1)then
1289 !  van der Waals correction with MLWFs related variables
1290    if(any(dtsets(1:ndtset_alloc)%vdw_xc==10).or.any(dtsets(1:ndtset_alloc)%vdw_xc==11))then
1291      do iat=1,mxvals%natom
1292        intarr(iat,:)=dtsets(:)%vdw_typfrag(iat)
1293      end do
1294      call prttagm(dprarr,intarr,iout,jdtset_,2,marr,mxvals%natom,narrm,ncid,ndtset_alloc,'vdw_typfrag','INT',0)
1295    end if !vdw_xc==10 or xc==11
1296  end if !prtwant>1
1297 
1298  intarr(1,:)=dtsets(:)%vdw_xc
1299  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'vdw_xc','INT',0)
1300 
1301  if(sum(dtsets(1:ndtset_alloc)%prtvdw) >1)then
1302    if(any(dtsets(1:ndtset_alloc)%vdw_xc<10))then
1303      dprarr(1,:)=dtsets(:)%vdw_df_threshold
1304      call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'vdw_df_threshold','ENE',0)
1305    end if
1306  end if
1307 
1308 !vel
1309  prtimg(:,:)=1
1310  do idtset=0,ndtset_alloc
1311    if(idtset/=0)then
1312      size1=dtsets(idtset)%natom
1313    else
1314      size1=mxvals%natom
1315    end if
1316    narrm(idtset)=3*size1
1317    do iimage=1,nimagem(idtset)
1318      if (narrm(idtset)>0) then
1319        dprarr_images(1:narrm(idtset),iimage,idtset)=&
1320 &       reshape(results_out(idtset)%vel(1:3,1:size1,iimage), (/ narrm(idtset) /) )
1321      end if
1322    end do
1323  end do
1324  call prttagm_images(dprarr_images,iout,jdtset_,2,marr,narrm,ncid,ndtset_alloc,'vel','DPR',&
1325 & mxvals%nimage,nimagem,ndtset,prtimg,strimg)
1326 
1327 !vel_cell
1328 !At present, vel_cell does not depend on image... but this might change in the future.
1329  prtimg(:,:)=1
1330  if (.true.) then
1331 !  if(mxvals%nimage==1)then
1332    do idtset=0,ndtset_alloc
1333      dprarr(1:9,idtset)= reshape(results_out(idtset)%vel_cell(:,:,1),(/9/))
1334    end do
1335    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,9,narrm,ncid,ndtset_alloc,'vel_cell','DPR',0)
1336 !  else
1337 !  do idtset=1,ndtset_alloc       ! specific size for each dataset
1338 !  nimagem(idtset)=dtsets(idtset)%nimage
1339 !  narrm(idtset)=9
1340 !  do iimage=1,dtsets(idtset)%nimage
1341 !  if (narrm(idtset)>0) then
1342 !  dprarr_images(1:narrm(idtset),iimage,idtset)=&
1343 !  &         reshape(results_out(idtset)%vel_cell(1:3,1:3,iimage),&
1344 !  &         (/ narrm(idtset) /) )
1345 !  end if
1346 !  end do
1347 !  end do
1348 !  call prttagm_images(dprarr_images,iout,jdtset_,&
1349 !  &   marr,narrm,ncid,ndtset_alloc,'vel_cell',&
1350 !  &   mxvals%nimage,nimagem,ndtset,prtimg,strimg)
1351  end if
1352 
1353  dprarr(1,:)=dtsets(:)%vis
1354  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'vis','DPR',0)
1355 
1356  dprarr(1,:)=dtsets(:)%vprtrb(1)
1357  dprarr(2,:)=dtsets(:)%vprtrb(2)
1358  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'vprtrb','ENE',0)
1359 
1360 
1361 !###########################################################
1362 !### 03. Print all the input variables (W)
1363 !##
1364  intarr(1,:)=dtsets(:)%wfinit
1365  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'wfinit','INT',0)
1366 
1367  dprarr(1,:)=dtsets(:)%wfmix
1368  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'wfmix','DPR',0)
1369 
1370  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%wfk_task
1371  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'wfk_task','INT',0)
1372 
1373  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%wfoptalg
1374  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'wfoptalg','INT',0,firstchar=firstchar_gpu)
1375 
1376 !wtatcon
1377  narr=3*natom*dtsets(1)%nconeq ! default size for all datasets
1378  do idtset=0,ndtset_alloc       ! specific size for each dataset
1379    if(idtset/=0)then
1380      narrm(idtset)=3*dtsets(idtset)%natom*dtsets(idtset)%nconeq
1381      if (narrm(idtset)>0)&
1382 &     dprarr(1:narrm(idtset),idtset)=&
1383 &     reshape(dtsets(idtset)%wtatcon(1:3,1:dtsets(idtset)%natom,1:dtsets(idtset)%nconeq),(/ narrm(idtset) /) )
1384    else
1385      narrm(idtset)=3*mxvals%natom*mxvals%nconeq
1386      if (narrm(idtset)>0)&
1387 &     dprarr(1:narrm(idtset),idtset)=&
1388 &     reshape(dtsets(idtset)%wtatcon(1:3,1:mxvals%natom,1:mxvals%nconeq),(/ narrm(idtset) /) )
1389    end if
1390  end do
1391 
1392  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'wtatcon','DPR',&
1393    multivals%natom+multivals%nconeq)
1394 
1395 !wtk
1396  if (allocated(dtsets(0)%wtk)) then
1397    tnkpt=0
1398    dprarr(:,0)=1
1399    narr=dtsets(1)%nkpt ! default size for all datasets
1400    if(prtvol_glob==0 .and. narr>nkpt_max)then
1401      narr=nkpt_max
1402      tnkpt=1
1403    end if
1404    do idtset=1,ndtset_alloc       ! specific size for each dataset
1405      narrm(idtset)=dtsets(idtset)%nkpt
1406      if (narrm(idtset)>0) dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%wtk(1:narrm(idtset))+tol12
1407 
1408      if(prtvol_glob==0 .and. narrm(idtset)>nkpt_max)then
1409        narrm(idtset)=nkpt_max
1410        tnkpt=1
1411      end if
1412    end do
1413    call prttagm(dprarr,intarr,iout,jdtset_,4,marr,narr,narrm,ncid,ndtset_alloc,'wtk','DPR',multivals%nkpt)
1414    if(tnkpt==1) write(iout,'(23x,a,i3,a)' ) 'outvars : Printing only first ',nkpt_max,' k-points.'
1415  end if
1416 
1417 !WVL - wavelets variables
1418  if (any(dtsets(:)%usewvl==1)) then
1419    intarr(1,:)=dtsets(:)%wvl_bigdft_comp
1420    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'wvl_bigdft_comp','INT',0)
1421    dprarr(1,:)=dtsets(:)%wvl_crmult
1422    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'wvl_crmult','DPR',0)
1423    dprarr(1,:)=dtsets(:)%wvl_frmult
1424    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'wvl_frmult','DPR',0)
1425    dprarr(1,:)=dtsets(:)%wvl_hgrid
1426    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'wvl_hgrid','DPR',0)
1427    intarr(1,:)=dtsets(:)%wvl_nprccg
1428    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'wvl_nprccg','INT',0)
1429  end if
1430 
1431 !Wannier90 interface related variables
1432  if(sum(dtsets(1:ndtset_alloc)%prtwant) >1)then
1433    intarr(1,:)=dtsets(:)%w90iniprj
1434    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'w90iniprj','INT',0)
1435    intarr(1,:)=dtsets(:)%w90prtunk
1436    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'w90prtunk','INT',0)
1437  end if !prtwant>1
1438 
1439 !###########################################################
1440 !### 03. Print all the input variables (X)
1441 !##
1442 
1443 !xangst
1444  prtimg(:,:)=1
1445  do idtset=0,ndtset_alloc
1446    if(idtset/=0)then
1447      size1=dtsets(idtset)%natom
1448    else
1449      size1=mxvals%natom
1450    end if
1451    narrm(idtset)=3*size1
1452    do iimage=1,nimagem(idtset)
1453      if (narrm(idtset)>0) then
1454        dprarr_images(1:narrm(idtset),iimage,idtset)=reshape(xangst_(1:3,1:size1,iimage,idtset), (/narrm(idtset)/))
1455      end if
1456    end do
1457  end do
1458 
1459  call prttagm_images(dprarr_images,iout,jdtset_,-2,marr,narrm,ncid,ndtset_alloc,'xangst','DPR',&
1460 & mxvals%nimage,nimagem,ndtset,prtimg,strimg)
1461 
1462 !xcart
1463  prtimg(:,:)=1
1464  do idtset=0,ndtset_alloc
1465    if(idtset/=0)then
1466      size1=dtsets(idtset)%natom
1467    else
1468      size1=mxvals%natom
1469    end if
1470    narrm(idtset)=3*size1
1471    do iimage=1,nimagem(idtset)
1472      if (narrm(idtset)>0) then
1473        dprarr_images(1:narrm(idtset),iimage,idtset)=reshape(xcart_(1:3,1:size1,iimage,idtset), (/ narrm(idtset) /) )
1474      end if
1475    end do
1476  end do
1477 
1478  call prttagm_images(dprarr_images,iout,jdtset_,-2,marr,narrm,ncid,ndtset_alloc,'xcart','DPR',&
1479 & mxvals%nimage,nimagem,ndtset,prtimg,strimg)
1480 
1481  dprarr(1,:)=dtsets(:)%xc_denpos
1482  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'xc_denpos','DPR',0)
1483 
1484  dprarr(1,:)=dtsets(:)%xc_tb09_c
1485  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'xc_tb09_c','DPR',0)
1486 
1487  intarr(1,:)=dtsets(:)%x1rdm
1488  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'x1rdm','INT',0)
1489 
1490 !xred
1491  prtimg(:,:)=1
1492  do idtset=0,ndtset_alloc
1493    if(idtset/=0)then
1494      size2=dtsets(idtset)%natom
1495    else
1496      size2=mxvals%natom
1497    end if
1498    narrm(idtset)=3*size2
1499    do iimage=1,nimagem(idtset)
1500      if (narrm(idtset)>0) then
1501        dprarr_images(1:narrm(idtset),iimage,idtset)=&
1502 &       reshape(results_out(idtset)%xred(:,1:size2,iimage), (/ narrm(idtset) /) )
1503      end if
1504    end do
1505  end do
1506  call prttagm_images(dprarr_images,iout,jdtset_,-2,marr,narrm,ncid,ndtset_alloc,'xred','DPR',&
1507 & mxvals%nimage,nimagem,ndtset,prtimg,strimg,forceprint=2)
1508 
1509 !xredsph_extra
1510  do idtset=0,ndtset_alloc
1511    if(idtset/=0)then
1512      size2=dtsets(idtset)%natsph_extra
1513    else
1514      size2=0
1515    end if
1516    narrm(idtset)=3*size2
1517    if (narrm(idtset)>0) then
1518      dprarr(1:narrm(idtset),idtset)= reshape(dtsets(idtset)%xredsph_extra(:,1:size2), (/ narrm(idtset) /) )
1519    end if
1520  end do
1521  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'xredsph_extra','DPR',1)
1522 
1523 
1524 !###########################################################
1525 !### 03. Print all the input variables (Y)
1526 !##
1527 
1528 !###########################################################
1529 !### 03. Print all the input variables (Z)
1530 !##
1531 
1532  dprarr(1,:)=dtsets(:)%zcut
1533  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'zcut','ENE',0)
1534 
1535 !zeemanfield
1536  dprarr(1,:)=dtsets(:)%zeemanfield(1)
1537  dprarr(2,:)=dtsets(:)%zeemanfield(2)
1538  dprarr(3,:)=dtsets(:)%zeemanfield(3)
1539  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'zeemanfield','BFI',0)
1540 
1541 !ziontypat   ! After all, should always echo this value
1542  if(sum(dtsets(:)%ntypalch)>0)then
1543    narr=ntypat                    ! default size for all datasets
1544    do idtset=0,ndtset_alloc       ! specific size for each dataset
1545      narrm(idtset)=dtsets(idtset)%ntypat
1546      if(idtset==0)narrm(idtset)=mxvals%ntypat
1547      if (narrm(idtset)>0) then
1548        dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%ziontypat(1:narrm(idtset))
1549      end if
1550    end do
1551    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
1552 &   narrm,ncid,ndtset_alloc,'ziontypat','DPR',multivals%ntypat,forceprint=2)
1553  end if
1554 
1555  do idtset=0,ndtset_alloc
1556    dprarr(1:npsp,idtset)=dtsets(idtset)%znucl(1:npsp)
1557  end do
1558  call prttagm(dprarr,intarr,iout,jdtset_,4,marr,npsp,narrm,ncid,ndtset_alloc,'znucl','DPR',0,forceprint=2)
1559 
1560 !###########################################################
1561 !## Deallocation for generic arrays, and for n-z variables
1562 
1563  ABI_FREE(dprarr)
1564  ABI_FREE(intarr)
1565  ABI_FREE(narrm)
1566  ABI_FREE(nimagem)
1567  ABI_FREE(dprarr_images)
1568  ABI_FREE(prtimg)
1569  ABI_FREE(xangst_)
1570  ABI_FREE(xcart_)
1571 
1572 contains

ABINIT/prtocc [ Functions ]

[ Top ] [ Functions ]

NAME

 prtocc

FUNCTION

 Print the content of occ.
 Due to the need to distinguish between different k-points and
 different spin polarisations, prttagm.f cannot be used.
 So, need a dedicated routine.

INPUTS

  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables
  iout=unit number for echoed output
  jdtset_(0:ndtset_alloc)=list of dataset indices.
  ndtset_alloc=govern second dimension of intarr and dprarr
  prtvol_glob= if 0, minimal output volume, if 1, no restriction.
  results_out(0:ndtset_alloc)=<type results_out_type>contains the results
   needed for outvars, including occ, an evolving variable

OUTPUT

  (only writing)

SOURCE

1600 subroutine prtocc(dtsets,iout,jdtset_,mxvals,ndtset_alloc,nimagem,prtvol_glob,results_out,strimg)
1601 
1602 !Arguments ------------------------------------
1603 !scalars
1604  integer,intent(in) :: iout,ndtset_alloc,prtvol_glob
1605 !arrays
1606  integer,intent(in) :: jdtset_(0:ndtset_alloc)
1607  integer,intent(in) :: nimagem(0:ndtset_alloc)
1608  type(ab_dimensions),intent(in) :: mxvals
1609  type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc)
1610  type(results_out_type),intent(in) :: results_out(0:ndtset_alloc)
1611  character(len=8),intent(in) :: strimg(mxvals%nimage)
1612 
1613 !Local variables-------------------------------
1614  character(len=*), parameter :: f_occ    ="(1x,a16,1x,(t22,6f10.6))"
1615  character(len=*), parameter :: f_occa   ="(1x,a16,a,1x,(t22,6f10.6))"
1616  character(len=*), parameter :: token='occ'
1617 !scalars
1618  integer,parameter :: nkpt_max=50
1619  integer :: generic,iban,idtset,ikpsp,ikpt,isppol,jdtset,multi,multi_nband
1620  integer :: multi_nimage
1621  integer :: multi_nkpt,multi_nsppol,multi_occopt,nban,nkpt,nkpt_eff
1622  integer :: multi_tsmear
1623  integer :: print,tnkpt
1624  logical, allocatable :: test_multiimages(:)
1625  character(len=4) :: appen
1626  character(len=16) :: keywd
1627  character(len=500) :: message
1628 
1629 ! *************************************************************************
1630 
1631  if(ndtset_alloc<1)then
1632    write(message, '(a,i0,a)' )' ndtset_alloc=',ndtset_alloc,', while it should be >= 1.'
1633    ABI_BUG(message)
1634  end if
1635 
1636  if(ndtset_alloc>9999)then
1637    write(message, '(a,i0,a)' )' ndtset_alloc=',ndtset_alloc,', while it must be lower than 100.'
1638    ABI_BUG(message)
1639  end if
1640 
1641 !It is important to take iscf into account, since when it is -2, occupation numbers must be ignored
1642 
1643  multi_occopt=0
1644  if(ndtset_alloc>1)then
1645    do idtset=1,ndtset_alloc
1646      if(dtsets(1)%occopt/=dtsets(idtset)%occopt .and. dtsets(idtset)%iscf/=-2 )multi_occopt=1
1647    end do
1648  end if
1649 
1650  multi_tsmear=0
1651  if(ndtset_alloc>1)then
1652    do idtset=1,ndtset_alloc
1653      if(dtsets(1)%tsmear/=dtsets(idtset)%tsmear .and. dtsets(idtset)%iscf/=-2 )multi_tsmear=1
1654    end do
1655  end if
1656 
1657  multi_nkpt=0
1658  if(ndtset_alloc>1)then
1659    do idtset=1,ndtset_alloc
1660      if(dtsets(1)%nkpt/=dtsets(idtset)%nkpt .and. dtsets(idtset)%iscf/=-2 )multi_nkpt=1
1661    end do
1662  end if
1663  if(multi_nkpt==0)nkpt=dtsets(1)%nkpt
1664 
1665  multi_nsppol=0
1666  if(ndtset_alloc>1)then
1667    do idtset=1,ndtset_alloc
1668      if(dtsets(1)%nsppol/=dtsets(idtset)%nsppol .and. dtsets(idtset)%iscf/=-2 )multi_nsppol=1
1669    end do
1670  end if
1671 
1672  if(multi_nsppol==0 .and. multi_nkpt==0)then
1673    multi_nband=0
1674    if(ndtset_alloc>1)then
1675      do idtset=1,ndtset_alloc
1676        if(dtsets(idtset)%iscf/=-2)then
1677          do ikpsp=1,dtsets(1)%nkpt*dtsets(1)%nsppol
1678            if(dtsets(1)%nband(ikpsp)/=dtsets(idtset)%nband(ikpsp))multi_nband=1
1679          end do
1680        end if
1681      end do
1682    end if
1683  else
1684    multi_nband=1
1685  end if
1686 
1687  multi_nimage=0
1688  if(ndtset_alloc>1)then
1689    do idtset=1,ndtset_alloc
1690      if(dtsets(1)%nimage/=dtsets(idtset)%nimage .and. dtsets(idtset)%iscf/=-2 )multi_nimage=1
1691    end do
1692  end if
1693 
1694 !DEBUG
1695 ! write(std_out,*)' prtocc : 2, multi_nimage= ',multi_nimage
1696 !ENDDEBUG
1697 
1698 !Test whether for this variable, the content of different images differ.
1699 !test_multiimages(idtset)=.false. if, for that dataset, the content for different
1700 !images is identical.
1701  ABI_MALLOC(test_multiimages,(0:ndtset_alloc))
1702  test_multiimages=.false.
1703  do idtset=1,ndtset_alloc
1704    if(nimagem(idtset)>1)then
1705      nban=sum(dtsets(idtset)%nband(1:dtsets(idtset)%nsppol*dtsets(idtset)%nkpt))
1706      do iban=1,nban
1707        if(sum(abs( results_out(idtset)%occ(iban,2:nimagem(idtset))- results_out(idtset)%occ(iban,1)))>tol12)then
1708          test_multiimages(idtset)=.true.
1709        end if
1710      end do
1711    end if
1712  end do
1713  if(nimagem(0)==0)test_multiimages(0)=.true.
1714 
1715 !DEBUG
1716 ! write(std_out,*)' prtocc : 3, test_multiimages= ',test_multiimages
1717 ! write(std_out,*)' prtocc : multi_occopt, multi_nband, multi_nimage=',multi_occopt, multi_nband, multi_nimage
1718 ! write(std_out,*)' prtocc : test_multiimages(1:ndtset_alloc)=',test_multiimages(1:ndtset_alloc)
1719 ! write(std_out,*)' prtocc : any(test_multiimages(1:ndtset_alloc))=',any(test_multiimages(1:ndtset_alloc))
1720 !ENDDEBUG
1721 
1722 !There is a possibility of a single generic occupation-number set (common to all datasets) if
1723 !multi_occopt==0 and multi_nband==0  and (multi_nimage==0  or the content of the different images is always the same)
1724 !This might occur even if occupation numbers differ for different images.
1725  multi=1
1726  if(multi_occopt==0 .and. multi_nband==0 .and. (multi_nimage==0 .or. .not. any(test_multiimages(1:ndtset_alloc)))) then
1727    nban=sum(dtsets(1)%nband(1:dtsets(1)%nsppol*dtsets(1)%nkpt))
1728    multi=0
1729    if(ndtset_alloc>1)then
1730      do idtset=1,ndtset_alloc
1731        if(dtsets(idtset)%iscf/=-2)then
1732 !        nban counts all bands and kpoints and spins: see above
1733          do iimage=1,nimagem(idtset)
1734            if(iimage==1 .or. test_multiimages(idtset))then
1735              do iban=1,nban
1736 !              Use of tol8, because the format for multi=1 is f16.6, so will not
1737 !              discriminate between relative values, or absolute values that
1738 !              agree within more than 6 digits
1739                if( abs(results_out(1)%occ(iban,iimage)-results_out(idtset)%occ(iban,iimage)) > tol8) multi=1
1740              end do
1741            end if
1742          end do
1743        end if
1744      end do
1745    end if
1746  end if
1747 
1748 ! write(std_out,*)' prtocc : 4, multi= ',multi
1749 
1750 !At this stage, if multi==1, the occ must be printed
1751 !if multi==0, then it might be that we have the default values.
1752 !Since the default is all zeros, it only happens when iscf=-2
1753 !Also initialize the number of a idtset that can be used as generic
1754 !(this might not be the case for idtset=1 !)
1755 
1756  generic=0
1757  print=0
1758  do idtset=1,ndtset_alloc
1759    if(dtsets(idtset)%iscf/=-2)then
1760      print=1
1761      generic=idtset
1762    end if
1763  end do
1764 
1765 ! write(std_out,*)' prtocc : 5, print= ',print
1766 
1767 !Now, print occ in the generic occupation-number set case (occ is independent of the dtset).
1768  if(print==1 .and. multi==0)then
1769 !  Might restrict the number of k points to be printed
1770    tnkpt=0
1771    nkpt_eff=dtsets(1)%nkpt
1772    if(prtvol_glob==0 .and. nkpt_eff>nkpt_max)then
1773      nkpt_eff=nkpt_max
1774      tnkpt=1
1775    end if
1776 
1777 ! write(std_out,*)' prtocc : 6, do-loop over iimage '
1778 
1779    do iimage=1,nimagem(generic)
1780      if(iimage==1 .or. test_multiimages(generic) )then
1781        keywd=token//trim(strimg(iimage))
1782 !      The quantity of data to be output vary with occopt
1783        if(dtsets(generic)%occopt>=2)then
1784          iban=1
1785          do isppol=1,dtsets(generic)%nsppol
1786            do ikpt=1,nkpt_eff
1787              ikpsp=ikpt+dtsets(generic)%nkpt*(isppol-1)
1788              nban=dtsets(generic)%nband(ikpsp)
1789              if(ikpsp==1)then
1790                write(iout, '(1x,a16,1x,(t22,6f10.6))' )&
1791 &               trim(keywd),results_out(generic)%occ(iban:iban+nban-1,iimage)
1792              else
1793                write(iout, '((t22,6f10.6))' )results_out(generic)%occ(iban:iban+nban-1,iimage)
1794              end if
1795              iban=iban+nban
1796            end do
1797            if(tnkpt==1) write(iout,'(23x,a)' ) 'prtocc : prtvol=0, do not print more k-points.'
1798          end do
1799        else
1800 !        The number of bands is identical for all k points and spin
1801          nban=dtsets(generic)%nband(1)
1802          write(iout, '(1x,a16,1x,(t22,6f10.6))' )trim(keywd),results_out(generic)%occ(1:nban,iimage)
1803 !        if occopt==1, the occ might differ with the spin
1804          if(dtsets(generic)%nsppol/=1)then
1805            write(iout,'((t22,6f10.6))')results_out(generic)%occ(nban*dtsets(generic)%nkpt+1:&
1806 &           nban*dtsets(generic)%nkpt+nban,iimage)
1807          end if
1808        end if
1809      end if
1810    end do
1811  end if
1812 
1813 ! write(std_out,*)' prtocc : 7, finished do-loop over iimage '
1814 
1815 !Now, print occ in the other cases (occ depends on the dataset)
1816  if(print==1 .and. multi==1)then
1817    do idtset=1,ndtset_alloc
1818 !    Might restrict the number of k points to be printed
1819      tnkpt=0
1820      nkpt_eff=dtsets(idtset)%nkpt
1821      if(prtvol_glob==0 .and. nkpt_eff>nkpt_max)then
1822        nkpt_eff=nkpt_max
1823        tnkpt=1
1824      end if
1825      if(dtsets(idtset)%iscf/=-2)then
1826        jdtset=jdtset_(idtset)
1827        call appdig(jdtset,'',appen)
1828        do iimage=1,nimagem(idtset)
1829          if(iimage==1 .or. test_multiimages(idtset) )then
1830            keywd=trim(token)//trim(strimg(iimage))
1831 !          The quantity of data to be output vary with occopt
1832            if(dtsets(idtset)%occopt>=2)then
1833              iban=1
1834              do isppol=1,dtsets(idtset)%nsppol
1835                do ikpt=1,nkpt_eff
1836                  ikpsp=ikpt+dtsets(idtset)%nkpt*(isppol-1)
1837                  nban=dtsets(idtset)%nband(ikpsp)
1838                  if(ikpsp==1)then
1839                    write(iout, '(1x,a16,a,1x,(t22,6f10.6))' )&
1840 &                   trim(keywd),appen,results_out(idtset)%occ(iban:iban+nban-1,iimage)
1841                  else
1842                    write(iout, '((t22,6f10.6))' )results_out(idtset)%occ(iban:iban+nban-1,iimage)
1843                  end if
1844                  iban=iban+nban
1845                end do
1846                if(tnkpt==1) write(iout,'(23x,a)' ) 'prtocc : prtvol=0, do not print more k-points.'
1847              end do
1848            else
1849 !            The number of bands is identical for all k points and spin
1850              nban=dtsets(idtset)%nband(1)
1851              write(iout, '(1x,a16,a,1x,(t22,6f10.6))' )&
1852 &             trim(keywd),appen,results_out(idtset)%occ(1:nban,iimage)
1853 !            if occopt==1, the occ might differ with the spin
1854              if(dtsets(idtset)%nsppol/=1)then
1855                write(iout, '((t22,6f10.6))' ) &
1856 &               results_out(idtset)%occ(nban*dtsets(idtset)%nkpt+1:nban*dtsets(idtset)%nkpt+nban,iimage)
1857              end if
1858            end if
1859          end if
1860        enddo
1861      end if
1862 !    Endloop on idtset
1863    end do
1864  end if
1865 
1866  ABI_FREE(test_multiimages)
1867 
1868 end subroutine prtocc