TABLE OF CONTENTS


ABINIT/m_results_gs [ Modules ]

[ Top ] [ Modules ]

NAME

  m_results_gs

FUNCTION

  This module provides the definition of the results_gs_type
  used to store results from GS calculations.

COPYRIGHT

 Copyright (C) 2011-2018 ABINIT group (MT)
 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

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_results_gs
29 
30  use defs_basis
31  use m_abicore
32  use m_xmpi
33  use m_energies
34  use m_errors
35  use m_nctk
36 #ifdef HAVE_NETCDF
37  use netcdf
38 #endif
39 
40  use m_io_tools,  only : file_exists
41  use m_fstrings,  only : sjoin
42 
43  implicit none
44 
45  private

m_results_gs/copy_results_gs [ Functions ]

[ Top ] [ m_results_gs ] [ Functions ]

NAME

  copy_results_gs

FUNCTION

  Copy a results_gs datastructure into another

INPUTS

  results_gs_in=<type(results_gs_type)>=input results_gs datastructure

OUTPUT

  results_gs_out=<type(results_gs_type)>=output results_gs datastructure

PARENTS

      m_results_img

CHILDREN

      energies_copy,energies_ncwrite

SOURCE

587 subroutine copy_results_gs(results_gs_in,results_gs_out)
588 
589 
590 !This section has been created automatically by the script Abilint (TD).
591 !Do not modify the following lines by hand.
592 #undef ABI_FUNC
593 #define ABI_FUNC 'copy_results_gs'
594 !End of the abilint section
595 
596  implicit none
597 
598 !Arguments ------------------------------------
599 !arrays
600  type(results_gs_type),intent(in) :: results_gs_in
601  type(results_gs_type),intent(inout) :: results_gs_out !vz_i
602 
603 !Local variables-------------------------------
604 !scalars
605  integer :: natom_in,natom_out,ngrvdw_in,nsppol_in,nsppol_out
606 
607 !************************************************************************
608 
609  !@results_gs_type
610 
611  natom_in =results_gs_in%natom
612  natom_out=results_gs_out%natom
613  ngrvdw_in =results_gs_in%ngrvdw
614  nsppol_in =results_gs_in%nsppol
615  nsppol_out =results_gs_out%nsppol
616 
617  if (natom_in>natom_out) then
618    if (allocated(results_gs_out%fcart))   then
619      ABI_DEALLOCATE(results_gs_out%fcart)
620    end if
621    if (allocated(results_gs_out%fred))    then
622      ABI_DEALLOCATE(results_gs_out%fred)
623    end if
624    if (allocated(results_gs_out%gresid))  then
625      ABI_DEALLOCATE(results_gs_out%gresid)
626    end if
627    if (allocated(results_gs_out%grewtn))  then
628      ABI_DEALLOCATE(results_gs_out%grewtn)
629    end if
630   if (allocated(results_gs_out%grchempottn))  then
631      ABI_DEALLOCATE(results_gs_out%grchempottn)
632    end if
633    if (allocated(results_gs_out%grvdw))  then
634      ABI_DEALLOCATE(results_gs_out%grvdw)
635    end if
636    if (allocated(results_gs_out%grxc))    then
637      ABI_DEALLOCATE(results_gs_out%grxc)
638    end if
639    if (allocated(results_gs_out%synlgr))  then
640      ABI_DEALLOCATE(results_gs_out%synlgr)
641    end if
642 
643    if (allocated(results_gs_in%fcart))   then
644      ABI_ALLOCATE(results_gs_out%fcart,(3,natom_in))
645    end if
646    if (allocated(results_gs_in%fred))    then
647      ABI_ALLOCATE(results_gs_out%fred,(3,natom_in))
648    end if
649    if (allocated(results_gs_in%gresid))  then
650      ABI_ALLOCATE(results_gs_out%gresid,(3,natom_in))
651    end if
652    if (allocated(results_gs_in%grewtn))  then
653      ABI_ALLOCATE(results_gs_out%grewtn,(3,natom_in))
654    end if
655    if (allocated(results_gs_in%grchempottn))  then
656      ABI_ALLOCATE(results_gs_out%grchempottn,(3,natom_in))
657    end if
658    if (allocated(results_gs_in%grvdw))  then
659      ABI_ALLOCATE(results_gs_out%grvdw,(3,ngrvdw_in))
660    end if
661    if (allocated(results_gs_in%grxc))    then
662      ABI_ALLOCATE(results_gs_out%grxc,(3,natom_in))
663    end if
664    if (allocated(results_gs_in%synlgr))  then
665      ABI_ALLOCATE(results_gs_out%synlgr,(3,natom_in))
666    end if
667  end if
668 
669  if (nsppol_in>nsppol_out) then
670    if (allocated(results_gs_out%gaps))   then
671      ABI_DEALLOCATE(results_gs_out%gaps)
672    end if
673    if (allocated(results_gs_in%gaps))    then
674      ABI_ALLOCATE(results_gs_out%gaps,(3,nsppol_in))
675    end if
676  endif
677 
678 
679  results_gs_out%natom  =results_gs_in%natom
680  results_gs_out%ngrvdw =results_gs_in%ngrvdw
681  results_gs_out%nsppol =results_gs_in%nsppol
682  results_gs_out%deltae =results_gs_in%deltae
683  results_gs_out%diffor =results_gs_in%diffor
684  results_gs_out%entropy=results_gs_in%entropy
685  results_gs_out%etotal =results_gs_in%etotal
686  results_gs_out%fermie =results_gs_in%fermie
687  results_gs_out%residm =results_gs_in%residm
688  results_gs_out%res2   =results_gs_in%res2
689  results_gs_out%vxcavg =results_gs_in%vxcavg
690 
691  call energies_copy(results_gs_in%energies,results_gs_out%energies)
692 
693  results_gs_out%pel(:)=results_gs_in%pel(:)
694  results_gs_out%strten(:)=results_gs_in%strten(:)
695 
696  if (allocated(results_gs_in%fcart))  results_gs_out%fcart(:,1:natom_in) =results_gs_in%fcart(:,1:natom_in)
697  if (allocated(results_gs_in%fred))   results_gs_out%fred(:,1:natom_in)  =results_gs_in%fred(:,1:natom_in)
698  if (allocated(results_gs_in%gaps))   results_gs_out%gaps(:,1:nsppol_in) =results_gs_in%gaps(:,1:nsppol_in)
699  if (allocated(results_gs_in%gresid)) results_gs_out%gresid(:,1:natom_in)=results_gs_in%gresid(:,1:natom_in)
700  if (allocated(results_gs_in%grewtn)) results_gs_out%grewtn(:,1:natom_in)=results_gs_in%grewtn(:,1:natom_in)
701  if (allocated(results_gs_in%grchempottn))&
702 &  results_gs_out%grchempottn(:,1:natom_in)=results_gs_in%grchempottn(:,1:natom_in)
703  if (allocated(results_gs_in%grxc))   results_gs_out%grxc(:,1:natom_in)  =results_gs_in%grxc(:,1:natom_in)
704  if (allocated(results_gs_in%synlgr)) results_gs_out%synlgr(:,1:natom_in)=results_gs_in%synlgr(:,1:natom_in)
705  if (allocated(results_gs_in%grvdw).and.ngrvdw_in>0) then
706    results_gs_out%grvdw(:,1:ngrvdw_in)=results_gs_in%grvdw(:,1:ngrvdw_in)
707  end if
708 
709 end subroutine copy_results_gs

