TABLE OF CONTENTS


ABINIT/m_outvar_i_n [ Modules ]

[ Top ] [ Modules ]

NAME

  m_outvar_i_n

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_i_n
27 
28  use defs_basis
29  use defs_abitypes
30  use m_errors
31  use m_abicore
32  use m_errors
33  use m_results_out
34  use m_errors
35 
36  use m_parser,  only : prttagm, prttagm_images
37 
38  implicit none
39 
40  private

ABINIT/outvar_i_n [ Functions ]

[ Top ] [ Functions ]

NAME

 outvar_i_n

FUNCTION

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

INPUTS

  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables
  iout=unit number for echoed output
  jdtset_(0:ndtset_alloc)=actual index of the dataset (equal to dtsets(:)%jdtset)
  marr=maximum number of numbers in an array (might need to be increased ... !)
  multivals= <type ab_dimensions>  either 0 or 1 , depending whether the
     dimension has different values for different datasets
  mxvals= <type ab_dimensions>
     maximum size of some arrays along all datasets, including
         lpawu      =maximal value of input lpawu for all the datasets
         gw_nqlwl   =maximal value of input gw_nqlwl for all the datasets
         mband      =maximum number of bands
         natom      =maximal value of input natom for all the datasets
         natpawu    =maximal value of number of atoms on which +U is applied for all the datasets
         natsph     =maximal value of input natsph for all the datasets
         natvshift  =maximal value of input natvshift for all the datasets
         nconeq     =maximal value of input nconeq for all the datasets
         nimage     =maximal value of input nimage for all the datasets
         nkptgw     =maximal value of input nkptgw for all the datasets
         nkpthf     =maximal value of input nkpthf for all the datasets
         nkpt       =maximal value of input nkpt for all the datasets
         nnos       =maximal value of input nnos for all the datasets
         nqptdm     =maximal value of input nqptdm for all the datasets
         nspinor    =maximal value of input nspinor for all the datasets
         nsppol     =maximal value of input nsppol for all the datasets
         nsym       =maximum number of symmetries
         ntypat     =maximum number of type of atoms
         nzchempot  =maximal value of input nzchempot for all the datasets
  ndtset=number of datasets
  ndtset_alloc=number of datasets, corrected for allocation of at least
      one data set. Use for most dimensioned arrays.
  npsp=number of pseudopotentials
  nqptdm=the number of q vectors provided by the user to calculate DM in GW
  prtvol_glob= if 0, minimal output volume, if 1, no restriction.
  response_(0:ndtset_alloc)= 1 if response variables must be output, 0 otherwise,
   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)

PARENTS

      outvars

CHILDREN

      prttagm,prttagm_images

