TABLE OF CONTENTS
- ABINIT/m_results_gs
- m_results_gs/copy_results_gs
- m_results_gs/destroy_results_gs
- m_results_gs/destroy_results_gs_array
- m_results_gs/init_results_gs
- m_results_gs/init_results_gs_array
- m_results_gs/results_gs_ncwrite
- m_results_gs/results_gs_type
- m_results_gs/results_gs_yaml_write
ABINIT/m_results_gs [ 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-2022 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 .
SOURCE
18 #if defined HAVE_CONFIG_H 19 #include "config.h" 20 #endif 21 22 #include "abi_common.h" 23 24 MODULE m_results_gs 25 26 use defs_basis 27 use m_abicore 28 use m_xmpi 29 use m_energies 30 use m_errors 31 use m_yaml 32 use m_crystal 33 use m_stream_string 34 use m_pair_list 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 use m_numeric_tools, only : get_trace 43 use m_geometry, only : stress_voigt_to_mat 44 45 implicit none 46 47 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
SOURCE
559 subroutine copy_results_gs(results_gs_in,results_gs_out) 560 561 !Arguments ------------------------------------ 562 !arrays 563 class(results_gs_type),intent(in) :: results_gs_in 564 type(results_gs_type),intent(inout) :: results_gs_out !vz_i 565 566 !Local variables------------------------------- 567 !scalars 568 integer :: natom_in,natom_out,ngrvdw_in,nspden_in,nspden_out,nsppol_in,nsppol_out 569 570 !************************************************************************ 571 572 !@results_gs_type 573 574 natom_in =results_gs_in%natom 575 natom_out =results_gs_out%natom 576 ngrvdw_in =results_gs_in%ngrvdw 577 nspden_in =results_gs_in%nspden 578 nspden_out=results_gs_out%nspden 579 nsppol_in =results_gs_in%nsppol 580 nsppol_out=results_gs_out%nsppol 581 582 if (natom_in>natom_out) then 583 ABI_SFREE(results_gs_out%fcart) 584 ABI_SFREE(results_gs_out%gred) 585 ABI_SFREE(results_gs_out%grchempottn) 586 ABI_SFREE(results_gs_out%grcondft) 587 ABI_SFREE(results_gs_out%gresid) 588 ABI_SFREE(results_gs_out%grewtn) 589 ABI_SFREE(results_gs_out%grvdw) 590 ABI_SFREE(results_gs_out%grxc) 591 ABI_SFREE(results_gs_out%intgres) 592 ABI_SFREE(results_gs_out%synlgr) 593 594 if (allocated(results_gs_in%fcart)) then 595 ABI_MALLOC(results_gs_out%fcart,(3,natom_in)) 596 end if 597 if (allocated(results_gs_in%gred)) then 598 ABI_MALLOC(results_gs_out%gred,(3,natom_in)) 599 end if 600 if (allocated(results_gs_in%gresid)) then 601 ABI_MALLOC(results_gs_out%gresid,(3,natom_in)) 602 end if 603 if (allocated(results_gs_in%grchempottn)) then 604 ABI_MALLOC(results_gs_out%grchempottn,(3,natom_in)) 605 end if 606 if (allocated(results_gs_in%grcondft)) then 607 ABI_MALLOC(results_gs_out%grcondft,(3,natom_in)) 608 end if 609 if (allocated(results_gs_in%grewtn)) then 610 ABI_MALLOC(results_gs_out%grewtn,(3,natom_in)) 611 end if 612 if (allocated(results_gs_in%grvdw)) then 613 ABI_MALLOC(results_gs_out%grvdw,(3,ngrvdw_in)) 614 end if 615 if (allocated(results_gs_in%grxc)) then 616 ABI_MALLOC(results_gs_out%grxc,(3,natom_in)) 617 end if 618 if (allocated(results_gs_in%synlgr)) then 619 ABI_MALLOC(results_gs_out%synlgr,(3,natom_in)) 620 end if 621 end if 622 623 if (nsppol_in>nsppol_out) then 624 ABI_SFREE(results_gs_out%gaps) 625 if (allocated(results_gs_in%gaps)) then 626 ABI_MALLOC(results_gs_out%gaps,(3,nsppol_in)) 627 end if 628 endif 629 630 if (nspden_in>nspden_out .or. natom_in>natom_out) then 631 ABI_SFREE(results_gs_out%intgres) 632 if (allocated(results_gs_in%intgres)) then 633 ABI_MALLOC(results_gs_out%intgres,(max(nspden_in,nspden_out),max(natom_in,natom_out))) 634 end if 635 endif 636 637 638 results_gs_out%natom =results_gs_in%natom 639 results_gs_out%ngrvdw =results_gs_in%ngrvdw 640 results_gs_out%nspden =results_gs_in%nspden 641 results_gs_out%nsppol =results_gs_in%nsppol 642 results_gs_out%berryopt=results_gs_in%berryopt 643 results_gs_out%deltae =results_gs_in%deltae 644 results_gs_out%diffor =results_gs_in%diffor 645 results_gs_out%entropy=results_gs_in%entropy 646 results_gs_out%entropy_extfpmd=results_gs_in%entropy_extfpmd 647 results_gs_out%etotal =results_gs_in%etotal 648 results_gs_out%fermie =results_gs_in%fermie 649 results_gs_out%fermih =results_gs_in%fermih ! CP added for occopt 9 650 results_gs_out%nelect_extfpmd=results_gs_in%nelect_extfpmd 651 results_gs_out%residm =results_gs_in%residm 652 results_gs_out%res2 =results_gs_in%res2 653 results_gs_out%shiftfactor_extfpmd=results_gs_in%shiftfactor_extfpmd 654 results_gs_out%vxcavg =results_gs_in%vxcavg 655 656 call energies_copy(results_gs_in%energies,results_gs_out%energies) 657 658 results_gs_out%pel(:)=results_gs_in%pel(:) 659 results_gs_out%pion(:)=results_gs_in%pion(:) 660 results_gs_out%strten(:)=results_gs_in%strten(:) 661 662 if (allocated(results_gs_in%fcart)) results_gs_out%fcart(:,1:natom_in) =results_gs_in%fcart(:,1:natom_in) 663 if (allocated(results_gs_in%gred)) results_gs_out%gred(:,1:natom_in) =results_gs_in%gred(:,1:natom_in) 664 if (allocated(results_gs_in%gaps)) results_gs_out%gaps(:,1:nsppol_in) =results_gs_in%gaps(:,1:nsppol_in) 665 if (allocated(results_gs_in%grchempottn))& 666 & results_gs_out%grchempottn(:,1:natom_in)=results_gs_in%grchempottn(:,1:natom_in) 667 if (allocated(results_gs_in%grcondft)) results_gs_out%grcondft(:,1:natom_in)=results_gs_in%grcondft(:,1:natom_in) 668 if (allocated(results_gs_in%gresid)) results_gs_out%gresid(:,1:natom_in)=results_gs_in%gresid(:,1:natom_in) 669 if (allocated(results_gs_in%grewtn)) results_gs_out%grewtn(:,1:natom_in)=results_gs_in%grewtn(:,1:natom_in) 670 if (allocated(results_gs_in%grxc)) results_gs_out%grxc(:,1:natom_in) =results_gs_in%grxc(:,1:natom_in) 671 if (allocated(results_gs_in%intgres))results_gs_out%intgres(1:nspden_in,1:natom_in) =results_gs_in%intgres(1:nspden_in,1:natom_in) 672 if (allocated(results_gs_in%synlgr)) results_gs_out%synlgr(:,1:natom_in)=results_gs_in%synlgr(:,1:natom_in) 673 if (allocated(results_gs_in%grvdw).and.ngrvdw_in>0) then 674 results_gs_out%grvdw(:,1:ngrvdw_in)=results_gs_in%grvdw(:,1:ngrvdw_in) 675 end if 676 677 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
SOURCE
450 subroutine destroy_results_gs(results_gs) 451 452 !Arguments ------------------------------------ 453 !arrays 454 type(results_gs_type),intent(inout) :: results_gs 455 456 !************************************************************************ 457 458 !@results_gs_type 459 460 results_gs%natom =0 461 results_gs%ngrvdw=0 462 results_gs%nspden=0 463 results_gs%nsppol=0 464 results_gs%berryopt=0 465 466 ABI_SFREE(results_gs%fcart) 467 ABI_SFREE(results_gs%gred) 468 ABI_SFREE(results_gs%gaps) 469 ABI_SFREE(results_gs%grcondft) 470 ABI_SFREE(results_gs%gresid) 471 ABI_SFREE(results_gs%grewtn) 472 ABI_SFREE(results_gs%grchempottn) 473 ABI_SFREE(results_gs%grvdw) 474 ABI_SFREE(results_gs%grxc) 475 ABI_SFREE(results_gs%intgres) 476 ABI_SFREE(results_gs%synlgr) 477 478 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
SOURCE
499 subroutine destroy_results_gs_array(results_gs) 500 501 !Arguments ------------------------------------ 502 !arrays 503 type(results_gs_type),intent(inout) :: results_gs(:,:) 504 !Local variables------------------------------- 505 !scalars 506 integer :: ii,jj,results_gs_size1,results_gs_size2 507 508 !************************************************************************ 509 510 !@results_gs_type 511 512 results_gs_size1=size(results_gs,1) 513 results_gs_size2=size(results_gs,2) 514 if (results_gs_size1>0.and.results_gs_size2>0) then 515 516 do ii=1,results_gs_size2 517 do jj=1,results_gs_size1 518 results_gs(jj,ii)%natom =0 519 results_gs(jj,ii)%ngrvdw=0 520 results_gs(jj,ii)%nspden=0 521 results_gs(jj,ii)%nsppol=0 522 results_gs(jj,ii)%berryopt=0 523 524 ABI_SFREE(results_gs(jj,ii)%fcart) 525 ABI_SFREE(results_gs(jj,ii)%gred) 526 ABI_SFREE(results_gs(jj,ii)%gaps) 527 ABI_SFREE(results_gs(jj,ii)%grchempottn) 528 ABI_SFREE(results_gs(jj,ii)%grcondft) 529 ABI_SFREE(results_gs(jj,ii)%gresid) 530 ABI_SFREE(results_gs(jj,ii)%grewtn) 531 ABI_SFREE(results_gs(jj,ii)%grvdw) 532 ABI_SFREE(results_gs(jj,ii)%grxc) 533 ABI_SFREE(results_gs(jj,ii)%intgres) 534 ABI_SFREE(results_gs(jj,ii)%synlgr) 535 end do 536 end do 537 538 end if 539 540 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,gred,strten
OUTPUT
SIDE EFFECTS
results_gs=<type(results_gs_type)>=results_gs datastructure
SOURCE
253 subroutine init_results_gs(natom,nspden,nsppol,results_gs,only_part) 254 255 !Arguments ------------------------------------ 256 !scalars 257 integer,intent(in) :: natom,nspden,nsppol 258 logical,optional,intent(in) :: only_part 259 !arrays 260 type(results_gs_type),intent(inout) :: results_gs 261 !Local variables------------------------------- 262 !scalars 263 logical :: full_init 264 265 !************************************************************************ 266 267 !@results_gs_type 268 269 full_init=.true.;if (present(only_part)) full_init=(.not.only_part) 270 271 results_gs%berryopt=0 272 results_gs%natom =natom 273 results_gs%ngrvdw =0 274 results_gs%nspden =nspden 275 results_gs%nsppol =nsppol 276 277 results_gs%deltae =zero 278 results_gs%diffor =zero 279 results_gs%entropy=zero 280 results_gs%entropy_extfpmd=zero 281 results_gs%etotal =zero 282 results_gs%fermie =zero 283 results_gs%fermih =zero ! CP added for case occopt 9 284 results_gs%nelect_extfpmd=zero 285 results_gs%residm =zero 286 results_gs%res2 =zero 287 results_gs%shiftfactor_extfpmd=zero 288 results_gs%vxcavg =zero 289 290 call energies_init(results_gs%energies) 291 292 results_gs%strten=zero 293 ABI_MALLOC(results_gs%fcart,(3,natom)) 294 results_gs%fcart=zero 295 ABI_MALLOC(results_gs%gred,(3,natom)) 296 results_gs%gred =zero 297 ABI_MALLOC(results_gs%gaps,(3,nsppol)) 298 results_gs%gaps =zero 299 ABI_MALLOC(results_gs%intgres,(nspden,natom)) 300 results_gs%intgres=zero 301 302 if (full_init) then 303 results_gs%pel=zero 304 results_gs%pion=zero 305 306 ABI_MALLOC(results_gs%grchempottn,(3,natom)) 307 results_gs%grchempottn=zero 308 ABI_MALLOC(results_gs%grcondft,(3,natom)) 309 results_gs%grcondft=zero 310 ABI_MALLOC(results_gs%gresid,(3,natom)) 311 results_gs%gresid=zero 312 ABI_MALLOC(results_gs%grewtn,(3,natom)) 313 results_gs%grewtn=zero 314 ABI_MALLOC(results_gs%grvdw,(3,natom)) 315 results_gs%grvdw=zero 316 ABI_MALLOC(results_gs%grxc,(3,natom)) 317 results_gs%grxc =zero 318 ABI_MALLOC(results_gs%synlgr,(3,natom)) 319 results_gs%synlgr=zero 320 end if 321 322 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,gred,strten
OUTPUT
SIDE EFFECTS
results_gs(:)=<type(results_gs_type)>=results_gs datastructure 2Darray
SOURCE
348 subroutine init_results_gs_array(natom,nspden,nsppol,results_gs,only_part) 349 350 !Arguments ------------------------------------ 351 !scalars 352 integer,intent(in) :: natom,nspden,nsppol 353 logical,optional,intent(in) :: only_part 354 !arrays 355 type(results_gs_type),intent(inout) :: results_gs(:,:) 356 !Local variables------------------------------- 357 !scalars 358 integer :: ii,jj,results_gs_size1,results_gs_size2 359 logical :: full_init 360 !arrays 361 362 !************************************************************************ 363 364 !@results_gs_type 365 366 results_gs_size1=size(results_gs,1) 367 results_gs_size2=size(results_gs,2) 368 full_init=.true.;if (present(only_part)) full_init=(.not.only_part) 369 370 if (results_gs_size1>0.and.results_gs_size2>0) then 371 372 do ii=1,results_gs_size2 373 do jj=1,results_gs_size1 374 375 results_gs(jj,ii)%berryopt=0 376 results_gs(jj,ii)%natom =natom 377 results_gs(jj,ii)%ngrvdw =0 378 results_gs(jj,ii)%nspden =nspden 379 results_gs(jj,ii)%nsppol =nsppol 380 381 results_gs(jj,ii)%deltae =zero 382 results_gs(jj,ii)%diffor =zero 383 results_gs(jj,ii)%entropy=zero 384 results_gs(jj,ii)%entropy_extfpmd=zero 385 results_gs(jj,ii)%etotal =zero 386 results_gs(jj,ii)%fermie =zero 387 results_gs(jj,ii)%fermih =zero ! CP added for occopt 9 cases 388 results_gs(jj,ii)%nelect_extfpmd=zero 389 results_gs(jj,ii)%residm =zero 390 results_gs(jj,ii)%res2 =zero 391 results_gs(jj,ii)%shiftfactor_extfpmd=zero 392 results_gs(jj,ii)%vxcavg =zero 393 394 call energies_init(results_gs(jj,ii)%energies) 395 396 results_gs(jj,ii)%strten=zero 397 ABI_MALLOC(results_gs(jj,ii)%fcart,(3,natom)) 398 results_gs(jj,ii)%fcart=zero 399 ABI_MALLOC(results_gs(jj,ii)%gred,(3,natom)) 400 results_gs(jj,ii)%gred =zero 401 ABI_MALLOC(results_gs(jj,ii)%gaps,(3,nsppol)) 402 results_gs(jj,ii)%gaps =zero 403 ABI_MALLOC(results_gs(jj,ii)%intgres,(nspden,natom)) 404 results_gs(jj,ii)%intgres =zero 405 406 if (full_init) then 407 results_gs(jj,ii)%pel=zero 408 results_gs(jj,ii)%pion=zero 409 ABI_MALLOC(results_gs(jj,ii)%grchempottn,(3,natom)) 410 results_gs(jj,ii)%grchempottn=zero 411 ABI_MALLOC(results_gs(jj,ii)%grcondft,(3,natom)) 412 results_gs(jj,ii)%grcondft=zero 413 ABI_MALLOC(results_gs(jj,ii)%gresid,(3,natom)) 414 results_gs(jj,ii)%gresid=zero 415 ABI_MALLOC(results_gs(jj,ii)%grewtn,(3,natom)) 416 results_gs(jj,ii)%grewtn=zero 417 ABI_MALLOC(results_gs(jj,ii)%grxc,(3,natom)) 418 results_gs(jj,ii)%grxc =zero 419 ABI_MALLOC(results_gs(jj,ii)%grvdw,(3,results_gs(jj,ii)%ngrvdw)) 420 results_gs(jj,ii)%grvdw =zero 421 ABI_MALLOC(results_gs(jj,ii)%synlgr,(3,natom)) 422 results_gs(jj,ii)%synlgr=zero 423 end if 424 425 end do 426 end do 427 end if 428 429 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
SOURCE
696 integer function results_gs_ncwrite(res, ncid, ecut, pawecutdg) result(ncerr) 697 698 !Arguments ------------------------------------ 699 !scalars 700 integer,intent(in) :: ncid 701 real(dp),intent(in) :: ecut,pawecutdg 702 class(results_gs_type),intent(in) :: res 703 704 !Local variables------------------------------- 705 !scalars 706 #ifdef HAVE_NETCDF 707 708 ! ************************************************************************* 709 710 ! ============================================== 711 ! === Write the dimensions specified by ETSF === 712 ! ============================================== 713 !FIXME: do not handle k_dependent = 1 714 ncerr = nctk_def_dims(ncid, [nctkdim_t("six", 6), nctkdim_t("number_of_cartesian_dimensions", 3), & 715 nctkdim_t("number_of_atoms", res%natom), nctkdim_t("number_of_spins", res%nsppol)], defmode=.True.) 716 NCF_CHECK(ncerr) 717 718 ! Define variables. 719 ! scalars passed in input (not belonging to results_gs) as well as scalars defined in results_gs 720 ! CP modified 721 !ncerr = nctk_def_dpscalars(ncid, [character(len=nctk_slen) :: & 722 ! "ecut", "pawecutdg", "deltae", "diffor", "entropy", "etotal", "fermie", "residm", "res2"]) 723 ncerr = nctk_def_dpscalars(ncid, [character(len=nctk_slen) :: & 724 "ecut", "pawecutdg", "deltae", "diffor", "entropy", "entropy_extfpmd", "etotal", "fermie", "fermih",& 725 & "nelect_extfpmd", "residm", "res2", "shiftfactor_extfpmd"]) ! CP added fermih 726 ! End CP modified 727 NCF_CHECK(ncerr) 728 729 ! arrays 730 ! 731 ! Note: unlike gred, this array has been corrected by enforcing 732 ! the translational symmetry, namely that the sum of force on all atoms is zero. 733 734 ncerr = nctk_def_arrays(ncid, [& 735 nctkarr_t('cartesian_forces', "dp", "number_of_cartesian_dimensions, number_of_atoms"),& 736 nctkarr_t('cartesian_stress_tensor', "dp", 'six')]) 737 NCF_CHECK(ncerr) 738 739 ! In case of a Berry phase calculation output the polarization 740 if (res%berryopt/=0) then 741 ncerr = nctk_def_arrays(ncid, [& 742 nctkarr_t('reduced_electronic_polarization', "dp", "number_of_cartesian_dimensions"),& 743 nctkarr_t('reduced_ionic_polarization', "dp", "number_of_cartesian_dimensions")]) 744 NCF_CHECK(ncerr) 745 end if 746 747 ! Write data. 748 ! Write variables 749 ! CP modified 750 ! ncerr = nctk_write_dpscalars(ncid, [character(len=nctk_slen) :: & 751 !& 'ecut', 'pawecutdg', 'deltae', 'diffor', 'entropy', 'etotal', 'fermie', 'residm', 'res2'],& 752 !& [ecut, pawecutdg, res%deltae, res%diffor, res%entropy, res%etotal, res%fermie, res%residm, res%res2],& 753 !& datamode=.True.) 754 ncerr = nctk_write_dpscalars(ncid, [character(len=nctk_slen) :: & 755 & 'ecut', 'pawecutdg', 'deltae', 'diffor', 'entropy', 'entropy_extfpmd', 'etotal', 'fermie', 'fermih',& 756 & 'nelect_extfpmd', 'residm', 'res2', 'shiftfactor_extfpmd'],& 757 & [ecut, pawecutdg, res%deltae, res%diffor, res%entropy, res%entropy_extfpmd, res%etotal, res%fermie, res%fermih,& 758 & res%nelect_extfpmd, res%residm, res%res2, res%shiftfactor_extfpmd],& 759 & datamode=.True.) 760 ! End CP modified 761 NCF_CHECK(ncerr) 762 763 NCF_CHECK(nctk_set_datamode(ncid)) 764 NCF_CHECK(nf90_put_var(ncid, vid("cartesian_forces"), res%fcart)) 765 NCF_CHECK(nf90_put_var(ncid, vid("cartesian_stress_tensor"), res%strten)) 766 767 if (res%berryopt/=0) then 768 NCF_CHECK(nf90_put_var(ncid, vid("reduced_electronic_polarization"), res%pel)) 769 NCF_CHECK(nf90_put_var(ncid, vid("reduced_ionic_polarization"), res%pion)) 770 end if 771 772 ! Add energies 773 call energies_ncwrite(res%energies, ncid) 774 775 #else 776 ABI_ERROR("netcdf support is not activated.") 777 #endif 778 779 contains 780 integer function vid(vname) 781 782 character(len=*),intent(in) :: vname 783 vid = nctk_idname(ncid, vname) 784 end function vid 785 786 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
61 type, public :: results_gs_type 62 63 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 64 ! declared in another part of ABINIT, that might need to take into account your modification. 65 66 ! Integer scalar 67 68 integer :: natom 69 ! The number of atoms for this dataset 70 71 integer :: nspden 72 ! The number of spin-density components for this dataset 73 74 integer :: nsppol 75 ! The number of spin channels for this dataset 76 77 integer :: ngrvdw 78 ! Size of grvdw array 79 ! Can be 0 (not allocated) or natom 80 81 integer :: berryopt 82 ! Store the Berry phase option to know whether to use pel and pion (0 means no, otherwise yes) 83 84 ! Real (real(dp)) scalars 85 86 real(dp) :: deltae 87 ! change in energy (Hartree) 88 89 real(dp) :: diffor 90 ! maximal absolute value of changes in the components of force 91 92 real(dp) :: nelect_extfpmd 93 ! Contribution of the Extended FPMD model to the number of electrons for high temperature simulations 94 95 ! All the energies are in Hartree, obtained "per unit cell". 96 type(energies_type) :: energies 97 !!! real(dp) :: eei ! local pseudopotential energy (Hartree) 98 !!! real(dp) :: eeig ! sum of eigenvalue energy (Hartree) 99 !!! real(dp) :: eew ! Ewald energy (Hartree) 100 !!! real(dp) :: ehart ! Hartree part of total energy (Hartree) 101 !!! real(dp) :: eii ! pseudopotential core-core energy 102 !!! real(dp) :: ek ! kinetic energy (Hartree) 103 !!! real(dp) :: enefield ! the term of the energy functional that depends 104 !!! ! explicitely on the electric field 105 !!! ! enefield = -ucvol*E*P 106 !!! real(dp) :: enl ! nonlocal pseudopotential energy (Hartree) 107 real(dp) :: entropy ! entropy (Hartree) 108 !!! real(dp) :: enxc ! exchange-correlation energy (Hartree) 109 !!! real(dp) :: enxcdc ! exchange-correlation double-counting energy (Hartree) 110 !!! real(dp) :: epaw ! PAW spherical energy (Hartree) 111 !!! real(dp) :: epawdc ! PAW spherical double-counting energy (Hartree) 112 real(dp) :: entropy_extfpmd ! Entropy contribution of the Extended FPMD model 113 ! for high temperature simulations 114 real(dp) :: etotal ! total energy (Hartree) 115 ! for fixed occupation numbers (occopt==0,1,or 2): 116 ! etotal=ek+ehart+enxc+eei+eew+eii+enl+PAW_spherical_part 117 ! for varying occupation numbers (occopt>=3): 118 ! etotal=ek+ehart+enxc+eei+eew+eii+enl - tsmear*entropy +PAW_spherical_part 119 real(dp) :: fermie ! Fermi energy (Hartree) 120 real(dp) :: fermih ! Fermi energy (Hartree) for excited holes in case occopt 9 (CP added) 121 real(dp) :: residm ! maximum value for the residual over all bands, all k points, 122 ! and all spins (Hartree or Hartree**2, to be checked !) 123 real(dp) :: res2 ! density/potential residual (squared) 124 real(dp) :: vxcavg ! Average of the exchange-correlation energy. The average 125 ! of the local psp pot and the Hartree pot is set to zero (due 126 ! to the usual problem at G=0 for Coulombic system, so vxcavg 127 ! is also the average of the local part of the Hamiltonian 128 129 ! Real (real(dp)) arrays 130 131 real(dp), allocatable :: fcart(:,:) 132 ! fcart(3,natom) 133 ! Cartesian forces (Hartree/Bohr) 134 ! Note: unlike gred, this array has been corrected by enforcing 135 ! the translational symmetry, namely that the sum of force 136 ! on all atoms is zero. 137 138 real(dp), allocatable :: gred(:,:) 139 ! gred(3,natom) 140 ! Forces in reduced coordinates (Hartree) 141 ! Actually, gradient of the total energy with respect 142 ! to change of reduced coordinates 143 144 real(dp),allocatable :: gaps(:,:) 145 ! gaps(3,nsppol) 146 ! gaps(1,:) : fundamental gap 147 ! gaps(2,:) : optical gap 148 ! gaps(3,:) : "status" for each channel: 149 ! 0.0dp if the gap was not computed (because there are only valence bands); 150 ! -1.0dp if the system (or spin-channel) is metallic; 151 ! 1.0dp if the gap was computed 152 153 real(dp), allocatable :: grchempottn(:,:) 154 ! grchempottn(3,natom) 155 ! Part of the gradient of the total energy (Hartree) with respect 156 ! to change of reduced coordinates, that comes from the spatially-varying chemical potential 157 158 real(dp), allocatable :: grcondft(:,:) 159 ! grcondft(nspden,natom) 160 ! Part of the gradient of the total energy (Hartree) with respect 161 ! to change of reduced coordinates, that comes from the constrained DFT contribution 162 163 real(dp), allocatable :: gresid(:,:) 164 ! gresid(3,natom) 165 ! Part of the gradient of the total energy (Hartree) with respect 166 ! to change of reduced coordinates, that comes from the residual 167 ! of the potential 168 169 real(dp), allocatable :: grewtn(:,:) 170 ! grewtn(3,natom) 171 ! Part of the gradient of the total energy (Hartree) with respect 172 ! to change of reduced coordinates, that comes from the Ewald energy 173 174 real(dp), allocatable :: grvdw(:,:) 175 ! grvdw(3,natom) 176 ! Part of the gradient of the total energy (Hartree) with respect 177 ! to change of reduced coordinates, that comes from 178 ! Van der Waals DFT-D dispersion (hartree) 179 ! ngrvdw can be 0 or natom 180 181 real(dp), allocatable :: grxc(:,:) 182 ! grxc(3,natom) 183 ! Part of the gradient of the total energy (Hartree) with respect 184 ! to change of reduced coordinates, that comes from the XC energy 185 186 real(dp), allocatable :: intgres(:,:) 187 ! intgres(nspden,natom) 188 ! Derivative of the total energy with respect to changes of constraints, in constrained DFT. 189 190 real(dp) :: pel(3) 191 ! ucvol times the electronic polarization in reduced coordinates 192 193 real(dp) :: pion(3) 194 ! ucvol times the ionic polarization in reduced coordinates 195 196 real(dp) :: shiftfactor_extfpmd 197 ! Energy shift factor of the Extended FPMD model for high temperature simulations 198 199 real(dp) :: strten(6) 200 ! Stress tensor in cartesian coordinates (Hartree/Bohr^3) 201 ! 6 unique components of this symmetric 3x3 tensor: 202 ! Given in order (1,1), (2,2), (3,3), (3,2), (3,1), (2,1). 203 204 real(dp), allocatable :: synlgr(:,:) 205 ! synlgr(3,natom) 206 ! Part of the gradient of the total energy (Hartree) with respect 207 ! to change of reduced coordinates, that comes from the non-local energy 208 ! The "sy" prefix refer to the fact that this gradient has been 209 ! symmetrized. 210 211 contains 212 213 procedure :: yaml_write => results_gs_yaml_write 214 ! Write the most important results in Yaml format. 215 216 end type results_gs_type 217 218 !public procedures. 219 public :: init_results_gs 220 public :: init_results_gs_array 221 public :: destroy_results_gs 222 public :: destroy_results_gs_array 223 public :: copy_results_gs 224 public :: results_gs_ncwrite
m_results_gs/results_gs_yaml_write [ Functions ]
[ Top ] [ m_results_gs ] [ Functions ]
NAME
results_gs_yaml_write
FUNCTION
Write results_gs in yaml format to unit.
INPUTS
results <type(results_gs_type)>=miscellaneous information about the system after ground state computation unit= unit of output file [cryst]: optional Crystal structure [info]: optional info for the final document [occopt]: optional Input variable occopt [with_conv]: optional True if the convergence dictionary with residuals and diffs should be written.
SOURCE
807 ! CP modified argument list 808 !subroutine results_gs_yaml_write(results, unit, cryst, with_conv, info) 809 subroutine results_gs_yaml_write(results, unit, cryst, info, occopt, with_conv) 810 ! End CP modified 811 812 class(results_gs_type),intent(in) :: results 813 integer,intent(in) :: unit 814 type(crystal_t),intent(in),optional :: cryst 815 integer, intent(in),optional :: occopt ! CP added for special output in the case occopt 9 816 logical,intent(in),optional :: with_conv 817 character(len=*),intent(in),optional :: info 818 819 !Local variables------------------------------- 820 integer,parameter :: width=10 821 integer :: ii 822 type(yamldoc_t) :: ydoc 823 !arrays 824 real(dp) :: strten(3,3), abc(3), fnorms(results%natom) 825 character(len=2) :: species(results%natom) 826 827 !************************************************************************ 828 829 if (unit == dev_null) return 830 831 if (present(info)) then 832 ydoc = yamldoc_open('ResultsGS', info=info, width=width) 833 else 834 ydoc = yamldoc_open('ResultsGS', width=width) 835 end if 836 837 ! Write lattice parameters 838 if (present(cryst)) then 839 call ydoc%add_real2d('lattice_vectors', cryst%rprimd, real_fmt="(f11.7)") 840 abc = [(sqrt(sum(cryst%rprimd(:, ii) ** 2)), ii=1,3)] 841 call ydoc%add_real1d('lattice_lengths', abc, real_fmt="(f10.5)") 842 call ydoc%add_real1d('lattice_angles', cryst%angdeg, real_fmt="(f7.3)", comment="degrees, (23, 13, 12)") 843 call ydoc%add_real('lattice_volume', cryst%ucvol + tol10, real_fmt="(es15.7)") 844 endif 845 846 ! Write convergence degree. 847 ! It seems there's a portability problem for residm computed with nstep = 0 and iscf -3 848 ! because one may get very small value e.g. 7.91-323. residm with nstep > 0 are OK though 849 ! so print zero if residm < tol30 or allow the caller not to write the convergence dict. 850 if(present(with_conv))then 851 if (with_conv) then 852 call ydoc%add_reals( & 853 "deltae, res2, residm, diffor", & 854 [results%deltae, results%res2, merge(results%residm, zero, results%residm > tol30), results%diffor], & 855 real_fmt="(es10.3)", dict_key="convergence") 856 else 857 call ydoc%set_keys_to_string("deltae, res2, residm, diffor", "null", dict_key="convergence") 858 end if 859 else 860 call ydoc%set_keys_to_string("deltae, res2, residm, diffor", "null", dict_key="convergence") 861 endif 862 863 ! Write energies. 864 ! CP modified 865 !call ydoc%add_reals("etotal, entropy, fermie", [results%etotal, results%entropy, results%fermie]) 866 ! CP add modificatopm for occopt 9: fermih 867 if (present(occopt))then 868 if (occopt == 9) then 869 call ydoc%add_reals("etotal, entropy, fermie, fermih", [results%etotal, results%entropy, results%fermie, results%fermih]) 870 else 871 call ydoc%add_reals("etotal, entropy, fermie", [results%etotal, results%entropy, results%fermie]) 872 endif 873 else 874 call ydoc%add_reals("etotal, entropy, fermie", [results%etotal, results%entropy, results%fermie]) 875 endif 876 ! End CP modified 877 878 ! Cartesian stress tensor and forces. 879 call stress_voigt_to_mat(results%strten, strten) 880 if (strten(1,1) /= MAGIC_UNDEF) then 881 call ydoc%add_real2d('cartesian_stress_tensor', strten, comment="hartree/bohr^3") 882 call ydoc%add_real('pressure_GPa', - get_trace(strten) * HaBohr3_GPa / three, real_fmt="(es12.4)") 883 else 884 call ydoc%set_keys_to_string("cartesian_stress_tensor, pressure_GPa", "null") 885 end if 886 887 if(present(cryst))then 888 species = [(cryst%symbol_iatom(ii), ii=1,cryst%natom)] 889 call ydoc%add_real2d('xred', cryst%xred, slist=species, real_fmt="(es12.4)") 890 endif 891 892 if (results%fcart(1,1) /= MAGIC_UNDEF) then 893 call ydoc%add_real2d('cartesian_forces', results%fcart, comment="hartree/bohr") 894 fnorms = [(sqrt(sum(results%fcart(:, ii) ** 2)), ii=1,results%natom)] 895 ! Write force statistics 896 call ydoc%add_reals('min, max, mean', & 897 values=[minval(fnorms), maxval(fnorms), sum(fnorms) / results%natom], dict_key="force_length_stats") 898 else 899 ! Set entries to null (python None) to facilitate life to the parsing routines! 900 call ydoc%add_string('cartesian_forces', "null") 901 call ydoc%set_keys_to_string("min, max, mean", "null", dict_key="force_length_stats") 902 end if 903 904 call ydoc%write_and_free(unit) 905 906 end subroutine results_gs_yaml_write