m_results_gs/destroy_results_gs [ Functions ]

[ Top ] [ m_results_gs ] [ Functions ]

NAME

  destroy_results_gs

FUNCTION

  Clean and destroy a results_gs datastructure

INPUTS

OUTPUT

SIDE EFFECTS

  results_gs(:)=<type(results_gs_type)>=results_gs datastructure

PARENTS

      m_results_img,mover_effpot

CHILDREN

      energies_copy,energies_ncwrite

SOURCE

422 subroutine destroy_results_gs(results_gs)
423 
424 
425 !This section has been created automatically by the script Abilint (TD).
426 !Do not modify the following lines by hand.
427 #undef ABI_FUNC
428 #define ABI_FUNC 'destroy_results_gs'
429 !End of the abilint section
430 
431  implicit none
432 
433 !Arguments ------------------------------------
434 !arrays
435  type(results_gs_type),intent(inout) :: results_gs
436 !Local variables-------------------------------
437 
438 !************************************************************************
439 
440  !@results_gs_type
441 
442  results_gs%natom =0
443  results_gs%ngrvdw=0
444  results_gs%nsppol=0
445  if (allocated(results_gs%fcart))   then
446    ABI_DEALLOCATE(results_gs%fcart)
447  end if
448  if (allocated(results_gs%fred))    then
449    ABI_DEALLOCATE(results_gs%fred)
450  end if
451  if (allocated(results_gs%gaps))    then
452    ABI_DEALLOCATE(results_gs%gaps)
453  end if
454  if (allocated(results_gs%gresid))  then
455    ABI_DEALLOCATE(results_gs%gresid)
456  end if
457  if (allocated(results_gs%grewtn))  then
458    ABI_DEALLOCATE(results_gs%grewtn)
459  end if
460  if (allocated(results_gs%grchempottn))  then
461    ABI_DEALLOCATE(results_gs%grchempottn)
462  end if
463  if (allocated(results_gs%grvdw))  then
464    ABI_DEALLOCATE(results_gs%grvdw)
465  end if
466  if (allocated(results_gs%grxc))    then
467    ABI_DEALLOCATE(results_gs%grxc)
468  end if
469  if (allocated(results_gs%synlgr))  then
470    ABI_DEALLOCATE(results_gs%synlgr)
471  end if
472 
473 end subroutine destroy_results_gs