SOURCE

 111 subroutine outvar_i_n (dtsets,iout,&
 112 & jdtset_,marr,multivals,mxvals,ncid,ndtset,ndtset_alloc,npsp,prtvol_glob,&
 113 & response_,results_out,strimg)
 114 
 115 
 116 !This section has been created automatically by the script Abilint (TD).
 117 !Do not modify the following lines by hand.
 118 #undef ABI_FUNC
 119 #define ABI_FUNC 'outvar_i_n'
 120 !End of the abilint section
 121 
 122  implicit none
 123 
 124 !Arguments ------------------------------------
 125 !scalars
 126  integer,intent(in) :: iout,marr,ndtset
 127  integer,intent(in) :: ndtset_alloc,prtvol_glob,ncid,npsp
 128 !arrays
 129  integer,intent(in) :: jdtset_(0:ndtset_alloc),response_(ndtset_alloc)
 130  type(ab_dimensions),intent(in) :: multivals,mxvals
 131  type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc)
 132  type(results_out_type),intent(in) :: results_out(0:ndtset_alloc)
 133  character(len=8),intent(in) :: strimg(mxvals%nimage)
 134 
 135 !Local variables-------------------------------
 136 !scalars
 137  integer,parameter :: nkpt_max=50
 138  integer :: allowed,allowed_sum,iatom,idtset,ii,iimage,ikpt,kptopt,narr
 139  integer :: multival,multi_natfix,multi_natfixx,multi_natfixy,multi_natfixz
 140  integer :: multi_atsph,multi_occopt
 141  integer :: natfix,natfixx,natfixy,natfixz,natom
 142  integer :: ndtset_kptopt,nimage,nqpt,nkpt_eff
 143  integer :: ntypalch,ntypat,size1,size2,tnkpt
 144  real(dp) :: kpoint
 145  character(len=1) :: firstchar_gpu
 146 !arrays
 147  integer,allocatable :: iatfixio_(:,:),iatfixx_(:,:),iatfixy_(:,:)
 148  integer,allocatable :: iatfixz_(:,:),intarr(:,:),istwfk_2(:,:)
 149  integer,allocatable :: jdtset_kptopt(:),natfix_(:),natfixx_(:),natfixy_(:)
 150  integer,allocatable :: natfixz_(:)
 151  integer,allocatable :: narrm(:)
 152  integer,allocatable :: nimagem(:),prtimg(:,:)
 153  real(dp),allocatable :: dprarr(:,:),dprarr_images(:,:,:)
 154 
 155 ! *************************************************************************
 156 
 157 !###########################################################
 158 !### 01. Initial allocations and initialisations.
 159 
 160  DBG_ENTER("COLL")
 161 
 162  ABI_ALLOCATE(dprarr,(marr,0:ndtset_alloc))
 163  ABI_ALLOCATE(dprarr_images,(marr,mxvals%nimage,0:ndtset_alloc))
 164  ABI_ALLOCATE(intarr,(marr,0:ndtset_alloc))
 165  ABI_ALLOCATE(narrm,(0:ndtset_alloc))
 166  ABI_ALLOCATE(nimagem,(0:ndtset_alloc))
 167  ABI_ALLOCATE(prtimg,(mxvals%nimage,0:ndtset_alloc))
 168 
 169  do idtset=0,ndtset_alloc
 170    nimagem(idtset)=dtsets(idtset)%nimage
 171  end do
 172 
 173  firstchar_gpu=' '; if (maxval(dtsets(1:ndtset_alloc)%use_gpu_cuda)>0) firstchar_gpu='-'
 174 
 175 !if(multivals%natom==0)natom=dtsets(1)%natom
 176  natom=dtsets(1)%natom
 177 !if(multivals%nimage==0)nimage=dtsets(1)%nimage
 178  nimage=dtsets(1)%nimage
 179 !if(multivals%ntypalch==0)ntypalch=dtsets(1)%ntypalch
 180  ntypalch=dtsets(1)%ntypalch
 181 !if(multivals%ntypat==0)ntypat=dtsets(1)%ntypat
 182  ntypat=dtsets(1)%ntypat
 183 
 184 !###########################################################
 185 !### 02. Specific treatment for partially fixed atoms. Also compute multi_occopt for nband
 186 
 187 !Must treat separately the translation of iatfix from the internal
 188 !representation to the input/output representation
 189  ABI_ALLOCATE(natfix_,(0:ndtset_alloc))
 190  ABI_ALLOCATE(iatfixio_,(mxvals%natom,0:ndtset_alloc))
 191  ABI_ALLOCATE(natfixx_,(0:ndtset_alloc))
 192  ABI_ALLOCATE(iatfixx_,(mxvals%natom,0:ndtset_alloc))
 193  ABI_ALLOCATE(natfixy_,(0:ndtset_alloc))
 194  ABI_ALLOCATE(iatfixy_,(mxvals%natom,0:ndtset_alloc))
 195  ABI_ALLOCATE(natfixz_,(0:ndtset_alloc))
 196  ABI_ALLOCATE(iatfixz_,(mxvals%natom,0:ndtset_alloc))
 197  natfix_(0:ndtset_alloc)=0 ; iatfixio_(:,0:ndtset_alloc)=0
 198  natfixx_(0:ndtset_alloc)=0 ; iatfixx_(:,0:ndtset_alloc)=0
 199  natfixy_(0:ndtset_alloc)=0 ; iatfixy_(:,0:ndtset_alloc)=0
 200  natfixz_(0:ndtset_alloc)=0 ; iatfixz_(:,0:ndtset_alloc)=0
 201  natfixx=-1; natfixy=-1; natfixz=-1;nimage=-1;
 202  do idtset=1,ndtset_alloc
 203 !  DEBUG
 204 !  write(std_out,*)' outvar_i_n : iatfix_ for idtset= ',idtset
 205 !  ENDDEBUG
 206    do iatom=1,dtsets(idtset)%natom
 207 !    First look whether the atom is fixed along the three directions
 208      if( dtsets(idtset)%iatfix(1,iatom)+ &
 209 &     dtsets(idtset)%iatfix(2,iatom)+ &
 210 &     dtsets(idtset)%iatfix(3,iatom)   ==3 )then
 211        natfix_(idtset)=natfix_(idtset)+1
 212 !      DEBUG
 213 !      write(std_out,*)' outvar_i_n: iatom,natfix_(idtset)=',iatom,natfix_(idtset)
 214 !      ENDDEBUG
 215        iatfixio_(natfix_(idtset),idtset)=iatom
 216      else
 217 !      Now examine each direction, one at a time
 218        if( dtsets(idtset)%iatfix(1,iatom) ==1)then
 219          natfixx_(idtset)=natfixx_(idtset)+1
 220          iatfixx_(natfixx_(idtset),idtset)=iatom
 221        end if
 222        if( dtsets(idtset)%iatfix(2,iatom) ==1)then
 223          natfixy_(idtset)=natfixy_(idtset)+1
 224          iatfixy_(natfixy_(idtset),idtset)=iatom
 225        end if
 226        if( dtsets(idtset)%iatfix(3,iatom) ==1)then
 227          natfixz_(idtset)=natfixz_(idtset)+1
 228          iatfixz_(natfixz_(idtset),idtset)=iatom
 229        end if
 230      end if
 231    end do
 232 !  DEBUG
 233 !  write(std_out,*)' natfix ...'
 234 !  write(std_out,*)natfix_(idtset),natfixx_(idtset),natfixy_(idtset),natfixz_(idtset)
 235 !  ENDDEBUG
 236  end do
 237 
 238  multi_natfix=0
 239  if(ndtset_alloc>1)then
 240    do idtset=1,ndtset_alloc
 241      if(natfix_(1)/=natfix_(idtset))multi_natfix=1
 242    end do
 243  end if
 244  if(multi_natfix==0)natfix=natfix_(1)
 245 
 246  multi_natfixx=0
 247  if(ndtset_alloc>1)then
 248    do idtset=1,ndtset_alloc
 249      if(natfixx_(1)/=natfixx_(idtset))multi_natfixx=1
 250    end do
 251  end if
 252  if(multi_natfixx==0)natfixx=natfixx_(1)
 253 
 254  multi_natfixy=0
 255  if(ndtset_alloc>1)then
 256    do idtset=1,ndtset_alloc
 257      if(natfixy_(1)/=natfixy_(idtset))multi_natfixy=1
 258    end do
 259  end if
 260  if(multi_natfixy==0)natfixy=natfixy_(1)
 261 
 262  multi_natfixz=0
 263  if(ndtset_alloc>1)then
 264    do idtset=1,ndtset_alloc
 265      if(natfixz_(1)/=natfixz_(idtset))multi_natfixz=1
 266    end do
 267  end if
 268  if(multi_natfixz==0)natfixz=natfixz_(1)
 269 
 270 
 271  multi_occopt=0
 272  if(ndtset_alloc>1)then
 273    do idtset=1,ndtset_alloc
 274      if(dtsets(1)%occopt/=dtsets(idtset)%occopt)multi_occopt=1
 275    end do
 276  end if
 277 
 278 !write(ab_out,*)' outvar_i_n : I '
 279 !call flush(ab_out)
 280 !###########################################################
 281 !### 03. Print all the input variables (I)
 282 !##
 283 
 284 !iatfix
 285  narr=natfix                    ! default size for all datasets
 286  do idtset=0,ndtset_alloc       ! especific size for each dataset
 287    narrm(idtset)=natfix_(idtset)
 288    if(idtset==0)narrm(idtset)=mxvals%natom
 289    intarr(1:narrm(idtset),idtset)=iatfixio_(1:narrm(idtset),idtset)
 290  end do
 291  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 292 & narrm,ncid,ndtset_alloc,'iatfix','INT',multi_natfix)
 293 
 294 !iatfixx
 295  narr=natfixx                   ! default size for all datasets
 296  do idtset=0,ndtset_alloc       ! especific size for each dataset
 297    narrm(idtset)=natfixx_(idtset)
 298    if(idtset==0)narrm(idtset)=mxvals%natom
 299    intarr(1:narrm(idtset),idtset)=iatfixx_(1:narrm(idtset),idtset)
 300  end do
 301  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 302 & narrm,ncid,ndtset_alloc,'iatfixx','INT',multi_natfixx)
 303 
 304 !iatfixy
 305  narr=natfixy                   ! default size for all datasets
 306  do idtset=0,ndtset_alloc       ! especific size for each dataset
 307    narrm(idtset)=natfixy_(idtset)
 308    if(idtset==0)narrm(idtset)=mxvals%natom
 309    intarr(1:narrm(idtset),idtset)=iatfixy_(1:narrm(idtset),idtset)
 310  end do
 311  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 312 & narrm,ncid,ndtset_alloc,'iatfixy','INT',multi_natfixy)
 313 
 314 !iatfixz
 315  narr=natfixz                   ! default size for all datasets
 316  do idtset=0,ndtset_alloc       ! especific size for each dataset
 317    narrm(idtset)=natfixz_(idtset)
 318    if(idtset==0)narrm(idtset)=mxvals%natom
 319    intarr(1:narrm(idtset),idtset)=iatfixz_(1:narrm(idtset),idtset)
 320  end do
 321  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 322 & narrm,ncid,ndtset_alloc,'iatfixz','INT',multi_natfixz)
 323 
 324 !iatsph
 325  multi_atsph=1
 326  narr=dtsets(1)%natsph          ! default size for all datasets
 327  do idtset=0,ndtset_alloc       ! specific size for each dataset
 328    narrm(idtset)=dtsets(idtset)%natsph
 329    if(idtset==0)narrm(idtset)=mxvals%natsph
 330 !  Need to be printed only if there is some occurence of prtdos==3 or pawfatbnd
 331    if (narrm(idtset)>0) then
 332      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%iatsph(1:narrm(idtset))
 333    end if
 334    if(dtsets(idtset)%prtdos==3.or.dtsets(idtset)%pawfatbnd>0)then
 335      narrm(idtset)=dtsets(idtset)%natsph
 336    else
 337      narrm(idtset)=0
 338    end if
 339  end do
 340  if (ndtset_alloc==1.and.sum(narrm(1:ndtset_alloc))==1) multi_atsph=0
 341 
 342  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 343 & narrm,ncid,ndtset_alloc,'iatsph','INT',multi_atsph) ! Emulating the case of multiple narr
 344 
 345 
 346  intarr(1,:)=dtsets(:)%iboxcut
 347  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'iboxcut','INT',0)
 348 
 349  intarr(1,:)=dtsets(:)%icoulomb
 350  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'icoulomb','INT',0)
 351 
 352  intarr(1,:)=dtsets(:)%icutcoul
 353  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'icutcoul','INT',0)
 354 
 355  intarr(1,:)=dtsets(:)%ieig2rf
 356  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ieig2rf','INT',0)
 357 
 358  intarr(1,:)=dtsets(:)%imgmov
 359  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'imgmov','INT',0)
 360 
 361  intarr(1,:)=dtsets(:)%imgwfstor
 362  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'imgwfstor','INT',0)
 363 
 364  intarr(1,:)=dtsets(:)%inclvkb
 365  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'inclvkb','INT',0)
 366 
 367  intarr(1,:)=dtsets(:)%intxc
 368  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'intxc','INT',0)
 369 
 370  intarr(1,:)=dtsets(:)%ionmov
 371  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ionmov','INT',0)
 372 
 373  intarr(1,:)=dtsets(:)%iprcel
 374  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'iprcel','INT',0)
 375 
 376  intarr(1,:)=dtsets(:)%iprcfc
 377  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'iprcfc','INT',0)
 378 
 379  intarr(1,:)=dtsets(:)%irandom
 380  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irandom','INT',0)
 381 
 382  intarr(1,:)=dtsets(:)%irdbscoup
 383  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdbscoup','INT',0)
 384 
 385  intarr(1,:)=dtsets(:)%irdbseig
 386  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdbseig','INT',0)
 387 
 388  intarr(1,:)=dtsets(:)%irdbsreso
 389  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdbsreso','INT',0)
 390 
 391  intarr(1,:)=dtsets(:)%irdddk
 392  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdddk','INT',0)
 393 
 394  intarr(1,:)=dtsets(:)%irdden
 395  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdden','INT',0)
 396 
 397  intarr(1,:)=dtsets(:)%irdefmas
 398  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdefmas','INT',0)
 399 
 400  intarr(1,:)=dtsets(:)%irdhaydock
 401  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdhaydock','INT',0)
 402 
 403  intarr(1,:)=dtsets(:)%irdpawden
 404  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdpawden','INT',0)
 405 
 406  intarr(1,:)=dtsets(:)%irdqps
 407  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdqps','INT',0)
 408 
 409  intarr(1,:)=dtsets(:)%irdscr
 410  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdscr','INT',0)
 411 
 412  intarr(1,:)=dtsets(:)%irdsuscep
 413  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdsuscep','INT',0)
 414 
 415  intarr(1,:)=dtsets(:)%irdvdw
 416  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdvdw','INT',0)
 417 
 418  intarr(1,:)=dtsets(:)%irdwfk
 419  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdwfk','INT',0)
 420 
 421  intarr(1,:)=dtsets(:)%irdwfkfine
 422  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdwfkfine','INT',0)
 423 
 424  intarr(1,:)=dtsets(:)%irdwfq
 425  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'irdwfq','INT',0)
 426 
 427  intarr(1,:)=dtsets(:)%ird1den
 428  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ird1den','INT',0)
 429 
 430  intarr(1,:)=dtsets(:)%ird1wf
 431  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ird1wf','INT',0)
 432 
 433  intarr(1,:)=dtsets(:)%iscf
 434  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'iscf','INT',0)
 435 
 436  intarr(1,:)=dtsets(:)%isecur
 437  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'isecur','INT',0)
 438 
 439  intarr(1,:)=dtsets(:)%istatimg
 440  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'istatimg','INT',0)
 441 
 442  intarr(1,:)=dtsets(:)%istatr
 443  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'istatr','INT',0)
 444 
 445  intarr(1,:)=dtsets(:)%istatshft
 446  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'istatshft','INT',0)
 447 
 448 !istwfk (must first restore the default istwf=0 for non-allowed k points)
 449  ABI_ALLOCATE(istwfk_2,(mxvals%nkpt,0:ndtset_alloc))
 450  istwfk_2=0;allowed_sum=0
 451  do idtset=1,ndtset_alloc
 452    nqpt=dtsets(idtset)%nqpt
 453    do ikpt=1,dtsets(idtset)%nkpt
 454      allowed=1
 455      do ii=1,3
 456        kpoint=dtsets(idtset)%kpt(ii,ikpt)/dtsets(idtset)%kptnrm
 457        if(nqpt/=0.and.response_(idtset)==0)kpoint=kpoint+dtsets(idtset)%qptn(ii)
 458        if(abs(kpoint)>1.d-10.and.abs(kpoint-0.5_dp)>1.e-10_dp )allowed=0
 459      end do
 460      allowed_sum=allowed_sum+allowed
 461      if(allowed==1)istwfk_2(ikpt,idtset)=dtsets(idtset)%istwfk(ikpt)
 462    end do
 463  end do
 464 
 465 !istwfk
 466  tnkpt=0
 467  intarr(1:marr,0:ndtset_alloc)=0 ! default value
 468  do idtset=1,ndtset_alloc
 469    nkpt_eff=dtsets(idtset)%nkpt
 470    if(prtvol_glob==0 .and. nkpt_eff>nkpt_max)then
 471      nkpt_eff=nkpt_max
 472      tnkpt=1
 473    end if
 474    if((multivals%nkpt/=0).and.(sum(istwfk_2(1:nkpt_eff,idtset))==0)) nkpt_eff=0
 475    narrm(idtset)=nkpt_eff
 476    intarr(1:narrm(idtset),idtset)=istwfk_2(1:narrm(idtset),idtset)
 477  end do
 478 
 479  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,nkpt_eff,&
 480 & narrm,ncid,ndtset_alloc,'istwfk','INT',multivals%nkpt)
 481 
 482  if(tnkpt==1 .and. sum(istwfk_2(1:nkpt_eff,1:ndtset_alloc))/=0 ) &
 483 & write(iout,'(23x,a,i3,a)' ) 'outvar_i_n : Printing only first ',nkpt_max,' k-points.'
 484  ABI_DEALLOCATE(istwfk_2)
 485 
 486 !ixc
 487  intarr(1,:)=dtsets(:)%ixc
 488  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ixc','INT',0)
 489 
 490 !ixcpositron
 491  intarr(1,:)=dtsets(:)%ixcpositron
 492  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ixcpositron','INT',0)
 493 
 494 !ixcrot
 495  intarr(1,:)=dtsets(:)%ixcrot
 496  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ixcrot','INT',0)
 497 
 498 !ixc_sigma
 499  intarr(1,:)=dtsets(:)%ixc_sigma
 500  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ixc_sigma','INT',0)
 501 
 502 !write(ab_out,*)' outvar_i_n : J '
 503 !call flush(ab_out)
 504 !###########################################################
 505 !### 03. Print all the input variables (J)
 506 !##
 507 
 508  if (ndtset > 0) write(iout,"(1x,a16,1x,(t22,10i5))") 'jdtset',jdtset_(1:ndtset)
 509 
 510  intarr(1,:)=dtsets(:)%jellslab
 511  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'jellslab','INT',0)
 512 
 513  intarr(1,:)=dtsets(:)%jfielddir(1)
 514  intarr(2,:)=dtsets(:)%jfielddir(2)
 515  intarr(3,:)=dtsets(:)%jfielddir(3)
 516  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'jfielddir','INT',0)
 517 
 518 !jpawu
 519  prtimg(:,:)=1
 520  do idtset=0,ndtset_alloc
 521    narrm(idtset)=dtsets(idtset)%ntypat
 522    if (idtset==0) narrm(idtset)=mxvals%ntypat
 523    do iimage=1,nimagem(idtset)
 524      if (narrm(idtset)>0) then
 525        dprarr_images(1:narrm(idtset),iimage,idtset)=dtsets(idtset)%jpawu(1:narrm(idtset),iimage)
 526      end if
 527    end do
 528  end do
 529  call prttagm_images(dprarr_images,iout,jdtset_,1,marr,narrm,&
 530 & ncid,ndtset_alloc,'jpawu','ENE',mxvals%nimage,nimagem,ndtset,prtimg,strimg)
 531 
 532 
 533 !write(ab_out,*)' outvar_i_n : K '
 534 !call flush(ab_out)
 535 !###########################################################
 536 !### 03. Print all the input variables (K)
 537 !##
 538 
 539 !kberry
 540  narr=3*dtsets(1)%nberry ! default size for all datasets
 541  do idtset=0,ndtset_alloc       ! especific size for each dataset
 542    if(idtset/=0)then
 543      narrm(idtset)=3*dtsets(idtset)%nberry
 544      if (narrm(idtset)>0)&
 545 &     intarr(1:narrm(idtset),idtset)=&
 546 &     reshape(dtsets(idtset)%kberry(1:3,1:dtsets(idtset)%nberry), [narrm(idtset)] )
 547    else
 548      narrm(idtset)=3*mxvals%nberry
 549      if (narrm(idtset)>0)&
 550 &     intarr(1:narrm(idtset),idtset)=&
 551 &     reshape(dtsets(idtset)%kberry(1:3,1:mxvals%nberry), [narrm(idtset)] )
 552    end if
 553  end do
 554  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'kberry','INT',multivals%nberry)
 555 
 556 !kpt
 557  tnkpt=0
 558  dprarr(:,0)=0
 559  narr=3*dtsets(1)%nkpt            ! default size for all datasets
 560  if(prtvol_glob==0 .and. narr>3*nkpt_max)then
 561    narr=3*nkpt_max
 562    tnkpt=1
 563  end if
 564 
 565  do idtset=1,ndtset_alloc       ! specific size for each dataset
 566    narrm(idtset)=3*dtsets(idtset)%nkpt
 567    if (narrm(idtset)>0) then
 568      dprarr(1:narrm(idtset),idtset)=reshape(&
 569 &     dtsets(idtset)%kpt(1:3,1:dtsets(idtset)%nkpt), [narrm(idtset)] )
 570    end if
 571 
 572    if(prtvol_glob==0 .and. narrm(idtset)>3*nkpt_max)then
 573      narrm(idtset)=3*nkpt_max
 574      tnkpt=1
 575    end if
 576 
 577  end do
 578  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr, narrm,ncid,ndtset_alloc,'kpt','DPR',multivals%nkpt)
 579 
 580  if(tnkpt==1) write(iout,'(23x,a,i3,a)' ) 'outvar_i_n : Printing only first ',nkpt_max,' k-points.'
 581 
 582  !intarr(1,:)=dtsets(:)%kptbounds
 583  !call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'kptbounds','INT',0)
 584 
 585 !kptgw
 586  narr=3*dtsets(1)%nkptgw ! default size for all datasets
 587  do idtset=0,ndtset_alloc       ! specific size for each dataset
 588    if(idtset/=0)then
 589      narrm(idtset)=3*dtsets(idtset)%nkptgw
 590      if (narrm(idtset)>0)&
 591 &     dprarr(1:narrm(idtset),idtset)=&
 592 &     reshape(dtsets(idtset)%kptgw(1:3,1:dtsets(idtset)%nkptgw), [narrm(idtset)])
 593    else
 594      narrm(idtset)=mxvals%nkptgw
 595      if (narrm(idtset)>0)&
 596 &     dprarr(1:narrm(idtset),idtset)=&
 597 &     reshape(dtsets(idtset)%kptgw(1:3,1:mxvals%nkptgw), [narrm(idtset)] )
 598    end if
 599  end do
 600  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'kptgw','DPR',multivals%nkptgw)
 601 
 602 
 603 !kptns_hf
 604  if(sum(dtsets(1:ndtset_alloc)%usefock)/=0)then
 605    tnkpt=0
 606    dprarr(:,0)=0
 607    do idtset=1,ndtset_alloc       ! specific size for each dataset
 608      if(dtsets(idtset)%usefock/=0)then
 609        narrm(idtset)=3*dtsets(idtset)%nkpthf
 610        narr=narrm(idtset)
 611        if (narrm(idtset)>0) then
 612          dprarr(1:narrm(idtset),idtset)=reshape(&
 613 &         dtsets(idtset)%kptns_hf(1:3,1:dtsets(idtset)%nkpthf), [narrm(idtset)] )
 614        end if
 615      else
 616        narrm(idtset)=0
 617      end if
 618      if(prtvol_glob==0 .and. narrm(idtset)>3*nkpt_max)then
 619        narrm(idtset)=3*nkpt_max
 620        narr=narrm(idtset)
 621        tnkpt=1
 622      end if
 623    end do
 624    call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,narrm,ncid,ndtset_alloc,'kptns_hf','DPR',multivals%nkpthf)
 625    if(tnkpt==1) write(iout,'(23x,a,i3,a)' ) 'outvar_i_n : Printing only first ',nkpt_max,' k-points.'
 626  end if
 627 
 628  dprarr(1,:)=dtsets(:)%kptnrm
 629  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'kptnrm','DPR',0)
 630 
 631  intarr(1,:)=dtsets(:)%kptopt
 632  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'kptopt','INT',0)
 633 
 634 !kptrlatt
 635  if(sum((dtsets(1:ndtset_alloc)%kptopt)**2)/=0)then
 636    ndtset_kptopt=0
 637    intarr(1:9,0)=reshape( dtsets(0)%kptrlatt, [9] )
 638    ABI_ALLOCATE(jdtset_kptopt,(0:ndtset_alloc))
 639 !  Define the set of datasets for which kptopt>0
 640    do idtset=1,ndtset_alloc
 641      kptopt=dtsets(idtset)%kptopt
 642      if(kptopt>0)then
 643        ndtset_kptopt=ndtset_kptopt+1
 644        jdtset_kptopt(ndtset_kptopt)=jdtset_(idtset)
 645        intarr(1:9,ndtset_kptopt)=reshape( dtsets(idtset)%kptrlatt , [9] )
 646      end if
 647    end do
 648    if(ndtset_kptopt>0)then
 649      call prttagm(dprarr,intarr,iout,jdtset_kptopt,6,marr,9,narrm,ncid,ndtset_kptopt,'kptrlatt','INT',0)
 650    end if
 651    ABI_DEALLOCATE(jdtset_kptopt)
 652  end if
 653 
 654  dprarr(1,:)=dtsets(:)%kptrlen
 655  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'kptrlen','DPR',0)
 656 
 657  intarr(1,:)=dtsets(:)%kssform
 658  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'kssform','INT',0)
 659 
 660 !write(ab_out,*)' outvar_i_n : L '
 661 !call flush(ab_out)
 662 !###########################################################
 663 !### 03. Print all the input variables (L)
 664 !##
 665 
 666 !lexexch
 667  narr=mxvals%ntypat             ! default size for all datasets
 668  do idtset=0,ndtset_alloc       ! specific size for each dataset
 669    narrm(idtset)=dtsets(idtset)%ntypat
 670    if(idtset==0)narrm(idtset)=mxvals%ntypat
 671    if (narrm(idtset)>0) then
 672      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%lexexch(1:narrm(idtset))
 673    end if
 674  end do
 675  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,&
 676 & narrm,ncid,ndtset_alloc,'lexexch','INT',multivals%ntypat)
 677 
 678 !ldaminushalf
 679  do idtset=0,ndtset_alloc       ! specific size for each dataset
 680    narrm(idtset)=dtsets(idtset)%ntypat
 681    if(idtset==0)narrm(idtset)=mxvals%ntypat
 682    if (narrm(idtset)>0) then
 683      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%ldaminushalf(1:narrm(idtset))
 684    end if
 685  end do
 686  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,&
 687 & narrm,ncid,ndtset_alloc,'ldaminushalf','INT',multivals%ntypat)
 688 
 689 !localrdwf
 690  intarr(1,:)=dtsets(:)%localrdwf
 691  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'localrdwf','INT',0)
 692 
 693 !lpawu
 694  narr=mxvals%ntypat                    ! default size for all datasets
 695  do idtset=0,ndtset_alloc       ! specific size for each dataset
 696    narrm(idtset)=dtsets(idtset)%ntypat
 697    if(idtset==0)narrm(idtset)=mxvals%ntypat
 698    if (narrm(idtset)>0) then
 699      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%lpawu(1:narrm(idtset))
 700    end if
 701  end do
 702  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,&
 703 & narrm,ncid,ndtset_alloc,'lpawu','INT',multivals%ntypat)
 704 
 705 #if defined HAVE_LOTF
 706  if (any(dtsets(:)%ionmov==23)) then
 707    intarr(1,:)=dtsets(:)%lotf_classic
 708    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'lotf_classic','INT',0)
 709    intarr(1,:)=dtsets(:)%lotf_nitex
 710    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'lotf_nitex','INT',0)
 711    intarr(1,:)=dtsets(:)%lotf_nneigx
 712    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'lotf_nneigx','INT',0)
 713    intarr(1,:)=dtsets(:)%lotf_version
 714    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'lotf_version','INT',0)
 715  end if
 716 #endif
 717 
 718 
 719 !write(ab_out,*)' outvar_i_n : M '
 720 !call flush(ab_out)
 721 !###########################################################
 722 !### 03. Print all the input variables (M)
 723 !##
 724 
 725  intarr(1,:)=dtsets(:)%macro_uj
 726  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'macro_uj','INT',0)
 727 
 728  dprarr(1,:)=dtsets(:)%magcon_lambda
 729  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'magcon_lambda','DPR',0)
 730 
 731  intarr(1,:)=dtsets(:)%magconon
 732  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'magconon','INT',0)
 733 
 734  dprarr(1,:)=dtsets(:)%maxestep
 735  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'maxestep','ENE',0)
 736 
 737  intarr(1,:)=dtsets(:)%max_ncpus
 738  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'max_ncpus','INT',0)
 739 
 740  intarr(1,:)=dtsets(:)%maxnsym
 741  call prttagm(dprarr,intarr,iout,jdtset_,4,marr,1,narrm,ncid,ndtset_alloc,'maxnsym','INT',0)
 742 
 743  dprarr(1,:)=dtsets(:)%mbpt_sciss
 744  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'mbpt_sciss','ENE',0)
 745 
 746  dprarr(1,:)=dtsets(:)%mdf_epsinf
 747  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'mdf_epsinf','DPR',0)
 748 
 749  dprarr(1,:)=dtsets(:)%mdtemp(1)
 750  dprarr(2,:)=dtsets(:)%mdtemp(2)
 751  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'mdtemp','DPR',0)
 752 
 753  dprarr(1,:)=dtsets(:)%mdwall
 754  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'mdwall','LEN',0)
 755 
 756  intarr(1,:)=dtsets(:)%mem_test
 757  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'mem_test','INT',0)
 758 
 759  dprarr(1,:)=dtsets(:)%mep_mxstep
 760  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'mep_mxstep','LEN',0)
 761 
 762  intarr(1,:)=dtsets(:)%mep_solver
 763  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'mep_solver','INT',0)
 764 
 765  intarr(1,:)=dtsets(:)%mffmem
 766  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'mffmem','INT',0)
 767 
 768 !mixalch
 769  prtimg(:,:)=1
 770  do idtset=0,ndtset_alloc
 771    if(idtset/=0)then
 772      size1=dtsets(idtset)%npspalch ; size2=dtsets(idtset)%ntypalch
 773    else
 774      size1=npsp ; size2=mxvals%ntypat
 775    end if
 776    narrm(idtset)=size1*size2
 777    do iimage=1,nimagem(idtset)
 778      if (narrm(idtset)>0) then
 779        dprarr_images(1:narrm(idtset),iimage,idtset)=&
 780 &       reshape(results_out(idtset)%mixalch(1:size1,1:size2,iimage), [narrm(idtset)] )
 781      end if
 782    end do
 783  end do
 784  call prttagm_images(dprarr_images,iout,jdtset_,1,marr,narrm,ncid,ndtset_alloc,'mixalch','DPR',&
 785 & mxvals%nimage,nimagem,ndtset,prtimg,strimg)
 786 
 787 !mixesimgf
 788  dprarr(1:marr,0)=zero              ! default value
 789  narr=mxvals%nimage                 ! default size for all datasets
 790  if(any(abs(dtsets(1:ndtset_alloc)%imgmov-6)==0))then
 791    multival=multivals%nimage
 792    do idtset=1,ndtset_alloc           ! specific size and array for each dataset
 793      narrm(idtset)=dtsets(idtset)%nimage
 794      dprarr(1:narrm(idtset),idtset)=dtsets(idtset)%mixesimgf(1:narrm(idtset))
 795    end do
 796  else
 797    multival=0
 798    narrm(1:ndtset_alloc)=narr
 799    dprarr(1:marr,1:ndtset_alloc)=zero  
 800  endif
 801  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 802 & narrm,ncid,ndtset_alloc,'mixesimgf','DPR',multival)
 803 
 804  intarr(1,:)=dtsets(:)%mkmem
 805  call prttagm(dprarr,intarr,iout,jdtset_,5,marr,1,narrm,ncid,ndtset_alloc,'mkmem','INT',0)
 806 
 807  if(sum(response_(:))/=0)then
 808    intarr(1,:)=dtsets(:)%mkqmem
 809    call prttagm(dprarr,intarr,iout,jdtset_,5,marr,1,narrm,ncid,ndtset_alloc,'mkqmem','INT',0)
 810  end if
 811 
 812  if(sum(response_(:))/=0)then
 813    intarr(1,:)=dtsets(:)%mk1mem
 814    call prttagm(dprarr,intarr,iout,jdtset_,5,marr,1,narrm,ncid,ndtset_alloc,'mk1mem','INT',0)
 815  end if
 816 
 817  intarr(1,:)=dtsets(:)%mqgrid
 818  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'mqgrid','INT',0)
 819 
 820  intarr(1,:)=dtsets(:)%mqgriddg
 821  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'mqgriddg','INT',0)
 822 
 823 !###########################################################
 824 !### 03. Print all the input variables (N)
 825 !##
 826 
 827  intarr(1,:)=natfix_(:)
 828  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'natfix','INT',0)
 829 
 830  intarr(1,:)=natfixx_(:)
 831  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'natfixx','INT',0)
 832 
 833  intarr(1,:)=natfixy_(:)
 834  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'natfixy','INT',0)
 835 
 836  intarr(1,:)=natfixz_(:)
 837  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'natfixz','INT',0)
 838 
 839  intarr(1,:)=dtsets(0:ndtset_alloc)%natom
 840  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'natom','INT',0,forceprint=2)
 841 
 842 !natsph
 843 !Need to be printed only if there is some occurence of prtdos==3 or
 844 !pawfatbnd>0
 845  narr=1                      ! default size for all datasets
 846  do idtset=0,ndtset_alloc       ! especific size for each dataset
 847    narrm(idtset)=1
 848    intarr(1,idtset)=dtsets(idtset)%natsph
 849 
 850    if(dtsets(idtset)%prtdos==3.or.dtsets(idtset)%pawfatbnd>0)then
 851      narrm(idtset)=1
 852    else
 853      narrm(idtset)=0
 854    end if
 855  end do
 856  if (ndtset_alloc==1.and.sum(narrm(1:ndtset_alloc))==1) multi_atsph=0
 857  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,narr,&
 858 & narrm,ncid,ndtset_alloc,'natsph','INT',multi_atsph) ! Emulating multiple size for narrm
 859 
 860 !natsph_extra
 861  intarr(1,:)=dtsets(0:ndtset_alloc)%natsph_extra
 862  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'natsph_extra','INT',0)
 863 
 864  if(dtsets(1)%occopt==2)then
 865    narr=dtsets(1)%nkpt*dtsets(1)%nsppol                      ! default size for all datasets
 866  else
 867    narr=1
 868  end if
 869  do idtset=0,ndtset_alloc       ! specific size for each dataset
 870    if(dtsets(idtset)%occopt==2)then
 871      narrm(idtset)=dtsets(idtset)%nkpt*dtsets(idtset)%nsppol
 872    else
 873      narrm(idtset)=1
 874    end if
 875 
 876    if (narrm(idtset)>0) then
 877      intarr(1:narrm(idtset),idtset)=dtsets(idtset)%nband(1:narrm(idtset))
 878    end if
 879  end do
 880 
 881  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,&
 882 & narrm,ncid,ndtset_alloc,'nband','INT',multivals%nkpt+multivals%nsppol+multi_occopt)
 883 
 884  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%natvshift
 885  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'natvshift','INT',0)
 886 
 887  if(sum(dtsets(1:ndtset_alloc)%usefock)/=0)then
 888    intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%nbandhf
 889    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nbandhf','INT',0)
 890  end if
 891 
 892  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%nbandkss
 893  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nbandkss','INT',0)
 894 
 895  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%nbdblock
 896  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nbdblock','INT',0)
 897 
 898  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%nbdbuf
 899  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nbdbuf','INT',0)
 900 
 901  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%nberry
 902  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nberry','INT',0)
 903 
 904  intarr(1,:)=dtsets(:)%nc_xccc_gspace
 905  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nc_xccc_gspace','INT',0) !, firstchar="-")
 906 
 907  intarr(1,:)=dtsets(:)%nconeq
 908  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nconeq','INT',0)
 909 
 910  intarr(1,:)=dtsets(:)%nctime
 911  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nctime','INT',0)
 912 
 913 !ndtset
 914  if(ndtset>0)then
 915    intarr(1,:)=ndtset
 916    intarr(1,0)=0
 917    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ndtset','INT',0)
 918  end if
 919 
 920  intarr(1,:)=dtsets(:)%ndivsm
 921  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ndivsm','INT',0)
 922 
 923  !intarr(1,:)=dtsets(:)%nkpath
 924  !call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nkpath','INT',0)
 925 
 926  intarr(1,:)=dtsets(:)%ndynimage
 927  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ndynimage','INT',0)
 928 
 929  intarr(1,:)=dtsets(:)%neb_algo
 930  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'neb_algo','INT',0)
 931 
 932  dprarr(1,:)=dtsets(:)%neb_spring(1)
 933  dprarr(2,:)=dtsets(:)%neb_spring(2)
 934  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,2,narrm,ncid,ndtset_alloc,'neb_spring','DPR',0)
 935 
 936  intarr(1,:)=dtsets(:)%nfreqim
 937  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nfreqim','INT',0)
 938 
 939  intarr(1,:)=dtsets(:)%nfreqre
 940  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nfreqre','INT',0)
 941 
 942  intarr(1,:)=dtsets(:)%nfreqsp
 943  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nfreqsp','INT',0)
 944 
 945  intarr(1,:)=dtsets(:)%ngfft(1)
 946  intarr(2,:)=dtsets(:)%ngfft(2)
 947  intarr(3,:)=dtsets(:)%ngfft(3)
 948  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'ngfft','INT',0)
 949 
 950  intarr(1,:)=dtsets(:)%ngfftdg(1)
 951  intarr(2,:)=dtsets(:)%ngfftdg(2)
 952  intarr(3,:)=dtsets(:)%ngfftdg(3)
 953  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,3,narrm,ncid,ndtset_alloc,'ngfftdg','INT',0)
 954 
 955  intarr(1,:)=dtsets(:)%nimage
 956  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nimage','INT',0)
 957 
 958  intarr(1,:)=dtsets(:)%nkpt
 959  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nkpt','INT',0)
 960 
 961  intarr(1,:)=dtsets(:)%nkptgw
 962  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nkptgw','INT',0)
 963 
 964  if(sum(dtsets(1:ndtset_alloc)%usefock)/=0)then
 965    intarr(1,:)=dtsets(:)%nkpthf
 966 !  do idtset=1,ndtset_alloc       ! specific size for each dataset
 967 !    if(dtsets(idtset)%usefock/=0)then
 968 !      intarr(1,idtset)=dtsets(idtset)%nkpthf
 969 !    else
 970 !      intarr(1,idtset)=0
 971 !    endif
 972 !  enddo
 973    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nkpthf','INT',0)
 974  end if
 975 
 976  intarr(1,:)=dtsets(:)%nline
 977  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nline','INT',0)
 978 
 979  intarr(1,:)=dtsets(:)%nloalg(1)
 980  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nloc_alg','INT',0)
 981 
 982  intarr(1,:)=dtsets(:)%nloalg(2)*(dtsets(:)%nloalg(3)+1)
 983  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nloc_mem','INT',0)
 984 
 985  intarr(1,:)=dtsets(:)%nnos
 986  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nnos','INT',0)
 987 
 988  intarr(1,:)=dtsets(:)%nnsclo
 989  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nnsclo','INT',0)
 990 
 991  intarr(1,:)=dtsets(:)%nnsclohf
 992  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nnsclohf','INT',0)
 993 
 994  intarr(1,:)=dtsets(:)%nomegasf
 995  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nomegasf','INT',0)
 996 
 997  intarr(1,:)=dtsets(:)%nomegasi
 998  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nomegasi','INT',0)
 999 
