TABLE OF CONTENTS


ABINIT/m_outvar_a_h [ Modules ]

[ Top ] [ Modules ]

NAME

  m_outvar_a_h

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 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 .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_outvar_a_h
27 
28  use defs_basis
29  use defs_abitypes
30  use m_abicore
31  use m_results_out
32 
33  use m_parser,  only : prttagm, prttagm_images
34 
35  implicit none
36 
37  private

ABINIT/outvar_a_h [ Functions ]

[ Top ] [ Functions ]

NAME

 outvar_a_h

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
  dmatpuflag=flag controlling the use of an initial density matrix in PAW+U (max. value over datasets)
  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
         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
         nkpt       =maximal value of input nkpt for all the datasets
         nkptgw     =maximal value of input nkptgw 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.
   for different datasets
  results_out(0:ndtset_alloc)=<type results_out_type>contains the results
   needed for outvars, including evolving variables

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.

PARENTS

      outvars

CHILDREN

      prttagm,prttagm_images

SOURCE

 108 subroutine outvar_a_h (choice,dmatpuflag,dtsets,iout,&
 109 & jdtset_,marr,multivals,mxvals,ncid,ndtset,ndtset_alloc,&
 110 & results_out,strimg)
 111 
 112 
 113 !This section has been created automatically by the script Abilint (TD).
 114 !Do not modify the following lines by hand.
 115 #undef ABI_FUNC
 116 #define ABI_FUNC 'outvar_a_h'
 117 !End of the abilint section
 118 
 119  implicit none
 120 
 121 !Arguments ------------------------------------
 122 !scalars
 123  integer,intent(in) :: choice,dmatpuflag,iout,marr,ndtset
 124  integer,intent(in) :: ndtset_alloc,ncid
 125 !arrays
 126  integer,intent(in) :: jdtset_(0:ndtset_alloc)
 127  type(ab_dimensions),intent(in) :: multivals,mxvals
 128  type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc)
 129  type(results_out_type),intent(in) :: results_out(0:ndtset_alloc)
 130  character(len=8),intent(in) :: strimg(mxvals%nimage)
 131 
 132 !Local variables-------------------------------
 133 !scalars
 134  integer,parameter :: nkpt_max=50
 135  integer :: defo,idtset,ii,iimage,ga_n_rules,nn
 136  integer :: lpawu1,narr,mxnsp
 137  integer :: natom,nimfrqs,nimage
 138  integer :: ntypalch,ntypat,size1,size2,tmpimg0
 139  logical :: compute_static_images
 140  real(dp) :: cpus
 141  character(len=1) :: firstchar_fftalg,firstchar_gpu
 142  character(len=14) :: str_hyb
 143 !arrays
 144  integer,allocatable :: narrm(:)
 145  integer,allocatable :: nimagem(:),prtimg(:,:)
 146  integer,allocatable :: intarr(:,:)
 147  real(dp),allocatable :: dprarr(:,:),dprarr_images(:,:,:)
 148 
 149 ! *************************************************************************
 150 
 151 !###########################################################
 152 !### 01. Initial allocations and initialisations.
 153 
 154  ABI_ALLOCATE(dprarr,(marr,0:ndtset_alloc))
 155  ABI_ALLOCATE(dprarr_images,(marr,mxvals%nimage,0:ndtset_alloc))
 156  ABI_ALLOCATE(intarr,(marr,0:ndtset_alloc))
 157  ABI_ALLOCATE(narrm,(0:ndtset_alloc))
 158  ABI_ALLOCATE(nimagem,(0:ndtset_alloc))
 159  ABI_ALLOCATE(prtimg,(mxvals%nimage,0:ndtset_alloc))
 160 
 161  do idtset=0,ndtset_alloc
 162    nimagem(idtset)=dtsets(idtset)%nimage
 163  end do
 164 
 165  firstchar_gpu=' ';if (maxval(dtsets(1:ndtset_alloc)%use_gpu_cuda)>0) firstchar_gpu='-'
 166 
 167 !if(multivals%ga_n_rules==0)ga_n_rules=dtsets(1)%ga_n_rules
 168  ga_n_rules=dtsets(1)%ga_n_rules
 169 !if(multivals%natom==0)natom=dtsets(1)%natom
 170  natom=dtsets(1)%natom
 171 !if(multivals%nimage==0)nimage=dtsets(1)%nimage
 172  nimage=dtsets(1)%nimage
 173 
 174  nimfrqs=dtsets(1)%cd_customnimfrqs
 175 !if(multivals%ntypalch==0)ntypalch=dtsets(1)%ntypalch
 176  ntypalch=dtsets(1)%ntypalch
 177 !if(multivals%ntypat==0)ntypat=dtsets(1)%ntypat
 178  ntypat=dtsets(1)%ntypat
 179 
 180 !write(ab_out,*)' outvar_a_h : A '
 181 !call flush(ab_out)
 182 !###########################################################
 183 !### 03. Print all the input variables (A)
 184 !##
 185 
 186  intarr(1,:)=dtsets(:)%iomode
 187  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'iomode','INT',0,firstchar="-")
 188 
 189  intarr(1,:)=dtsets(:)%accuracy
 190  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'accuracy','INT',0)
 191 
 192 !acell
 193  prtimg(:,:)=1
 194  do idtset=0,ndtset_alloc
 195    narrm(idtset)=3
 196    do iimage=1,nimagem(idtset)
 197      if (narrm(idtset)>0) then
 198        dprarr_images(1:narrm(idtset),iimage,idtset)=results_out(idtset)%acell(1:3,iimage)
 199      end if
 200    end do
 201  end do
 202  call prttagm_images(dprarr_images,iout,jdtset_,2,marr,narrm,ncid,ndtset_alloc,'acell','LEN',&
 203 & mxvals%nimage,nimagem,ndtset,prtimg,strimg)
 204 
 205 !adpimd and adpimd_gamma
 206  intarr(1,:)=dtsets(:)%adpimd
 207  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'adpimd','INT',0)
 208 
 209  dprarr(1,:)=dtsets(:)%adpimd_gamma
 210  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'adpimd_gamma','DPR',0)
 211 
 212 !algalch
 213  narr=ntypalch                      ! default size for all datasets
 214  do idtset=0,ndtset_alloc       ! especific size for each dataset
 215    narrm(idtset)=dtsets(idtset)%ntypalch
 216    if(idtset==0)narrm(idtset)=mxvals%ntypat
 217    intarr(1:narrm(idtset),idtset)=dtsets(idtset)%algalch(1:narrm(idtset))
 218  end do
 219  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 220 & narrm,ncid,ndtset_alloc,'algalch','INT',multivals%ntypalch)
 221 
 222 !amu
 223  prtimg(:,:)=1
 224  do idtset=0,ndtset_alloc
 225    if(idtset/=0)then
 226      size1=dtsets(idtset)%ntypat
 227    else
 228      size1=mxvals%ntypat
 229    end if
 230    narrm(idtset)=size1
 231    do iimage=1,nimagem(idtset)
 232      if (narrm(idtset)>0) then
 233        dprarr_images(1:narrm(idtset),iimage,idtset)=results_out(idtset)%amu(1:size1,iimage)
 234      end if
 235    end do
 236  end do
 237  call prttagm_images(dprarr_images,iout,jdtset_,1,marr,narrm,ncid,ndtset_alloc,'amu','DPR',&
 238 & mxvals%nimage,nimagem,ndtset,prtimg,strimg,forceprint=2)
 239 
 240  intarr(1,:)=dtsets(:)%asr
 241  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'asr','INT',0)
 242 
 243 !atvshift
 244  if(mxvals%natpawu>0)then
 245    narr=dtsets(1)%natvshift*dtsets(1)%nsppol*mxvals%natom ! default size for all datasets
 246    do idtset=0,ndtset_alloc       ! especific size for each dataset
 247      if(idtset/=0)then
 248        narrm(idtset)=dtsets(idtset)%natvshift*dtsets(idtset)%nsppol*mxvals%natom
 249        if(narrm(idtset)/=0)&
 250 &       dprarr(1:narrm(idtset),idtset)=&
 251 &       reshape(dtsets(idtset)%atvshift(1:dtsets(idtset)%natvshift,&
 252 &       1:dtsets(idtset)%nsppol,1:mxvals%natom),&
 253 &       (/ narrm(idtset) /) )
 254      else
 255        narrm(idtset)=mxvals%natvshift*mxvals%nsppol*mxvals%natom
 256        if(narrm(idtset)/=0)&
 257 &       dprarr(1:narrm(idtset),idtset)=&
 258 &       reshape(dtsets(idtset)%atvshift(1:mxvals%natvshift,&
 259 &       1:mxvals%nsppol,1:mxvals%natom),&
 260 &       (/ narrm(idtset) /) )
 261      end if
 262    end do
 263    call prttagm(dprarr,intarr,iout,jdtset_,5,marr,narr,&
 264 &   narrm,ncid,ndtset_alloc,'atvshift','DPR',&
 265 &   multivals%natvshift+multivals%nsppol+multivals%natom)
 266  end if
 267 
 268  intarr(1,:)=dtsets(:)%autoparal
 269  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'autoparal','INT',0)
 270 
 271  intarr(1,:)=dtsets(:)%auxc_ixc
 272  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'auxc_ixc','INT',0)
 273 
 274  dprarr(1,:)=dtsets(:)%auxc_scal
 275  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'auxc_scal','DPR',0)
 276 
 277  intarr(1,:)=dtsets(:)%awtr
 278  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'awtr','INT',0)
 279 
 280 !write(ab_out,*)' outvar_a_h : B '
 281 !call flush(ab_out)
 282 !###########################################################
 283 !### 03. Print all the input variables (B)
 284 !##
 285 
 286  intarr(1,:)=dtsets(:)%bandpp
 287  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bandpp','INT',0)
 288 
 289  intarr(1,:)=dtsets(:)%bdberry(1)
 290  intarr(2,:)=dtsets(:)%bdberry(2)
 291  intarr(3,:)=dtsets(:)%bdberry(3)
 292  intarr(4,:)=dtsets(:)%bdberry(4)
 293  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,4,narrm,ncid,ndtset_alloc,'bdberry','INT',0)
 294 
 295  intarr(1,:)=dtsets(:)%bdeigrf
 296  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bdeigrf','INT',0)
 297 
 298 !bdgw
 299  narr=2*dtsets(1)%nkptgw*dtsets(1)%nsppol ! default size for all datasets
 300  do idtset=0,ndtset_alloc        ! especific size for each dataset
 301    if(idtset/=0)then
 302      narrm(idtset)=2*dtsets(idtset)%nkptgw*dtsets(idtset)%nsppol
 303      if (narrm(idtset)>0)&
 304 &     intarr(1:narrm(idtset),idtset)=&
 305 &     reshape(dtsets(idtset)%bdgw(1:2,1:dtsets(idtset)%nkptgw,1:dtsets(idtset)%nsppol),&
 306 &     (/ narrm(idtset) /) )
 307    else
 308      narrm(idtset)=2*mxvals%nkptgw*mxvals%nsppol
 309      if (narrm(idtset)>0)&
 310 &     intarr(1:narrm(idtset),idtset)=&
 311 &     reshape(dtsets(idtset)%bdgw(1:2,1:mxvals%nkptgw,1:mxvals%nsppol),&
 312 &     (/ narrm(idtset) /) )
 313    end if
 314  end do
 315  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,&
 316 & narrm,ncid,ndtset_alloc,'bdgw','INT',multivals%nkptgw+multivals%nsppol)
 317 
 318  intarr(1,:)=dtsets(:)%berryopt
 319  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'berryopt','INT',0)
 320 
 321  intarr(1,:)=dtsets(:)%berrysav
 322  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'berrysav','INT',0)
 323 
 324  intarr(1,:)=dtsets(:)%berrystep
 325  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'berrystep','INT',0)
 326 
 327  dprarr(1,:)=dtsets(:)%bfield(1)
 328  dprarr(2,:)=dtsets(:)%bfield(2)
 329  dprarr(3,:)=dtsets(:)%bfield(3)
 330  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'bfield','DPR',0)
 331 
 332  dprarr(1,:)=dtsets(:)%bmass
 333  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'bmass','DPR',0)
 334 
 335  dprarr(1,:)=dtsets(:)%boxcenter(1)
 336  dprarr(2,:)=dtsets(:)%boxcenter(2)
 337  dprarr(3,:)=dtsets(:)%boxcenter(3)
 338  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'boxcenter','DPR',0)
 339 
 340  dprarr(1,:)=dtsets(:)%boxcutmin
 341  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'boxcutmin','DPR',0)
 342 
 343  intarr(1,:)=dtsets(:)%bs_algorithm
 344  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_algorithm','INT',0)
 345 
 346  intarr(1,:)=dtsets(:)%bs_calctype
 347  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_calctype','INT',0)
 348 
 349  intarr(1,:)=dtsets(:)%bs_coulomb_term
 350  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_coulomb_term','INT',0)
 351 
 352  intarr(1,:)=dtsets(:)%bs_coupling
 353  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_coupling','INT',0)
 354 
 355  do idtset=0,ndtset_alloc
 356    dprarr(1:2,idtset)=dtsets(idtset)%bs_eh_cutoff(1:2)
 357  end do
 358  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'bs_eh_cutoff','ENE',0)
 359 
 360  intarr(1,:)=dtsets(:)%bs_exchange_term
 361  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_exchange_term','INT',0)
 362 
 363  do idtset=0,ndtset_alloc
 364    dprarr(1:3,idtset)=dtsets(idtset)%bs_freq_mesh(1:3)
 365  end do
 366  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'bs_freq_mesh','ENE',0)
 367 
 368  intarr(1,:)=dtsets(:)%bs_haydock_niter
 369  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_haydock_niter','INT',0)
 370 
 371  do idtset=0,ndtset_alloc
 372    dprarr(1:2,idtset)=dtsets(idtset)%bs_haydock_tol(:)
 373  end do
 374  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'bs_haydock_tol','DPR',0)
 375 
 376  intarr(1,:)=dtsets(:)%bs_hayd_term
 377  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_hayd_term','INT',0)
 378 
 379  do idtset=0,ndtset_alloc
 380    intarr(1:3,idtset)=dtsets(idtset)%bs_interp_kmult(1:3)
 381  end do
 382  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'bs_interp_kmult','INT',0)
 383 
 384  !dprarr(1,:)=dtsets(:)%bs_interp_m3_width
 385  !call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_interp_m3_width','DPR',0)
 386 
 387  intarr(1,:)=dtsets(:)%bs_interp_method
 388  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_interp_method','INT',0)
 389 
 390  intarr(1,:)=dtsets(:)%bs_interp_mode
 391  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_interp_mode','INT',0)
 392 
 393  intarr(1,:)=dtsets(:)%bs_interp_prep
 394  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_interp_prep','INT',0)
 395 
 396  intarr(1,:)=dtsets(:)%bs_interp_rl_nb
 397  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_interp_rl_nb','INT',0)
 398 
 399 !bs_loband
 400  narr=dtsets(1)%nsppol ! default size for all datasets
 401  intarr = 0
 402  do idtset=0,ndtset_alloc        ! especific size for each dataset
 403    if(idtset/=0)then
 404      narrm(idtset)=dtsets(idtset)%nsppol
 405    else
 406      narrm(idtset)=mxvals%nsppol
 407    end if
 408    intarr(1:narrm(idtset),idtset)=dtsets(idtset)%bs_loband(1:narrm(idtset))
 409  end do
 410 
 411  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,narrm,ncid,ndtset_alloc,'bs_loband','INT',multivals%nsppol)
 412 
 413  intarr(1,:)=dtsets(:)%bs_nstates
 414  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'bs_nstates','INT',0)
 415 
 416  intarr(1,:)=dtsets(:)%builtintest
 417  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'builtintest','INT',0)
 418 
 419  dprarr(1,:)=dtsets(:)%bxctmindg
 420  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'bxctmindg','DPR',0)
 421 
 422 !###########################################################
 423 !### 03. Print all the input variables (C)
 424 !##
 425 
 426  if (ANY(dtsets(:)%cd_customnimfrqs/=0)) then
 427    intarr(1,:)=dtsets(:)%cd_customnimfrqs
 428    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'cd_customnimfrqs','INT',0)
 429  end if
 430 
 431  intarr(1,:)=dtsets(:)%cd_frqim_method
 432  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'cd_frqim_method','INT',0)
 433 
 434  intarr(1,:)=dtsets(:)%cd_full_grid
 435  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'cd_full_grid','INT',0)
 436 
 437  dprarr(1,:)=dtsets(:)%cd_halfway_freq
 438  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'cd_halfway_freq','ENE',0)
 439 
 440 !cd_imfrqs
 441  narr=mxvals%nimfrqs            ! default size for all datasets
 442  do idtset=0,ndtset_alloc       ! specific size for each dataset
 443    narrm(idtset)=dtsets(idtset)%cd_customnimfrqs
 444    if(idtset==0)narrm(idtset)=mxvals%nimfrqs
 445    if (narrm(idtset)>0) then
 446      dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%cd_imfrqs(1:narrm(idtset))
 447    end if
 448  end do
 449  call prttagm(dprarr,intarr,iout,jdtset_,6,marr,narr,narrm,ncid,ndtset_alloc,'cd_imfrqs','ENE',multivals%nimfrqs)
 450 
 451  dprarr(1,:)=dtsets(:)%cd_max_freq
 452  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'cd_max_freq','ENE',0)
 453 
 454  if (ANY(dtsets(:)%cd_subset_freq(1)/=0)) then
 455    intarr(1,:)=dtsets(:)%cd_subset_freq(1)
 456    intarr(2,:)=dtsets(:)%cd_subset_freq(2)
 457    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'cd_subset_freq','INT',0)
 458  end if
 459 
 460  dprarr(1,:)=dtsets(:)%charge
 461  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'charge','DPR',0)
 462 
 463 !chempot
 464  narr=3*mxvals%nzchempot*mxvals%ntypat ! default size for all datasets
 465  if(narr/=0)then
 466    do idtset=0,ndtset_alloc       ! specific size for each dataset
 467      if(idtset/=0)then
 468        narrm(idtset)=3*dtsets(idtset)%nzchempot*dtsets(idtset)%ntypat
 469        if(narrm(idtset)/=0)&
 470 &       dprarr(1:narrm(idtset),idtset)=&
 471 &       reshape(dtsets(idtset)%chempot(1:3,1:dtsets(idtset)%nzchempot,&
 472 &       1:dtsets(idtset)%ntypat),&
 473 &       (/ narrm(idtset) /) )
 474      else
 475        narrm(idtset)=3*mxvals%nzchempot*mxvals%ntypat
 476        if(narrm(idtset)/=0)&
 477 &       dprarr(1:narrm(idtset),idtset)=&
 478 &       reshape(dtsets(idtset)%chempot(1:3,1:mxvals%nzchempot,1:mxvals%ntypat),&
 479 &       (/ narrm(idtset) /) )
 480      end if
 481    end do
 482    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 483 &   narrm,ncid,ndtset_alloc,'chempot','DPR',1)
 484  end if
 485 
 486  intarr(1,:)=dtsets(:)%chkdilatmx
 487  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'chkdilatmx','INT',0)
 488 
 489  intarr(1,:)=dtsets(:)%chkexit
 490  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'chkexit','INT',0)
 491 
 492  intarr(1,:)=dtsets(:)%chkprim
 493  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'chkprim','INT',0)
 494 
 495  intarr(1,:)=dtsets(:)%chksymbreak
 496  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'chksymbreak','INT',0)
 497 
 498  intarr(1,:)=dtsets(:)%chneut
 499  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'chneut','INT',0)
 500 
 501  intarr(1,:)=dtsets(:)%cineb_start
 502  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'cineb_start','INT',0)
 503 
 504  if(dtsets(1)%cpus>one)then
 505    cpus=dtsets(1)%cpus
 506    write(iout,'(1x,a16,1x,1p,t22,g10.2,t25,a)') 'cpus',cpus,'(seconds)'
 507    write(iout,'(1x,a16,1x,1p,t22,g10.2,t25,a)') 'cpum',cpus/60.0_dp,'(minutes)'
 508    write(iout,'(1x,a16,1x,1p,t22,g10.2,t25,a)') 'cpuh',cpus/3600.0_dp,'(hours)'
 509  end if
 510 
 511  do idtset=0, ndtset_alloc
 512    do ii = 1, ntypat
 513      dprarr(ii,idtset) = dtsets(idtset)%corecs(ii)
 514    end do ! end loop over ntypat
 515  end do ! end loop over datasets
 516  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,ntypat,narrm,ncid,ndtset_alloc,'corecs','DPR',0)
 517 
 518 !###########################################################
 519 !### 03. Print all the input variables (D)
 520 !##
 521 
 522  dprarr(1,:)=dtsets(:)%ddamp
 523  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ddamp','DPR',0)
 524 
 525  intarr(1,:)=dtsets(:)%delayperm
 526  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'delayperm','INT',0)
 527 
 528  intarr(1,:)=dtsets(:)%densfor_pred
 529  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'densfor_pred','INT',0)
 530 
 531 !densty
 532  narr=mxvals%ntypat              ! default size for all datasets
 533  do idtset=0,ndtset_alloc        ! especific size for each dataset
 534    narrm(idtset)=dtsets(idtset)%ntypat
 535    if(idtset==0)narrm(idtset)=mxvals%ntypat
 536 !  Only one component of densty is used until now
 537    dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%densty(1:narrm(idtset),1)
 538  end do
 539  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'densty','DPR',multivals%ntypat)
 540 
 541  dprarr(1,:)=dtsets(:)%dfield(1)
 542  dprarr(2,:)=dtsets(:)%dfield(2)
 543  dprarr(3,:)=dtsets(:)%dfield(3)
 544  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'dfield','DPR',0)
 545 
 546  dprarr(1,:)=dtsets(:)%dfpt_sciss
 547  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'dfpt_sciss','ENE',0)
 548 
 549 
 550  dprarr(1,:)=dtsets(:)%diecut
 551  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'diecut','ENE',0)
 552 
 553  dprarr(1,:)=dtsets(:)%diegap
 554  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'diegap','ENE',0)
 555 
 556  dprarr(1,:)=dtsets(:)%dielam
 557  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'dielam','DPR',0)
 558 
 559  dprarr(1,:)=dtsets(:)%dielng
 560  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'dielng','LEN',0)
 561 
 562  dprarr(1,:)=dtsets(:)%diemac
 563  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'diemac','DPR',0)
 564 
 565  dprarr(1,:)=dtsets(:)%diemix
 566  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'diemix','DPR',0)
 567 
 568  if (any(dtsets(1:ndtset_alloc)%diemixmag/=dtsets(1:ndtset_alloc)%diemix)) then
 569    dprarr(1,:)=dtsets(:)%diemixmag
 570    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'diemixmag','DPR',0)
 571  end if
 572 
 573  intarr(1,:)=dtsets(:)%diismemory
 574  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'diismemory','INT',0)
 575 
 576  dprarr(1,:)=dtsets(:)%dilatmx
 577  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'dilatmx','DPR',0)
 578 
 579  intarr(1,:)=dtsets(:)%dipdip
 580  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dipdip','INT',0)
 581 
 582 !dmatpawu
 583  if (dmatpuflag==1.and.mxvals%natpawu>0) then
 584    prtimg(:,:)=1
 585    do idtset=0,ndtset_alloc
 586      mxnsp=max(dtsets(idtset)%nsppol,dtsets(idtset)%nspinor)
 587      lpawu1=maxval(dtsets(idtset)%lpawu(:))
 588      narrm(idtset)=((2*lpawu1+1)**2)*mxnsp*dtsets(idtset)%natpawu
 589      do iimage=1,nimagem(idtset)
 590        if (narrm(idtset)>0) then
 591          dprarr_images(1:narrm(idtset),iimage,idtset)= &
 592 &         reshape(dtsets(idtset)%dmatpawu(&
 593 &         1:2*lpawu1+1,1:2*lpawu1+1,1:mxnsp,1:dtsets(idtset)%natpawu,iimage),&
 594 &         (/narrm(idtset)/))
 595        end if
 596      end do
 597    end do
 598    call prttagm_images(dprarr_images,iout,jdtset_,5,marr,narrm,&
 599 &   ncid,ndtset_alloc,'dmatpawu','DPR',mxvals%nimage,nimagem,ndtset,prtimg,strimg)
 600  end if
 601 
 602  intarr(1,:)=dtsets(:)%dmatpuopt
 603  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmatpuopt','INT',0)
 604 
 605  intarr(1,:)=dtsets(:)%dmatudiag
 606  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmatudiag','INT',0)
 607 
 608  intarr(1,:)=dtsets(:)%dmftbandf
 609  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmftbandf','INT',0)
 610 
 611  intarr(1,:)=dtsets(:)%dmftbandi
 612  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmftbandi','INT',0)
 613 
 614  intarr(1,:)=dtsets(:)%dmftcheck
 615  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmftcheck','INT',0)
 616 
 617  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_dc','INT',0)
 618 
 619  intarr(1,:)=dtsets(:)%dmft_iter
 620  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_iter','INT',0)
 621 
 622  dprarr(1,:)=dtsets(:)%dmft_mxsf
 623  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_mxsf','DPR',0)
 624 
 625  intarr(1,:)=dtsets(:)%dmft_nwli
 626  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_nwli','INT',0)
 627 
 628  intarr(1,:)=dtsets(:)%dmft_nwlo
 629  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_nwlo','INT',0)
 630 
 631  intarr(1,:)=dtsets(:)%dmft_read_occnd
 632  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_read_occnd','INT',0)
 633 
 634  intarr(1,:)=dtsets(:)%dmft_rslf
 635  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_rslf','INT',0)
 636 
 637  intarr(1,:)=dtsets(:)%dmft_solv
 638  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_solv','INT',0)
 639 
 640  dprarr(1,:)=dtsets(:)%dmft_tolfreq
 641  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_tolfreq','DPR',0)
 642 
 643  dprarr(1,:)=dtsets(:)%dmft_tollc
 644  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'dmft_tollc','DPR',0)
 645 
 646  dprarr(1,:)=dtsets(:)%dosdeltae
 647  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'dosdeltae','ENE',0)
 648 
 649  dprarr(1,:)=dtsets(:)%dtion
 650  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'dtion','DPR',0,forceprint=2)
 651 
 652 !dynimage
 653  intarr(1:marr,0)=1                 ! default value
 654  narr=nimage                        ! default size for all datasets
 655  do idtset=1,ndtset_alloc           ! specific size and array for each dataset
 656    narrm(idtset)=dtsets(idtset)%nimage
 657    intarr(1:narrm(idtset),idtset)=dtsets(idtset)%dynimage(1:narrm(idtset))
 658  end do
 659  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 660 & narrm,ncid,ndtset_alloc,'dynimage','INT',multivals%nimage)
 661 
 662 !Variables for nonlinear response
 663  intarr(1,:)=dtsets(:)%d3e_pert1_atpol(1)
 664  intarr(2,:)=dtsets(:)%d3e_pert1_atpol(2)
 665  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'d3e_pert1_atpol','INT',0)
 666 
 667  intarr(1,:)=dtsets(:)%d3e_pert1_dir(1)
 668  intarr(2,:)=dtsets(:)%d3e_pert1_dir(2)
 669  intarr(3,:)=dtsets(:)%d3e_pert1_dir(3)
 670  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'d3e_pert1_dir','INT',0)
 671 
 672  intarr(1,:)=dtsets(:)%d3e_pert1_elfd
 673  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'d3e_pert1_elfd','INT',0)
 674 
 675  intarr(1,:)=dtsets(:)%d3e_pert1_phon
 676  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'d3e_pert1_phon','INT',0)
 677 
 678  intarr(1,:)=dtsets(:)%d3e_pert2_atpol(1)
 679  intarr(2,:)=dtsets(:)%d3e_pert2_atpol(2)
 680  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'d3e_pert2_atpol','INT',0)
 681 
 682  intarr(1,:)=dtsets(:)%d3e_pert2_dir(1)
 683  intarr(2,:)=dtsets(:)%d3e_pert2_dir(2)
 684  intarr(3,:)=dtsets(:)%d3e_pert2_dir(3)
 685  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'d3e_pert2_dir','INT',0)
 686 
 687  intarr(1,:)=dtsets(:)%d3e_pert2_elfd
 688  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'d3e_pert2_elfd','INT',0)
 689 
 690  intarr(1,:)=dtsets(:)%d3e_pert2_phon
 691  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'d3e_pert2_phon','INT',0)
 692 
 693  intarr(1,:)=dtsets(:)%d3e_pert3_atpol(1)
 694  intarr(2,:)=dtsets(:)%d3e_pert3_atpol(2)
 695  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,2,narrm,ncid,ndtset_alloc,'d3e_pert3_atpol','INT',0)
 696 
 697  intarr(1,:)=dtsets(:)%d3e_pert3_dir(1)
 698  intarr(2,:)=dtsets(:)%d3e_pert3_dir(2)
 699  intarr(3,:)=dtsets(:)%d3e_pert3_dir(3)
 700  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'d3e_pert3_dir','INT',0)
 701 
 702  intarr(1,:)=dtsets(:)%d3e_pert3_elfd
 703  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'d3e_pert3_elfd','INT',0)
 704 
 705  intarr(1,:)=dtsets(:)%d3e_pert3_phon
 706  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'d3e_pert3_phon','INT',0)
 707 
 708 !###########################################################
 709 !### 03. Print all the input variables (E)
 710 !##
 711 
 712  dprarr(1,:)=dtsets(:)%ecut
 713  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ecut','ENE',0)
 714 
 715  dprarr(1,:)=dtsets(:)%ecuteps
 716  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ecuteps','ENE',0)
 717 
 718  dprarr(1,:)=dtsets(:)%ecutsigx
 719  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ecutsigx','ENE',0)
 720 
 721  dprarr(1,:)=dtsets(:)%ecutsm
 722  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ecutsm','ENE',0)
 723 
 724  dprarr(1,:)=dtsets(:)%ecutwfn
 725  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ecutwfn','ENE',0)
 726 
 727  dprarr(1,:)=dtsets(:)%effmass_free
 728  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'effmass_free','DPR',0)
 729 
 730  dprarr(1,:)=dtsets(:)%efield(1)
 731  dprarr(2,:)=dtsets(:)%efield(2)
 732  dprarr(3,:)=dtsets(:)%efield(3)
 733  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'efield','DPR',0)
 734 
 735  nn = size(dtsets(0)%einterp)
 736  do ii=1,nn
 737    dprarr(ii,:)=dtsets(:)%einterp(ii)
 738  end do
 739  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,nn,narrm,ncid,ndtset_alloc,'einterp','DPR',0)
 740 
 741  dprarr(1,:)=dtsets(:)%elph2_imagden
 742  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'elph2_imagden','ENE',0)
 743 
 744  intarr(1,:)=dtsets(:)%enunit
 745  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'enunit','INT',0)
 746 
 747 !eph variables
 748  intarr(1,:)=dtsets(:)%eph_intmeth
 749  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'eph_intmeth','INT',0)
 750 
 751  dprarr(1,:)=dtsets(:)%eph_extrael
 752  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'eph_extrael','DPR',0)
 753 
 754  dprarr(1,:)=dtsets(:)%eph_fermie
 755  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'eph_fermie','ENE',0)
 756 
 757  intarr(1,:)=dtsets(:)%eph_frohlichm
 758  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'eph_frohlichm','INT',0)
 759 
 760  dprarr(1,:)=dtsets(:)%eph_fsmear
 761  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'eph_fsmear','ENE',0)
 762 
 763  dprarr(1,:)=dtsets(:)%eph_fsewin
 764  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'eph_fsewin','ENE',0)
 765 
 766  dprarr(1,:)=dtsets(:)%eph_mustar
 767  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'eph_mustar','DPR',0)
 768 
 769  intarr(1,:)=dtsets(:)%eph_task
 770  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'eph_task','INT',0)
 771 
 772  do idtset=0,ndtset_alloc
 773    intarr(1:3,idtset)=dtsets(idtset)%eph_ngqpt_fine
 774  end do
 775  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'eph_ngqpt_fine','INT',0)
 776 
 777  intarr(1,:)=dtsets(:)%eph_transport
 778  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'eph_transport','INT',0)
 779 
 780  dprarr(1,:)=dtsets(:)%ph_wstep
 781  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_wstep','ENE',0)
 782 
 783  intarr(1,:)=dtsets(:)%ph_intmeth
 784  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_intmeth','INT',0)
 785 
 786  intarr(1,:)=dtsets(:)%ph_nqshift
 787  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_nqshift','INT',0)
 788 
 789  dprarr(1,:)=dtsets(:)%ph_smear
 790  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'ph_smear','ENE',0)
 791 
 792  do idtset=0,ndtset_alloc
 793    intarr(1:3,idtset)=dtsets(idtset)%ddb_ngqpt
 794  end do
 795  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'ddb_ngqpt','INT',0)
 796 
 797  dprarr(1,:)=dtsets(:)%ddb_shiftq(1)
 798  dprarr(2,:)=dtsets(:)%ddb_shiftq(2)
 799  dprarr(3,:)=dtsets(:)%ddb_shiftq(3)
 800  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'ddb_shiftq','DPR',0)
 801 
 802  intarr(1,:)=dtsets(:)%ph_ndivsm
 803  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ph_ndivsm','INT',0)
 804 
 805  intarr(1,:)=dtsets(:)%ph_nqpath
 806  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ph_nqpath','INT',0)
 807 
 808  do idtset=0,ndtset_alloc
 809    intarr(1:3,idtset)=dtsets(idtset)%ph_ngqpt
 810  end do
 811  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,3,narrm,ncid,ndtset_alloc,'ph_ngqpt','INT',0)
 812 !end e-ph variables
 813 
 814  dprarr(1,:)=dtsets(:)%eshift
 815  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'eshift','ENE',0)
 816 
 817  dprarr(1,:)=dtsets(:)%esmear
 818  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'esmear','ENE',0)
 819 
 820 !etotal
 821  if(choice==2)then
 822    prtimg(:,:)=1
 823    do idtset=0,ndtset_alloc       ! especific size for each dataset
 824      compute_static_images=(dtsets(idtset)%istatimg>0)
 825      narrm(idtset)=1
 826 
 827      if(dtsets(idtset)%iscf>=0 .or. dtsets(idtset)%iscf==-3)then
 828        do iimage=1,dtsets(idtset)%nimage
 829          if (narrm(idtset)>0) then
 830            dprarr_images(1:narrm(idtset),iimage,idtset)=&
 831 &           results_out(idtset)%etotal(iimage)
 832          end if
 833          if(.not.(dtsets(idtset)%dynimage(iimage)==1.or.compute_static_images))then
 834            prtimg(iimage,idtset)=0
 835          end if
 836        end do
 837      else
 838        narrm(idtset)=0
 839      end if
 840    end do
 841 !  This is a trick to force printing of etotal even if zero, still not destroying the value of nimagem(0).
 842    tmpimg0=nimagem(0)
 843    nimagem(0)=0
 844    call prttagm_images(dprarr_images,iout,jdtset_,2,marr,narrm,ncid,ndtset_alloc,'etotal','DPR',&
 845 &   mxvals%nimage,nimagem,ndtset,prtimg,strimg)
 846    nimagem(0)=tmpimg0
 847  end if
 848 
 849  dprarr(1,:)=dtsets(:)%exchmix
 850  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'exchmix','DPR',0)
 851 
 852  intarr(1,:)=dtsets(:)%exchn2n3d
 853  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'exchn2n3d','INT',0)
 854 
 855  intarr(1,:)=dtsets(:)%extrapwf
 856  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'extrapwf','INT',0)
 857 
 858 !###########################################################
 859 !### 03. Print all the input variables (F)
 860 !##
 861 
 862 !fcart
 863  if(choice==2)then
 864    prtimg(:,:)=1
 865    do idtset=0,ndtset_alloc       ! especific size for each dataset
 866      compute_static_images=(dtsets(idtset)%istatimg>0)
 867      size2=dtsets(idtset)%natom
 868      if(idtset==0)size2=0
 869      narrm(idtset)=3*size2
 870      if(dtsets(idtset)%iscf>=0 .or. idtset==0)then
 871        do iimage=1,dtsets(idtset)%nimage
 872          if (narrm(idtset)>0) then
 873            dprarr_images(1:narrm(idtset),iimage,idtset)=&
 874 &           reshape(results_out(idtset)%fcart(1:3,1:size2,iimage),&
 875 &           (/ narrm(idtset) /) )
 876          end if
 877          if(.not.(dtsets(idtset)%dynimage(iimage)==1.or.compute_static_images))then
 878            prtimg(iimage,idtset)=0
 879          end if
 880        end do
 881      else
 882        narrm(idtset)=0
 883      end if
 884    end do
 885 !  This is a trick to force printing of fcart even if zero, still not destroying the value of nimagem(0).
 886    tmpimg0=nimagem(0)
 887    nimagem(0)=0
 888    call prttagm_images(dprarr_images,iout,jdtset_,2,&
 889 &   marr,narrm,ncid,ndtset_alloc,'fcart','DPR',&
 890 &   mxvals%nimage,nimagem,ndtset,prtimg,strimg)
 891    nimagem(0)=tmpimg0
 892  end if
 893 
 894  dprarr(1,:)=dtsets(:)%fermie_nest
 895  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'fermie_nest','DPR',0)
 896 
 897  firstchar_fftalg = "_"
 898  intarr(1,:)=dtsets(:)%ngfft(7)
 899  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'fftalg','INT',0,&
 900 & firstchar="-",forceprint=3)
 901 
 902  intarr(1,:)=dtsets(:)%ngfft(8)
 903  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'fftcache','INT',0)
 904 
 905  intarr(1,:)=dtsets(:)%fftgw
 906  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'fftgw','INT',0)
 907 
 908  intarr(1,:)=dtsets(:)%fockoptmix
 909  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'fockoptmix','INT',0)
 910 
 911  dprarr(1,:)=dtsets(:)%focktoldfe
 912  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'focktoldfe','DPR',0)
 913 
 914  intarr(1,:)=dtsets(:)%fockdownsampling(1)
 915  intarr(2,:)=dtsets(:)%fockdownsampling(2)
 916  intarr(3,:)=dtsets(:)%fockdownsampling(3)
 917  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'fockdownsampling','INT',0)
 918 !DEBUG
 919 ! write(std_out,*)' outvar_a_h : echo fockdownsampling ...'
 920 ! write(std_out,*)' dtsets(0)%fockdownsampling(1:3)=',dtsets(0)%fockdownsampling(1:3)
 921 ! write(std_out,*)' dtsets(1)%fockdownsampling(1:3)=',dtsets(1)%fockdownsampling(1:3)
 922 ! write(std_out,*)' dtsets(2)%fockdownsampling(1:3)=',dtsets(2)%fockdownsampling(1:3)
 923 ! write(std_out,*)' dtsets(3)%fockdownsampling(1:3)=',dtsets(3)%fockdownsampling(1:3)
 924 ! call prttagm(dprarr,intarr,std_out,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'fockdownsampling','INT',0)
 925 !ENDDEBUG
 926 
 927  dprarr(1,:)=dtsets(:)%freqim_alpha
 928  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'freqim_alpha','DPR',0)
 929 
 930  dprarr(1,:)=dtsets(:)%freqremax
 931  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'freqremax','ENE',0)
 932 
 933  dprarr(1,:)=dtsets(:)%freqremin
 934  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'freqremin','ENE',0)
 935 
 936  dprarr(1,:)=dtsets(:)%freqspmax
 937  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'freqspmax','ENE',0)
 938 
 939  dprarr(1,:)=dtsets(:)%freqspmin
 940  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'freqspmin','ENE',0)
 941 
 942  dprarr(1,:)=dtsets(:)%friction
 943  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'friction','DPR',0)
 944 
 945  intarr(1,:)=dtsets(:)%frzfermi
 946  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'frzfermi','INT',0)
 947 
 948  dprarr(1,:)=dtsets(:)%fxcartfactor
 949  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'fxcartfactor','DPR',0)
 950 
 951 !f4of2_sla
 952  narr=mxvals%ntypat                    ! default size for all datasets
 953  do idtset=0,ndtset_alloc       ! specific size for each dataset
 954    narrm(idtset)=dtsets(idtset)%ntypat
 955    if(idtset==0)narrm(idtset)=mxvals%ntypat
 956    if (narrm(idtset)>0) then
 957      dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%f4of2_sla(1:narrm(idtset))
 958    end if
 959  end do
 960  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 961 & narrm,ncid,ndtset_alloc,'f4of2_sla','DPR',multivals%ntypat)
 962 
 963 !f6of2_sla
 964  narr=mxvals%ntypat                    ! default size for all datasets
 965  do idtset=0,ndtset_alloc       ! specific size for each dataset
 966    narrm(idtset)=dtsets(idtset)%ntypat
 967    if(idtset==0)narrm(idtset)=mxvals%ntypat
 968    if (narrm(idtset)>0) then
 969      dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%f6of2_sla(1:narrm(idtset))
 970    end if
 971  end do
 972  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 973 & narrm,ncid,ndtset_alloc,'f6of2_sla','DPR',multivals%ntypat)
 974 
 975 !###########################################################
 976 !### 03. Print all the input variables (G)
 977 !##
 978 
 979  intarr(1,:)=dtsets(:)%ga_algor
 980  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ga_algor','INT',0)
 981 
 982  intarr(1,:)=dtsets(:)%ga_fitness
 983  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ga_fitness','INT',0)
 984 
 985  intarr(1,:)=dtsets(:)%ga_n_rules
 986  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ga_n_rules','INT',0)
 987 
 988  dprarr(1,:)=dtsets(:)%ga_opt_percent
 989  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ga_opt_percent','DPR',0)
 990 
 991 !ga_rules
 992  narr=ga_n_rules                    ! default size for all datasets
 993  do idtset=0,ndtset_alloc       ! especific size for each dataset
 994    narrm(idtset)=dtsets(idtset)%ga_n_rules
 995    if(idtset==0)narrm(idtset)=mxvals%ga_n_rules
 996    intarr(1:narrm(idtset),idtset)=dtsets(idtset)%ga_rules(1:narrm(idtset))
 997  end do
 998  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'ga_rules','INT',multivals%ga_n_rules)
 999 