m_results_gs/destroy_results_gs_array [ Functions ]

[ Top ] [ m_results_gs ] [ Functions ]

NAME

  destroy_results_gs_array

FUNCTION

  Clean and destroy a 2D-array of results_gs datastructures

INPUTS

OUTPUT

SIDE EFFECTS

  results_gs(:)=<type(results_gs_type)>=results_gs datastructure 2D-array

PARENTS

CHILDREN

      energies_copy,energies_ncwrite

SOURCE

499 subroutine destroy_results_gs_array(results_gs)
500 
501 
502 !This section has been created automatically by the script Abilint (TD).
503 !Do not modify the following lines by hand.
504 #undef ABI_FUNC
505 #define ABI_FUNC 'destroy_results_gs_array'
506 !End of the abilint section
507 
508  implicit none
509 
510 !Arguments ------------------------------------
511 !arrays
512  type(results_gs_type),intent(inout) :: results_gs(:,:)
513 !Local variables-------------------------------
514 !scalars
515  integer :: ii,jj,results_gs_size1,results_gs_size2
516 
517 !************************************************************************
518 
519  !@results_gs_type
520 
521  results_gs_size1=size(results_gs,1)
522  results_gs_size2=size(results_gs,2)
523  if (results_gs_size1>0.and.results_gs_size2>0) then
524 
525    do ii=1,results_gs_size2
526      do jj=1,results_gs_size1
527        results_gs(jj,ii)%natom =0
528        results_gs(jj,ii)%ngrvdw=0
529        results_gs(jj,ii)%nsppol=0
530        if (allocated(results_gs(jj,ii)%fcart))   then
531          ABI_DEALLOCATE(results_gs(jj,ii)%fcart)
532        end if
533        if (allocated(results_gs(jj,ii)%fred))    then
534          ABI_DEALLOCATE(results_gs(jj,ii)%fred)
535        end if
536        if (allocated(results_gs(jj,ii)%gaps))   then
537          ABI_DEALLOCATE(results_gs(jj,ii)%gaps)
538        end if
539        if (allocated(results_gs(jj,ii)%gresid))  then
540          ABI_DEALLOCATE(results_gs(jj,ii)%gresid)
541        end if
542        if (allocated(results_gs(jj,ii)%grewtn))  then
543          ABI_DEALLOCATE(results_gs(jj,ii)%grewtn)
544        end if
545        if (allocated(results_gs(jj,ii)%grchempottn))  then
546          ABI_DEALLOCATE(results_gs(jj,ii)%grchempottn)
547        end if
548        if (allocated(results_gs(jj,ii)%grvdw))  then
549          ABI_DEALLOCATE(results_gs(jj,ii)%grvdw)
550        end if
551        if (allocated(results_gs(jj,ii)%grxc))    then
552          ABI_DEALLOCATE(results_gs(jj,ii)%grxc)
553        end if
554        if (allocated(results_gs(jj,ii)%synlgr))  then
555          ABI_DEALLOCATE(results_gs(jj,ii)%synlgr)
556        end if
557      end do
558    end do
559 
560  end if
561 
562 end subroutine destroy_results_gs_array