1000  intarr(1,:)=dtsets(:)%nomegasrd
1001  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nomegasrd','INT',0)
1002 
1003  intarr(1,:)=dtsets(:)%nonlinear_info
1004  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nonlinear_info','INT',0)
1005 
1006  dprarr(1,:)=dtsets(:)%noseinert
1007  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'noseinert','DPR',0)
1008 
1009  intarr(1,:)=dtsets(:)%npband
1010  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npband','INT',0,firstchar="-")
1011 
1012  intarr(1,:)=dtsets(:)%npfft
1013  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npfft','INT',0,firstchar="-")
1014 
1015  intarr(1,:)=dtsets(:)%nphf
1016  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nphf','INT',0,firstchar="-")
1017 
1018  intarr(1,:)=dtsets(:)%npimage
1019  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npimage','INT',0,firstchar="-")
1020 
1021  intarr(1,:)=dtsets(:)%npkpt
1022  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npkpt','INT',0,firstchar='-')
1023 
1024  intarr(1,:)=dtsets(:)%nppert
1025  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nppert','INT',0,firstchar="-")
1026 
1027  if(multivals%ntypat/=0 .or. (multivals%ntypat==0 .and. ntypat/=npsp) )then
1028    intarr(1,:)=dtsets(:)%npsp
1029    call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npsp','INT',0)
1030  end if
1031 
1032  intarr(1,:)=dtsets(:)%npspinor
1033  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npspinor','INT',0,firstchar="-")
1034 
1035  intarr(1,:)=dtsets(0:ndtset_alloc)%npulayit
1036  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npulayit','INT',0)
1037 
1038  intarr(1,:)=dtsets(:)%npweps
1039  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npweps','INT',0)
1040 
1041  intarr(1,:)=dtsets(:)%npwkss
1042  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npwkss','INT',0)
1043 
1044  intarr(1,:)=dtsets(:)%npwsigx
1045  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npwsigx','INT',0)
1046 
1047  intarr(1,:)=dtsets(:)%npwwfn
1048  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npwwfn','INT',0)
1049 
1050  intarr(1,:)=dtsets(:)%np_slk
1051  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'np_slk','INT',0,firstchar="-")
1052 
1053  intarr(1,:)=dtsets(:)%nqpt
1054  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nqpt','INT',0)
1055 
1056  intarr(1,:)=dtsets(:)%nqptdm
1057  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'nqptdm','INT',0)
1058 
1059  intarr(1,:)=dtsets(:)%npvel
1060  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'npvel','INT',0)
1061 
1062  intarr(1,:)=dtsets(:)%nscforder
1063  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nscforder','INT',0)
1064 
1065 !nshiftk
1066  if(sum((dtsets(1:ndtset_alloc)%kptopt)**2)/=0)then
1067    ndtset_kptopt=0
1068    intarr(1:1,0)=dtsets(0)%nshiftk
1069    ABI_ALLOCATE(jdtset_kptopt,(0:ndtset_alloc))
1070 !  Define the set of datasets for which kptopt>0
1071    do idtset=1,ndtset_alloc
1072      kptopt=dtsets(idtset)%kptopt
1073      if(kptopt>0)then
1074        ndtset_kptopt=ndtset_kptopt+1
1075        jdtset_kptopt(ndtset_kptopt)=jdtset_(idtset)
1076        intarr(1:1,ndtset_kptopt)=dtsets(idtset)%nshiftk
1077      end if
1078    end do
1079    if(ndtset_kptopt>0)then
1080      call prttagm(dprarr,intarr,iout,jdtset_kptopt,2,marr,1,narrm,ncid,ndtset_kptopt,'nshiftk','INT',0)
1081    end if
1082    ABI_DEALLOCATE(jdtset_kptopt)
1083  end if
1084 
1085  intarr(1,:)=dtsets(:)%nspden
1086  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nspden','INT',0)
1087 
1088  intarr(1,:)=dtsets(:)%nspinor
1089  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nspinor','INT',0)
1090 
1091  intarr(1,:)=dtsets(:)%nsppol
1092  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nsppol','INT',0)
1093 
1094  intarr(1,:)=dtsets(:)%nstep
1095  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nstep','INT',0)
1096 
1097  intarr(1,:)=dtsets(:)%nsym
1098  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nsym','INT',0)
1099 
1100  intarr(1,:)=dtsets(:)%ntime
1101  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ntime','INT',0)
1102 
1103  intarr(1,:)=dtsets(:)%ntimimage
1104  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ntimimage','INT',0)
1105 
1106  intarr(1,:)=dtsets(:)%ntypalch
1107  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ntypalch','INT',0)
1108 
1109  intarr(1,:)=dtsets(:)%ntypat
1110  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'ntypat','INT',0,forceprint=2)
1111 
1112 !nucdipmom
1113  dprarr(:,0)=0.0_dp
1114  narr=3*natom ! default size for all datasets
1115  do idtset=1,ndtset_alloc       ! specific size for each dataset
1116    narrm(idtset)=3*dtsets(idtset)%natom
1117    if (narrm(idtset)>0) then
1118      dprarr(1:narrm(idtset),idtset)=&
1119 &     reshape(dtsets(idtset)%nucdipmom(1:3,1:dtsets(idtset)%natom), [narrm(idtset)])
1120    end if
1121    if(sum(abs( dtsets(idtset)%nucdipmom(1:3,1:dtsets(idtset)%natom))) < tol12 ) narrm(idtset)=0
1122  end do
1123  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,narr,&
1124 & narrm,ncid,ndtset_alloc,'nucdipmom','DPR',multivals%natom)
1125 
1126  intarr(1,:)=dtsets(:)%nwfshist
1127  call prttagm(dprarr,intarr,iout,jdtset_,2,marr,1,narrm,ncid,ndtset_alloc,'nwfshist','INT',0)
1128 
1129  intarr(1,0:ndtset_alloc)=dtsets(0:ndtset_alloc)%nzchempot
1130  call prttagm(dprarr,intarr,iout,jdtset_,1,marr,1,narrm,ncid,ndtset_alloc,'nzchempot','INT',0)
1131 
1132 !###########################################################
1133 !## Deallocation for generic arrays, and for i-n variables
1134 
1135  ABI_DEALLOCATE(dprarr)
1136  ABI_DEALLOCATE(intarr)
1137  ABI_DEALLOCATE(narrm)
1138  ABI_DEALLOCATE(nimagem)
1139  ABI_DEALLOCATE(dprarr_images)
1140  ABI_DEALLOCATE(prtimg)
1141 
1142  ABI_DEALLOCATE(natfix_)
1143  ABI_DEALLOCATE(iatfixio_)
1144  ABI_DEALLOCATE(natfixx_)
1145  ABI_DEALLOCATE(iatfixx_)
1146  ABI_DEALLOCATE(natfixy_)
1147  ABI_DEALLOCATE(iatfixy_)
1148  ABI_DEALLOCATE(natfixz_)
1149  ABI_DEALLOCATE(iatfixz_)
1150 
1151  DBG_EXIT("COLL")
1152 
1153 end subroutine outvar_i_n