TABLE OF CONTENTS


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.

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 .

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

SIDE EFFECTS

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

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