m_results_gs/init_results_gs [ Functions ]

[ Top ] [ m_results_gs ] [ Functions ]

NAME

  init_results_gs

FUNCTION

  Init all (or part of) scalars and allocatables in a results_gs datastructure

INPUTS

  natom=number of atoms in cell
  nsppol=number of spin channels for this dataset
  only_part= --optional, default=false--
            if this flag is activated only the following parts of results_gs
            are initalized: all scalars, fcart,fred,strten

OUTPUT

SIDE EFFECTS

  results_gs=<type(results_gs_type)>=results_gs datastructure

PARENTS

      m_results_img,mover_effpot

CHILDREN

      energies_copy,energies_ncwrite

SOURCE

224 subroutine init_results_gs(natom,nsppol,results_gs,only_part)
225 
226 
227 !This section has been created automatically by the script Abilint (TD).
228 !Do not modify the following lines by hand.
229 #undef ABI_FUNC
230 #define ABI_FUNC 'init_results_gs'
231 !End of the abilint section
232 
233  implicit none
234 
235 !Arguments ------------------------------------
236 !scalars
237  integer,intent(in) :: natom,nsppol
238  logical,optional,intent(in) :: only_part
239 !arrays
240  type(results_gs_type),intent(inout) :: results_gs
241 !Local variables-------------------------------
242 !scalars
243  logical :: full_init
244 !arrays
245 
246 !************************************************************************
247 
248  !@results_gs_type
249 
250  full_init=.true.;if (present(only_part)) full_init=(.not.only_part)
251 
252  results_gs%natom  =natom
253  results_gs%ngrvdw =0
254  results_gs%nsppol =nsppol
255  results_gs%deltae =zero
256  results_gs%diffor =zero
257  results_gs%entropy=zero
258  results_gs%etotal =zero
259  results_gs%fermie =zero
260  results_gs%residm =zero
261  results_gs%res2   =zero
262  results_gs%vxcavg =zero
263 
264  call energies_init(results_gs%energies)
265 
266  results_gs%strten=zero
267  ABI_ALLOCATE(results_gs%fcart,(3,natom))
268  results_gs%fcart=zero
269  ABI_ALLOCATE(results_gs%fred,(3,natom))
270  results_gs%fred =zero
271  ABI_ALLOCATE(results_gs%gaps,(3,nsppol))
272  results_gs%gaps =zero
273  if (full_init) then
274    results_gs%pel=zero
275    ABI_ALLOCATE(results_gs%gresid,(3,natom))
276    results_gs%gresid=zero
277    ABI_ALLOCATE(results_gs%grewtn,(3,natom))
278    results_gs%grewtn=zero
279    ABI_ALLOCATE(results_gs%grchempottn,(3,natom))
280    results_gs%grchempottn=zero
281    ABI_ALLOCATE(results_gs%grxc,(3,natom))
282    results_gs%grxc  =zero
283    ABI_ALLOCATE(results_gs%synlgr,(3,natom))
284    results_gs%synlgr=zero
285    ABI_ALLOCATE(results_gs%grvdw,(3,results_gs%ngrvdw))
286  end if
287 
288 end subroutine init_results_gs