1000  intarr(1,:)=dtsets(:)%getbscoup
1001  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getbscoup','INT',0)
1002 
1003  intarr(1,:)=dtsets(:)%getbseig
1004  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getbseig','INT',0)
1005 
1006  intarr(1,:)=dtsets(:)%getbsreso
1007  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getbsreso','INT',0)
1008 
1009  intarr(1,:)=dtsets(:)%getcell
1010  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getcell','INT',0)
1011 
1012  intarr(1,:)=dtsets(:)%getddb
1013  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getddb','INT',0)
1014 
1015  intarr(1,:)=dtsets(:)%getddk
1016  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getddk','INT',0)
1017 
1018  intarr(1,:)=dtsets(:)%getdelfd
1019  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getdelfd','INT',0)
1020 
1021  intarr(1,:)=dtsets(:)%getdkdk
1022  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getdkdk','INT',0)
1023 
1024  intarr(1,:)=dtsets(:)%getdkde
1025  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getdkde','INT',0)
1026 
1027  intarr(1,:)=dtsets(:)%getden
1028  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getden','INT',0)
1029 
1030  intarr(1,:)=dtsets(:)%getefmas
1031  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getefmas','INT',0)
1032 
1033  intarr(1,:)=dtsets(:)%getgam_eig2nkq
1034  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getgam_eig2nkq','INT',0)
1035 
1036  intarr(1,:)=dtsets(:)%gethaydock
1037  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gethaydock','INT',0)
1038 
1039  intarr(1,:)=dtsets(:)%getocc
1040  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getocc','INT',0)
1041 
1042  intarr(1,:)=dtsets(:)%getpawden
1043  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getpawden','INT',0)
1044 
1045  intarr(1,:)=dtsets(:)%getqps
1046  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getqps','INT',0)
1047 
1048  intarr(1,:)=dtsets(:)%getscr
1049  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getscr','INT',0)
1050 
1051  intarr(1,:)=dtsets(:)%getsuscep
1052  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getsuscep','INT',0)
1053 
1054  intarr(1,:)=dtsets(:)%getvel
1055  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getvel','INT',0)
1056 
1057  intarr(1,:)=dtsets(:)%getwfk
1058  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getwfk','INT',0)
1059 
1060  intarr(1,:)=dtsets(:)%getwfkfine
1061  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getwfkfine','INT',0)
1062 
1063  intarr(1,:)=dtsets(:)%getwfq
1064  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getwfq','INT',0)
1065 
1066  intarr(1,:)=dtsets(:)%getxcart
1067  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getxcart','INT',0)
1068 
1069  intarr(1,:)=dtsets(:)%getxred
1070  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'getxred','INT',0)
1071 
1072  intarr(1,:)=dtsets(:)%get1den
1073  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'get1den','INT',0)
1074 
1075  intarr(1,:)=dtsets(:)%get1wf
1076  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'get1wf','INT',0)
1077 
1078  intarr(1,:)=dtsets(:)%goprecon
1079  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'goprecon','INT',0)
1080 
1081  dprarr(1,:)=dtsets(:)%goprecprm(1)
1082  dprarr(2,:)=dtsets(:)%goprecprm(2)
1083  dprarr(3,:)=dtsets(:)%goprecprm(3)
1084  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'goprecprm','DPR',0)
1085 
1086  intarr(1,:)=dtsets(:)%gpu_devices(1) ; intarr(2,:)=dtsets(:)%gpu_devices(2)
1087  intarr(3,:)=dtsets(:)%gpu_devices(3) ; intarr(4,:)=dtsets(:)%gpu_devices(4)
1088  intarr(5,:)=dtsets(:)%gpu_devices(5)
1089  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,5,narrm,ncid,ndtset_alloc,'gpu_devices','INT',0)
1090 
1091  intarr(1,:)=dtsets(:)%gpu_linalg_limit
1092  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gpu_linalg_limit','INT',0)
1093 
1094  intarr(1,:)=dtsets(:)%gwcalctyp
1095  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwcalctyp','INT',0)
1096 
1097  intarr(1,:)=dtsets(:)%gwcomp
1098  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwcomp','INT',0)
1099 
1100  dprarr(1,:)=dtsets(:)%gwencomp
1101  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwencomp','ENE',0)
1102 
1103  intarr(1,:)=dtsets(:)%gwgamma
1104  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwgamma','INT',0)
1105 
1106  intarr(1,:)=dtsets(:)%gwmem
1107  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwmem','INT',0)
1108 
1109  intarr(1,:)=dtsets(:)%gwpara
1110  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwpara','INT',0, firstchar="-")
1111 
1112  intarr(1,:)=dtsets(:)%gwrpacorr
1113  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gwrpacorr','INT',0)
1114 
1115 !gw_customnfreqsp
1116 !It actually overrides the content of nfreqsp (which is forbidden !) in dtset.
1117 !This is to be cleaned ...
1118  if (ANY(dtsets(:)%gw_customnfreqsp/=0)) then
1119    intarr(1,:)=dtsets(:)%gw_customnfreqsp
1120    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_customnfreqsp','INT',0)
1121  end if
1122 
1123 !gw_freqsp
1124 !This is to be cleaned ... See above ...
1125  narr=mxvals%nfreqsp ! default size for all datasets
1126  do idtset=0,ndtset_alloc       ! specific size for each dataset
1127    dprarr(1:narr,idtset)=zero
1128    narrm(idtset)=dtsets(idtset)%gw_customnfreqsp
1129    if(idtset==0)narrm(idtset)=mxvals%nfreqsp
1130    if (narrm(idtset)>0) then
1131      dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%gw_freqsp(1:narrm(idtset))
1132    end if
1133  end do
1134  call prttagm(dprarr,intarr,iout,jdtset_,6,marr,narr,narrm,ncid,ndtset_alloc,'gw_freqsp','ENE',multivals%nfreqsp)
1135 
1136  intarr(1,:)=dtsets(:)%gw_frqim_inzgrid
1137  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_frqim_inzgrid','INT',0)
1138 
1139  intarr(1,:)=dtsets(:)%gw_frqre_inzgrid
1140  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_frqre_inzgrid','INT',0)
1141 
1142  intarr(1,:)=dtsets(:)%gw_frqre_tangrid
1143  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_frqre_tangrid','INT',0)
1144 
1145  intarr(1,:)=dtsets(:)%gw_invalid_freq
1146  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_invalid_freq','INT',0)
1147 
1148  intarr(1,:)=dtsets(:)%gw_qprange
1149  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_qprange','INT',0)
1150 
1151  intarr(1,:)=dtsets(:)%gw_nqlwl
1152  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_nqlwl','INT',0)
1153 
1154  intarr(1,:)=dtsets(:)%gw_nstep
1155  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_nstep','INT',0)
1156 
1157 !gw_qlwl
1158  narr=3*dtsets(1)%gw_nqlwl ! default size for all datasets
1159  do idtset=0,ndtset_alloc       ! especific size for each dataset
1160    if(idtset/=0)then
1161      narrm(idtset)=3*dtsets(idtset)%gw_nqlwl
1162      if (narrm(idtset)>0)then
1163        dprarr(1:narrm(idtset),idtset)=&
1164 &       reshape(dtsets(idtset)%gw_qlwl(1:3,1:dtsets(idtset)%gw_nqlwl),(/ narrm(idtset) /) )
1165      end if
1166    else
1167      narrm(idtset)=3*mxvals%gw_nqlwl
1168      if (narrm(idtset)>0)then
1169        dprarr(1:narrm(idtset),idtset)=zero
1170        dprarr(1:3,idtset)=(/0.00001_dp, 0.00002_dp, 0.00003_dp/)
1171      end if
1172    end if
1173  end do
1174 
1175  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'gw_qlwl','DPR',multivals%gw_nqlwl)
1176 
1177  intarr(1,:)=dtsets(:)%gw_sctype
1178  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_sctype','INT',0)
1179 
1180 
1181  intarr(1,:)=dtsets(:)%gw_sigxcore
1182  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'gw_sigxcore','INT',0)
1183 
1184  dprarr(1,:)=dtsets(:)%gw_toldfeig
1185  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'gw_toldfeig','ENE',0)
1186 
1187  intarr(1,:)=dtsets(:)%hmcsst
1188  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'hmcsst','INT',0)
1189 
1190  intarr(1,:)=dtsets(:)%hmctt
1191  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'hmctt','INT',0)
1192 
1193 !DEBUG
1194 !write(std_out,*)' hyb_mixing=',dtsets(0:ndtset_alloc)%hyb_mixing
1195 !write(std_out,*)' hyb_mixing_sr=',dtsets(0:ndtset_alloc)%hyb_mixing_sr
1196 !write(std_out,*)' hyb_mixing_dft=',dtsets(0:ndtset_alloc)%hyb_range_dft
1197 !write(std_out,*)' hyb_mixing_fock=',dtsets(0:ndtset_alloc)%hyb_range_fock
1198 !ENDDEBUG
1199 
1200 !Special treatment of the default values for the hybrid functional parameters.
1201  do ii=1,4
1202    if(ii==1)dprarr(1,:)=dtsets(:)%hyb_mixing
1203    if(ii==2)dprarr(1,:)=dtsets(:)%hyb_mixing_sr
1204    if(ii==3)dprarr(1,:)=dtsets(:)%hyb_range_dft
1205    if(ii==4)dprarr(1,:)=dtsets(:)%hyb_range_fock
1206    defo=1
1207    do idtset=1,ndtset_alloc
1208      if(dprarr(1,idtset)<-tol8 .and. abs(dprarr(1,idtset)+999.0_dp)>tol8)defo=0
1209    end do
1210    if(defo==0)then
1211      do idtset=1,ndtset_alloc
1212 !      Change the sign of user defined input value
1213        if(dprarr(1,idtset)<-tol8 .and. abs(dprarr(1,idtset)+999.0_dp)>tol8)then
1214          dprarr(1,idtset)=abs(dprarr(1,idtset))
1215        end if
1216      end do
1217      if(ii==1)str_hyb='hyb_mixing'
1218      if(ii==2)str_hyb='hyb_mixing_sr'
1219      if(ii==3)str_hyb='hyb_range_dft'
1220      if(ii==4)str_hyb='hyb_range_fock'
1221      call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,str_hyb,'DPR',0)
1222    end if
1223  end do
1224 
1225 !###########################################################
1226 !## Deallocation for generic arrays, and for n-z variables
1227 
1228  ABI_DEALLOCATE(dprarr)
1229  ABI_DEALLOCATE(intarr)
1230  ABI_DEALLOCATE(narrm)
1231  ABI_DEALLOCATE(nimagem)
1232  ABI_DEALLOCATE(dprarr_images)
1233  ABI_DEALLOCATE(prtimg)
1234 
1235 end subroutine outvar_a_h