m_results_gs/init_results_gs_array [ Functions ]

[ Top ] [ m_results_gs ] [ Functions ]

NAME

  init_results_gs_array

FUNCTION

  Init all (or part of) scalars and allocatables in a 2D-array of results_gs datastructures

INPUTS

  natom=number of atoms in cell
  nsppol=number of spin channels for this dataset
  only_part= --optional, default=false--
            if this flag is activated only the following parts of results_gs
            are initalized: all scalars, fcart,fred,strten

OUTPUT

SIDE EFFECTS

  results_gs(:)=<type(results_gs_type)>=results_gs datastructure 2Darray

PARENTS

CHILDREN

      energies_copy,energies_ncwrite

SOURCE

319 subroutine init_results_gs_array(natom,nsppol,results_gs,only_part)
320 
321 
322 !This section has been created automatically by the script Abilint (TD).
323 !Do not modify the following lines by hand.
324 #undef ABI_FUNC
325 #define ABI_FUNC 'init_results_gs_array'
326 !End of the abilint section
327 
328  implicit none
329 
330 !Arguments ------------------------------------
331 !scalars
332  integer,intent(in) :: natom,nsppol
333  logical,optional,intent(in) :: only_part
334 !arrays
335  type(results_gs_type),intent(inout) :: results_gs(:,:)
336 !Local variables-------------------------------
337 !scalars
338  integer :: ii,jj,results_gs_size1,results_gs_size2
339  logical :: full_init
340 !arrays
341 
342 !************************************************************************
343 
344  !@results_gs_type
345 
346  results_gs_size1=size(results_gs,1)
347  results_gs_size2=size(results_gs,2)
348  full_init=.true.;if (present(only_part)) full_init=(.not.only_part)
349 
350  if (results_gs_size1>0.and.results_gs_size2>0) then
351 
352    do ii=1,results_gs_size2
353      do jj=1,results_gs_size1
354 
355        results_gs(jj,ii)%natom  =natom
356        results_gs(jj,ii)%ngrvdw =0
357        results_gs(jj,ii)%nsppol =nsppol
358        results_gs(jj,ii)%deltae =zero
359        results_gs(jj,ii)%diffor =zero
360        results_gs(jj,ii)%entropy=zero
361        results_gs(jj,ii)%etotal =zero
362        results_gs(jj,ii)%fermie =zero
363        results_gs(jj,ii)%residm =zero
364        results_gs(jj,ii)%res2   =zero
365        results_gs(jj,ii)%vxcavg =zero
366 
367        call energies_init(results_gs(jj,ii)%energies)
368 
369        results_gs(jj,ii)%strten=zero
370        ABI_ALLOCATE(results_gs(jj,ii)%fcart,(3,natom))
371        results_gs(jj,ii)%fcart=zero
372        ABI_ALLOCATE(results_gs(jj,ii)%fred,(3,natom))
373        results_gs(jj,ii)%fred =zero
374        ABI_ALLOCATE(results_gs(jj,ii)%gaps,(3,nsppol))
375        results_gs(jj,ii)%gaps =zero
376        if (full_init) then
377          results_gs(jj,ii)%pel=zero
378          ABI_ALLOCATE(results_gs(jj,ii)%gresid,(3,natom))
379          results_gs(jj,ii)%gresid=zero
380          ABI_ALLOCATE(results_gs(jj,ii)%grewtn,(3,natom))
381          results_gs(jj,ii)%grewtn=zero
382          ABI_ALLOCATE(results_gs(jj,ii)%grchempottn,(3,natom))
383          results_gs(jj,ii)%grchempottn=zero
384          ABI_ALLOCATE(results_gs(jj,ii)%grxc,(3,natom))
385          results_gs(jj,ii)%grxc  =zero
386          ABI_ALLOCATE(results_gs(jj,ii)%synlgr,(3,natom))
387          results_gs(jj,ii)%synlgr=zero
388          ABI_ALLOCATE(results_gs(jj,ii)%grvdw,(3,results_gs(jj,ii)%ngrvdw))
389        end if
390 
391      end do
392    end do
393  end if
394 
395 end subroutine init_results_gs_array

m_results_gs/results_gs_ncwrite [ Functions ]

[ Top ] [ m_results_gs ] [ Functions ]

NAME

 results_gs_ncwrite

FUNCTION

INPUTS

  ncid=NC file handle
  ecut, pawecutdg= Input cutoff energies in Ha.

OUTPUT

PARENTS

      m_results_gs

CHILDREN

      energies_ncwrite,results_gs_ncwrite

SOURCE

734 integer function results_gs_ncwrite(res,ncid,ecut,pawecutdg) result(ncerr)
735 
736 
737 !This section has been created automatically by the script Abilint (TD).
738 !Do not modify the following lines by hand.
739 #undef ABI_FUNC
740 #define ABI_FUNC 'results_gs_ncwrite'
741 !End of the abilint section
742 
743  implicit none
744 
745 !Arguments ------------------------------------
746 !scalars
747  integer,intent(in) :: ncid
748  real(dp),intent(in) :: ecut,pawecutdg
749  type(results_gs_type),intent(in) :: res
750 
751 !Local variables-------------------------------
752 !scalars
753 #ifdef HAVE_NETCDF
754 
755 ! *************************************************************************
756 
757  ! ==============================================
758  ! === Write the dimensions specified by ETSF ===
759  ! ==============================================
760  !FIXME: do not handle k_dependent = 1
761  ncerr = nctk_def_dims(ncid, [nctkdim_t("six", 6), nctkdim_t("number_of_cartesian_dimensions", 3), &
762    nctkdim_t("number_of_atoms", res%natom), nctkdim_t("number_of_spins", res%nsppol)], defmode=.True.)
763  NCF_CHECK(ncerr)
764 
765 ! Define variables.
766 ! scalars passed in input (not belonging to results_gs) as well as scalars defined in results_gs
767  ncerr = nctk_def_dpscalars(ncid, [character(len=nctk_slen) :: &
768    "ecut", "pawecutdg", "deltae", "diffor", "entropy", "etotal", "fermie", "residm", "res2"])
769  NCF_CHECK(ncerr)
770 
771  ! arrays
772  !
773  ! Note: unlike fred, this array has been corrected by enforcing
774  ! the translational symmetry, namely that the sum of force on all atoms is zero.
775 
776  ncerr = nctk_def_arrays(ncid, [&
777    nctkarr_t('cartesian_forces', "dp", "number_of_cartesian_dimensions, number_of_atoms"),&
778    nctkarr_t('cartesian_stress_tensor', "dp", 'six')])
779  NCF_CHECK(ncerr)
780 
781 ! Write data.
782 ! Write variables
783  ncerr = nctk_write_dpscalars(ncid, [character(len=nctk_slen) :: &
784 &  'ecut', 'pawecutdg', 'deltae', 'diffor', 'entropy', 'etotal', 'fermie', 'residm', 'res2'],&
785 &  [ecut, pawecutdg, res%deltae, res%diffor, res%entropy, res%etotal, res%fermie, res%residm, res%res2],&
786 &  datamode=.True.)
787  NCF_CHECK(ncerr)
788 
789  NCF_CHECK(nctk_set_datamode(ncid))
790  NCF_CHECK(nf90_put_var(ncid, vid("cartesian_forces"), res%fcart))
791  NCF_CHECK(nf90_put_var(ncid, vid("cartesian_stress_tensor"), res%strten))
792 
793 ! Add energies
794  call energies_ncwrite(res%energies, ncid)
795 
796 #else
797  MSG_ERROR("netcdf support is not activated.")
798 #endif
799 
800 contains
801  integer function vid(vname)
802 
803 
804 !This section has been created automatically by the script Abilint (TD).
805 !Do not modify the following lines by hand.
806 #undef ABI_FUNC
807 #define ABI_FUNC 'vid'
808 !End of the abilint section
809 
810    character(len=*),intent(in) :: vname
811    vid = nctk_idname(ncid, vname)
812  end function vid
813 
814 end function results_gs_ncwrite

m_results_gs/results_gs_type [ Types ]

[ Top ] [ m_results_gs ] [ Types ]

NAME

 results_gs_type

FUNCTION

 This structured datatype contains the results of a GS calculation :
 energy and its decomposition, forces and their decompositions, stresses
 and their decompositions

SOURCE

 59  type, public :: results_gs_type
 60 
 61 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
 62 ! declared in another part of ABINIT, that might need to take into account your modification.
 63 
 64 ! Integer scalar
 65 
 66   integer :: natom
 67    ! The number of atoms for this dataset
 68 
 69   integer :: nsppol
 70    ! The number of spin channels for this dataset
 71 
 72   integer :: ngrvdw
 73    ! Size of grvdw array
 74    ! Can be 0 (not allocated) or natom
 75 
 76 ! Real (real(dp)) scalars
 77 
 78   real(dp) :: deltae
 79    ! change in energy (Hartree)
 80 
 81   real(dp) :: diffor
 82    ! maximal absolute value of changes in the components of force
 83 
 84 ! All the energies are in Hartree, obtained "per unit cell".
 85   type(energies_type) :: energies
 86 !!!  real(dp) :: eei      ! local pseudopotential energy (Hartree)
 87 !!!  real(dp) :: eeig     ! sum of eigenvalue energy (Hartree)
 88 !!!  real(dp) :: eew      ! Ewald energy (Hartree)
 89 !!!  real(dp) :: ehart    ! Hartree part of total energy (Hartree)
 90 !!!  real(dp) :: eii      ! pseudopotential core-core energy
 91 !!!  real(dp) :: ek       ! kinetic energy (Hartree)
 92 !!!  real(dp) :: enefield ! the term of the energy functional that depends
 93 !!!                       ! explicitely on the electric field
 94 !!!                       ! enefield = -ucvol*E*P
 95 !!!  real(dp) :: enl      ! nonlocal pseudopotential energy (Hartree)
 96   real(dp) :: entropy  ! entropy (Hartree)
 97 !!!  real(dp) :: enxc     ! exchange-correlation energy (Hartree)
 98 !!!  real(dp) :: enxcdc   ! exchange-correlation double-counting energy (Hartree)
 99 !!!  real(dp) :: epaw     ! PAW spherical energy (Hartree)
100 !!!  real(dp) :: epawdc   ! PAW spherical double-counting energy (Hartree)
101   real(dp) :: etotal   ! total energy (Hartree)
102                        ! for fixed occupation numbers (occopt==0,1,or 2):
103                        !   etotal=ek+ehart+enxc+eei+eew+eii+enl+PAW_spherical_part
104                        ! for varying occupation numbers (occopt>=3):
105                        !   etotal=ek+ehart+enxc+eei+eew+eii+enl - tsmear*entropy +PAW_spherical_part
106   real(dp) :: fermie   ! Fermi energy (Hartree)
107   real(dp) :: residm   ! maximum value for the residual over all bands, all k points,
108                        !   and all spins (Hartree or Hartree**2, to be checked !)
109   real(dp) :: res2     ! density/potential residual (squared)
110   real(dp) :: vxcavg   ! Average of the exchange-correlation energy. The average
111                        ! of the local psp pot and the Hartree pot is set to zero (due
112                        ! to the usual problem at G=0 for Coulombic system, so vxcavg
113                        ! is also the average of the local part of the Hamiltonian
114 
115 ! Real (real(dp)) arrays
116 
117   real(dp), allocatable :: fcart(:,:)
118    ! fcart(3,natom)
119    ! Cartesian forces (Hartree/Bohr)
120    ! Note : unlike fred, this array has been corrected by enforcing
121    ! the translational symmetry, namely that the sum of force
122    ! on all atoms is zero.
123 
124   real(dp), allocatable :: fred(:,:)
125    ! fred(3,natom)
126    ! Forces in reduced coordinates (Hartree)
127    ! Actually, gradient of the total energy with respect
128    ! to change of reduced coordinates
129 
130   real(dp),allocatable :: gaps(:,:)
131    ! gaps(3,nsppol)
132    ! gaps(1,:) : fundamental gap
133    ! gaps(2,:) : optical gap
134    ! gaps(3,:) : "status" for each channel : 0.0dp if the gap was not computed
135    !   (because there are only valence bands) ; -1.0dp if the system (or spin-channel) is metallic ; 1.0dp if the
136    !   gap was computed
137 
138   real(dp), allocatable :: gresid(:,:)
139    ! gresid(3,natom)
140    ! Part of the gradient of the total energy (Hartree) with respect
141    ! to change of reduced coordinates, that comes from the residual
142    ! of the potential
143 
144   real(dp), allocatable :: grewtn(:,:)
145    ! grewtn(3,natom)
146    ! Part of the gradient of the total energy (Hartree) with respect
147    ! to change of reduced coordinates, that comes from the Ewald energy
148 
149   real(dp), allocatable :: grchempottn(:,:)
150    ! grchempottn(3,natom)
151    ! Part of the gradient of the total energy (Hartree) with respect
152    ! to change of reduced coordinates, that comes from the spatially-varying chemical potential
153 
154   real(dp), allocatable :: grvdw(:,:)
155    ! grvdw(3,ngrvdw)
156    ! Part of the gradient of the total energy (Hartree) with respect
157    ! to change of reduced coordinates, that comes from
158    ! Van der Waals DFT-D2 dispersion (hartree)
159    ! ngrvdw can be 0 or natom
160 
161   real(dp), allocatable :: grxc(:,:)
162    ! grxc(3,natom)
163    ! Part of the gradient of the total energy (Hartree) with respect
164    ! to change of reduced coordinates, that comes from the XC energy
165 
166   real(dp) :: pel(3)
167    ! ucvol times the electronic polarization in reduced coordinates
168 
169   real(dp) :: strten(6)
170    ! Stress tensor in cartesian coordinates (Hartree/Bohr^3)
171    ! 6 unique components of this symmetric 3x3 tensor:
172    ! Given in order (1,1), (2,2), (3,3), (3,2), (3,1), (2,1).
173 
174   real(dp), allocatable :: synlgr(:,:)
175    ! synlgr(3,natom)
176    ! Part of the gradient of the total energy (Hartree) with respect
177    ! to change of reduced coordinates, that comes from the non-local energy
178    ! The "sy" prefix refer to the fact that this gradient has been
179    ! symmetrized.
180 
181  end type results_gs_type
182 
183 !public procedures.
184  public :: init_results_gs
185  public :: init_results_gs_array
186  public :: destroy_results_gs
187  public :: destroy_results_gs_array
188  public :: copy_results_gs
189  public :: results_gs_ncwrite