TABLE OF CONTENTS


ABINIT/m_multibinit_dataset [ Modules ]

[ Top ] [ Modules ]

NAME

  m_multibinit_dataset

FUNCTION

  module with the type for the input of multibinit (should be clean)

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (AM)
  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

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_multibinit_dataset
28 
29  use defs_basis
30  use m_abicore
31  use m_errors
32 
33  use m_parser, only : intagm
34  use m_ddb,    only : DDB_QTOL
35 
36  implicit none
37 
38  private
39 
40  public :: multibinit_dtset_type
41  public :: multibinit_dtset_init
42  public :: multibinit_dtset_free
43  public :: outvars_multibinit
44  public :: invars10

m_multibinit_dataset/invars10 [ Functions ]

[ Top ] [ m_multibinit_dataset ] [ Functions ]

NAME

 invars10

FUNCTION

 Open input file for the multibinit code, then reads or echoes the input information.

INPUTS

 lenstr=actual length of string
 natom=number of atoms, needed for atifc
 string*(*)=string of characters containing all input variables and data

OUTPUT

 multibinit_dtset <type(multibinit_dtset_type)> = datatype with all the input variables

NOTES

 Should be executed by one processor only.

PARENTS

      multibinit

CHILDREN

SOURCE

 496 subroutine invars10(multibinit_dtset,lenstr,natom,string)
 497 
 498 
 499 !This section has been created automatically by the script Abilint (TD).
 500 !Do not modify the following lines by hand.
 501 #undef ABI_FUNC
 502 #define ABI_FUNC 'invars10'
 503 !End of the abilint section
 504 
 505  implicit none
 506 
 507 !Arguments -------------------------------
 508 !scalars
 509  integer,intent(in) :: lenstr,natom
 510  character(len=*),intent(in) :: string
 511  type(multibinit_dtset_type),intent(inout) :: multibinit_dtset
 512 
 513 !Local variables -------------------------
 514 !Dummy arguments for subroutine 'intagm' to parse input file
 515 !Set routine version number here:
 516 !scalars
 517  integer :: iatifc,ii,iph1,iph2,jdtset,jj,marr,tread
 518  character(len=500) :: message
 519 !arrays
 520  integer,allocatable :: intarr(:)
 521  real(dp),allocatable :: dprarr(:),work(:)
 522 
 523 !*********************************************************************
 524  marr=30
 525  ABI_ALLOCATE(intarr,(marr))
 526  ABI_ALLOCATE(dprarr,(marr))
 527 
 528  jdtset=1
 529 
 530 !copy natom to multibinit_dtset
 531  multibinit_dtset%natom=natom
 532 
 533 !=====================================================================
 534 !start reading in dimensions and non-dependent variables
 535 !=====================================================================
 536 
 537 !A
 538  multibinit_dtset%asr=2
 539  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'asr',tread,'INT')
 540  if(tread==1) multibinit_dtset%asr=intarr(1)
 541  if(multibinit_dtset%asr<-2.or.multibinit_dtset%asr>5)then
 542    write(message, '(a,i8,a,a,a,a,a)' )&
 543 &   'asr is',multibinit_dtset%asr,', but the only allowed values',ch10,&
 544 &   'are 0, 1, 2, 3, 4, 5, -1 or -2 .',ch10,&
 545 &   'Action: correct asr in your input file.'
 546 !  Note : negative values are allowed when the acoustic sum rule
 547 !  is to be applied after the analysis of IFCs
 548 !  3,4 are for rotational invariance (under development)
 549 !  5 is for hermitian imposition of the ASR
 550    MSG_ERROR(message)
 551  end if
 552 
 553 !B
 554  multibinit_dtset%brav=1
 555  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'brav',tread,'INT')
 556  if(tread==1) multibinit_dtset%brav=intarr(1)
 557  if(multibinit_dtset%brav/=1)then
 558    write(message, '(a,i8,a,a,a,a,a)' )&
 559 &   'brav is',multibinit_dtset%brav,', but the only allowed values',ch10,&
 560 &   'are 1 for multibinit (not implemented) .',ch10,&
 561 &   'Action: correct brav in your input file.'
 562    MSG_ERROR(message)
 563  end if
 564 
 565  multibinit_dtset%bmass=0
 566  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bmass',tread,'DPR')
 567  if(tread==1) multibinit_dtset%bmass=dprarr(1)
 568  if(multibinit_dtset%bmass<0)then
 569    write(message, '(a,f10.2,a,a,a,a,a)' )&
 570 &   'bmass is',multibinit_dtset%bmass,', but the only allowed values',ch10,&
 571 &   'is superior to 0.',ch10,&
 572 &   'Action: correct bmass in your input file.'
 573    MSG_ERROR(message)
 574  end if
 575 
 576 
 577 !C
 578  multibinit_dtset%chneut=0
 579  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chneut',tread,'INT')
 580  if(tread==1) multibinit_dtset%chneut=intarr(1)
 581  if(multibinit_dtset%chneut<0.or.multibinit_dtset%chneut>2)then
 582    write(message, '(a,i8,a,a,a,a,a)' )&
 583 &   'chneut is',multibinit_dtset%chneut,', but the only allowed values',ch10,&
 584 &   'are 0, 1 or 2 .',ch10,&
 585 &   'Action: correct chneut in your input file.'
 586    MSG_ERROR(message)
 587  end if
 588 
 589  multibinit_dtset%confinement=0
 590  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'confinement',tread,'INT')
 591  if(tread==1) multibinit_dtset%confinement=intarr(1)
 592  if(multibinit_dtset%confinement<0.or.multibinit_dtset%confinement>2)then
 593    write(message, '(a,i8,a,a,a,a,a)' )&
 594 &   'confinement is',multibinit_dtset%confinement,', but the only allowed values',ch10,&
 595 &   'are 0, 1 or 2 .',ch10,&
 596 &   'Action: correct confinement in your input file.'
 597    MSG_ERROR(message)
 598  end if
 599 
 600  multibinit_dtset%conf_power_disp=0
 601  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'conf_power_disp',tread,'INT')
 602  if(tread==1) multibinit_dtset%conf_power_disp=intarr(1)
 603  if(multibinit_dtset%conf_power_disp<0)then
 604    write(message, '(a,i8,a,a,a,a,a)' )&
 605 &   'conf_power_disp is',multibinit_dtset%conf_power_disp,', but the only allowed values',ch10,&
 606 &   'positive .',ch10,&
 607 &   'Action: correct conf_power_disp in your input file.'
 608    MSG_ERROR(message)
 609  end if
 610 
 611  multibinit_dtset%conf_power_strain=0
 612  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'conf_power_strain',tread,'INT')
 613  if(tread==1) multibinit_dtset%conf_power_strain=intarr(1)
 614  if(multibinit_dtset%conf_power_strain<0)then
 615    write(message, '(a,i8,a,a,a,a,a)' )&
 616 &   'conf_power_strain is',multibinit_dtset%conf_power_strain,', but the only allowed values',ch10,&
 617 &   'are positive .',ch10,&
 618 &   'Action: correct conf_power_strain in your input file.'
 619    MSG_ERROR(message)
 620  end if
 621 
 622  multibinit_dtset%conf_power_fact_disp=100
 623  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'conf_power_fact_disp',tread,'DPR')
 624  if(tread==1) multibinit_dtset%conf_power_fact_disp=dprarr(1)
 625 
 626  multibinit_dtset%conf_power_fact_strain=100
 627  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'conf_power_fact_strain',tread,'DPR')
 628  if(tread==1) multibinit_dtset%conf_power_fact_strain=dprarr(1)
 629 
 630 !D
 631  multibinit_dtset%dipdip=1
 632  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dipdip',tread,'INT')
 633  if(tread==1) multibinit_dtset%dipdip=intarr(1)
 634  if(multibinit_dtset%dipdip>1.or.multibinit_dtset%dipdip<0)then
 635    write(message, '(a,i8,a,a,a,a,a)' )&
 636 &   'dipdip is',multibinit_dtset%dipdip,', but the only allowed values',ch10,&
 637 &   'is 1.',ch10,&
 638 &   'Action: correct dipdip in your input file.'
 639    MSG_ERROR(message)
 640  end if
 641 
 642  multibinit_dtset%dipdip_prt=0
 643  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dipdip_prt',tread,'INT')
 644  if(tread==1) multibinit_dtset%dipdip_prt=intarr(1)
 645  if(multibinit_dtset%dipdip_prt<0.or.multibinit_dtset%dipdip_prt>1)then
 646    write(message, '(a,i8,a,a,a,a,a)' )&
 647 &   'dipdip_prt is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
 648     'are 0 or 1.',ch10,&
 649 &   'Action: correct dipdip_prt in your input file.'
 650    MSG_ERROR(message)
 651  end if
 652 
 653 
 654  multibinit_dtset%dtion=100
 655  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dtion',tread,'INT')
 656  if(tread==1) multibinit_dtset%dtion=intarr(1)
 657  if(multibinit_dtset%dtion<1)then
 658    write(message, '(a,i8,a,a,a,a,a)' )&
 659 &   'dtion is',multibinit_dtset%dtion,', but the only allowed values',ch10,&
 660 &   'is superior to 1.',ch10,&
 661 &   'Action: correct dtion in your input file.'
 662    MSG_ERROR(message)
 663  end if
 664 
 665 
 666  multibinit_dtset%delta_df= 1d-02
 667  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'delta_df',tread,'DPR')
 668  if(tread==1) multibinit_dtset%delta_df=dprarr(1)
 669  if(multibinit_dtset%delta_df<0)then
 670    write(message, '(a,es10.2,a,a,a,a,a)' )&
 671 &   'delta_df is',multibinit_dtset%delta_df,', but the only allowed values',ch10,&
 672 &   'are superior to 0  .',ch10,&
 673 &   'Action: correct delta_df in your input file.'
 674    MSG_ERROR(message)
 675  end if
 676 
 677 !E
 678  multibinit_dtset%energy_reference= zero
 679  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'energy_reference',tread,'DPR')
 680  if(tread==1) multibinit_dtset%energy_reference=dprarr(1)
 681 
 682  multibinit_dtset%enunit=0
 683  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'enunit',tread,'INT')
 684  if(tread==1) multibinit_dtset%enunit=intarr(1)
 685  if(multibinit_dtset%enunit<0.or.multibinit_dtset%enunit>2)then
 686    write(message, '(a,i0,a,a,a,a,a)' )&
 687 &   'enunit is',multibinit_dtset%enunit,', but the only allowed values',ch10,&
 688 &   'are 0, 1 or 2.',ch10,&
 689 &   'Action: correct enunit in your input file.'
 690    MSG_ERROR(message)
 691  end if
 692 
 693 !F
 694  multibinit_dtset%fit_option=0
 695  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_option',tread,'INT')
 696  if(tread==1) multibinit_dtset%fit_option=intarr(1)
 697  if(multibinit_dtset%fit_option<0.or.multibinit_dtset%fit_option>2)then
 698    write(message, '(a,i8,a,a,a,a,a)' )&
 699 &   'fit_option is',multibinit_dtset%fit_option,', but the only allowed values',ch10,&
 700 &   'are 0, 1 or 2 for multibinit.',ch10,&
 701 &   'Action: correct fit_option in your input file.'
 702    MSG_ERROR(message)
 703  end if
 704 
 705 
 706  multibinit_dtset%fit_ncoeff=0
 707  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_ncoeff',tread,'INT')
 708  if(tread==1) multibinit_dtset%fit_ncoeff=intarr(1)
 709  if(multibinit_dtset%fit_ncoeff<0)then
 710    write(message, '(a,i8,a,a,a,a,a)' )&
 711 &   'fit_ncoeff is',multibinit_dtset%fit_ncoeff,', but the only allowed values',ch10,&
 712 &   'are positives for multibinit.',ch10,&
 713 &   'Action: correct fit_ncoeff in your input file.'
 714    MSG_ERROR(message)
 715  end if
 716 
 717  multibinit_dtset%fit_nbancoeff=0
 718  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_nbancoeff',tread,'INT')
 719  if(tread==1) multibinit_dtset%fit_nbancoeff=intarr(1)
 720  if(multibinit_dtset%fit_nbancoeff<0)then
 721    write(message, '(a,i8,a,a,a,a,a)' )&
 722 &   'fit_nbancoeff is',multibinit_dtset%fit_nbancoeff,', but the only allowed values',ch10,&
 723 &   'are 0 or positive values for multibinit.',ch10,&
 724 &   'Action: correct fit_nbancoeff in your input file.'
 725    MSG_ERROR(message)
 726  end if
 727 
 728  multibinit_dtset%fit_nfixcoeff=0
 729  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_nfixcoeff',tread,'INT')
 730  if(tread==1) multibinit_dtset%fit_nfixcoeff=intarr(1)
 731  if(multibinit_dtset%fit_nfixcoeff<-2)then
 732    write(message, '(a,i8,a,a,a,a,a)' )&
 733 &   'fit_nfixcoeff is',multibinit_dtset%fit_nfixcoeff,', but the only allowed values',ch10,&
 734 &   'are -1 or positives for multibinit.',ch10,&
 735 &   'Action: correct fit_nfixcoeff in your input file.'
 736    MSG_ERROR(message)
 737  end if
 738 
 739  multibinit_dtset%ts_option=0
 740  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ts_option',tread,'INT')
 741  if(tread==1) multibinit_dtset%ts_option=intarr(1)
 742  if(multibinit_dtset%ts_option<0.or.multibinit_dtset%ts_option>1)then
 743    write(message, '(a,i8,a,a,a,a,a)' )&
 744 &   'ts_option is',multibinit_dtset%ts_option,', but the only allowed values',ch10,&
 745 &   'are positives for multibinit.',ch10,&
 746 &   'Action: correct ts_option in your input file.'
 747    MSG_ERROR(message)
 748  end if
 749 
 750  multibinit_dtset%ifcana=0
 751  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifcana',tread,'INT')
 752  if(tread==1) multibinit_dtset%ifcana=intarr(1)
 753  if(multibinit_dtset%ifcana<0.or.multibinit_dtset%ifcana>1)then
 754    write(message, '(a,i0,a,a,a,a,a)' )&
 755 &   'ifcana is',multibinit_dtset%ifcana,', but the only allowed values',ch10,&
 756 &   'are 0 or 1.',ch10,&
 757 &   'Action: correct ifcana in your input file.'
 758    MSG_ERROR(message)
 759  end if
 760 
 761  multibinit_dtset%ifcflag=1
 762  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifcflag',tread,'INT')
 763  if(tread==1) multibinit_dtset%ifcflag=intarr(1)
 764  if(multibinit_dtset%ifcflag<0.or.multibinit_dtset%ifcflag>1)then
 765    write(message, '(a,i0,a,a,a,a,a)' )&
 766 &   'ifcflag is',multibinit_dtset%ifcflag,', but the only allowed values',ch10,&
 767 &   'are 0 or 1.',ch10,&
 768 &   'Action: correct ifcflag in your input file.'
 769    MSG_ERROR(message)
 770  end if
 771 
 772  multibinit_dtset%prtsrlr=0
 773  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtsrlr',tread,'INT')
 774  if(tread==1) multibinit_dtset%prtsrlr=intarr(1)
 775  if(multibinit_dtset%prtsrlr<0.or.multibinit_dtset%prtsrlr>1)then
 776    write(message, '(a,i8,a,a,a,a,a)' )&
 777 &   'prtsrlr is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
 778 &   'are 0 or 1.',ch10,&
 779 &   'Action: correct prtsrlr in your input file.'
 780    MSG_ERROR(message)
 781  end if
 782 
 783  multibinit_dtset%ifcout=2000000 ! or -1 -> max number of ifc
 784  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifcout',tread,'INT')
 785  if(tread==1) multibinit_dtset%ifcout=intarr(1)
 786  if(multibinit_dtset%ifcout<-1)then
 787    write(message, '(a,i0,a,a,a)' )&
 788 &   'ifcout is',multibinit_dtset%ifcout,', which is lower than -1 (default = all ifc) .',ch10,&
 789 &   'Action: correct ifcout in your input file.'
 790    MSG_ERROR(message)
 791  end if
 792 
 793  multibinit_dtset%nctime=1
 794  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nctime',tread,'INT')
 795  if(tread==1) multibinit_dtset%nctime=intarr(1)
 796  if(multibinit_dtset%nctime<0)then
 797    write(message, '(a,i0,a,a,a)' )&
 798 &   'nctime is',multibinit_dtset%ntime,', which is lower than 0 .',ch10,&
 799 &   'Action: correct nctime in your input file.'
 800    MSG_ERROR(message)
 801  end if
 802 
 803 
 804  multibinit_dtset%ntime=200
 805  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntime',tread,'INT')
 806  if(tread==1) multibinit_dtset%ntime=intarr(1)
 807  if(multibinit_dtset%ntime<0)then
 808    write(message, '(a,i0,a,a,a)' )&
 809 &   'ntime is',multibinit_dtset%ntime,', which is lower than 0 .',ch10,&
 810 &   'Action: correct ntime in your input file.'
 811    MSG_ERROR(message)
 812  end if
 813 
 814  multibinit_dtset%dynamics=0
 815  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dynamics',tread,'INT')
 816  if(tread==1) multibinit_dtset%dynamics=intarr(1)
 817  if(multibinit_dtset%dynamics/=0.and.multibinit_dtset%dynamics/=6.and.&
 818 &   multibinit_dtset%dynamics/=12.and.multibinit_dtset%dynamics/=13.and.&
 819 &   multibinit_dtset%dynamics/=27.and.&
 820 &   multibinit_dtset%dynamics/=24.and.multibinit_dtset%dynamics/=25) then
 821    write(message, '(a,i8,a,a,a,a,a)' )&
 822 &   'dynamics is',multibinit_dtset%dynamics,', but the only allowed values',ch10,&
 823 &   'are 6,12,24,25 or  13 (see ionmov in abinit documentation).',ch10,&
 824 &   'Action: correct dynamics in your input file.'
 825    MSG_ERROR(message)
 826  end if
 827 
 828 !N
 829  multibinit_dtset%natifc=natom
 830  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natifc',tread,'INT')
 831  if(tread==1) multibinit_dtset%natifc=intarr(1)
 832  if(multibinit_dtset%natifc<0)then
 833    write(message, '(a,i0,a,a,a)' )&
 834 &   'natifc is',multibinit_dtset%natifc,', which is lower than 0 .',ch10,&
 835 &   'Action: correct natifc in your input file.'
 836    MSG_ERROR(message)
 837  end if
 838 
 839  if(multibinit_dtset%natifc>natom)then
 840    write(message, '(a,i0,a,a,a,i0,a,a,a)' )&
 841 &   'The number of atom ifc in the input files',multibinit_dtset%natifc,',',ch10,&
 842 &   'is larger than the number of atoms',natom,'.',ch10,&
 843 &   'Action: change natifc in the input file.'
 844    MSG_ERROR(message)
 845  end if
 846 
 847  multibinit_dtset%ncoeff=0
 848  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ncoeff',tread,'INT')
 849  if(tread==1) multibinit_dtset%ncoeff=intarr(1)
 850  if(multibinit_dtset%ncoeff<0)then
 851    write(message, '(a,i0,a,a,a)' )&
 852 &   'ncoeff is',multibinit_dtset%ncoeff,', which is lower than 0 .',ch10,&
 853 &   'Action: correct natifc in your input file.'
 854    MSG_ERROR(message)
 855  end if
 856 
 857  multibinit_dtset%ng2qpt(:)=0
 858  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ng2qpt',tread,'INT')
 859  if(tread==1) multibinit_dtset%ng2qpt(:)=intarr(1:3)
 860  do ii=1,3
 861    if(multibinit_dtset%ng2qpt(ii)<0)then
 862      write(message, '(a,i0,a,i0,a,a,a,i0,a)' )&
 863 &     'ng2qpt(',ii,') is',multibinit_dtset%ng2qpt(ii),', which is lower than 0 .',ch10,&
 864 &     'Action: correct ng2qpt(',ii,') in your input file.'
 865      MSG_ERROR(message)
 866    end if
 867  end do
 868 
 869  multibinit_dtset%ncell(:)= 1
 870  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ncell',tread,'INT')
 871  if(tread==1) multibinit_dtset%ncell(1:3)=intarr(1:3)
 872  do ii=1,3
 873    if(multibinit_dtset%ncell(ii)<0.or.multibinit_dtset%ncell(ii)>100)then
 874      write(message, '(a,i0,a,i0,3a,i0,a)' )&
 875 &     'ncell(',ii,') is ',multibinit_dtset%ncell(ii),', which is lower than 0 of superior than 50.',&
 876 &     ch10,'Action: correct ncell(',ii,') in your input file.'
 877      MSG_ERROR(message)
 878    end if
 879  end do
 880 
 881  multibinit_dtset%ngqpt(:)= 1
 882  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ngqpt',tread,'INT')
 883  if(tread==1) multibinit_dtset%ngqpt(1:3)=intarr(1:3)
 884  do ii=1,3
 885    if(multibinit_dtset%ngqpt(ii)<0)then
 886      write(message, '(a,i0,a,i0,a,a,a,i0,a)' )&
 887 &     'ngqpt(',ii,') is',multibinit_dtset%ngqpt(ii),', which is lower than 0 .',ch10,&
 888 &     'Action: correct ngqpt(',ii,') in your input file.'
 889      MSG_ERROR(message)
 890    end if
 891  end do
 892 
 893  multibinit_dtset%nph1l=1
 894  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nph1l',tread,'INT')
 895  if(tread==1) multibinit_dtset%nph1l=intarr(1)
 896  if(multibinit_dtset%nph1l<0)then
 897    write(message, '(a,i0,a,a,a)' )&
 898 &   'nph1l is',multibinit_dtset%nph1l,', which is lower than 0 .',ch10,&
 899 &   'Action: correct nph1l in your input file.'
 900    MSG_ERROR(message)
 901  end if
 902 
 903  multibinit_dtset%nph2l=0
 904  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nph2l',tread,'INT')
 905  if(tread==1) multibinit_dtset%nph2l=intarr(1)
 906  if(multibinit_dtset%nph2l<0)then
 907    write(message, '(a,i0,a,a,a)' )&
 908 &   'nph2l is',multibinit_dtset%nph2l,', which is lower than 0 .',ch10,&
 909 &   'Action: correct nph2l in your input file.'
 910    MSG_ERROR(message)
 911  end if
 912 
 913  multibinit_dtset%nqshft=1
 914  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqshft',tread,'INT')
 915  if(tread==1) multibinit_dtset%nqshft=intarr(1)
 916  if(multibinit_dtset%nqshft<0 .or. multibinit_dtset%nqshft==3 .or.&
 917 & multibinit_dtset%nqshft>=5 )then
 918    write(message, '(a,i0,a,a,a,a,a)' )&
 919 &   'nqshft is',multibinit_dtset%nqshft,', but the only allowed values',ch10,&
 920 &   'are 1, 2 or 4 .',ch10,&
 921 &   'Action: correct nqshft in your input file.'
 922    MSG_ERROR(message)
 923  end if
 924 
 925  multibinit_dtset%nnos=0
 926  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nnos',tread,'INT')
 927  if(tread==1) multibinit_dtset%nnos=intarr(1)
 928  if(multibinit_dtset%nnos<0)then
 929    write(message, '(a,i0,a,a,a)' )&
 930 &   'nnos is',multibinit_dtset%nnos,', which is lower than 0',ch10,&
 931 &   'Action: correct nnos in your input file.'
 932    MSG_ERROR(message)
 933  end if
 934 
 935 
 936  multibinit_dtset%nsphere=0
 937  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nsphere',tread,'INT')
 938  if(tread==1) multibinit_dtset%nsphere=intarr(1)
 939  if(multibinit_dtset%nsphere<0)then
 940    write(message, '(a,i0,a,a,a)' )&
 941 &   'nsphere is',multibinit_dtset%nsphere,', which is lower than 0',ch10,&
 942 &   'Action: correct nsphere in your input file.'
 943    MSG_ERROR(message)
 944  end if
 945 
 946 !O
 947  multibinit_dtset%optcell=0
 948  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optcell',tread,'INT')
 949  if(tread==1) multibinit_dtset%optcell=intarr(1)
 950  if(multibinit_dtset%optcell<0.or.multibinit_dtset%optcell>2)then
 951    write(message, '(a,i8,a,a,a,a,a)' )&
 952 &   'optcell is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
 953 &   'are 0, 1 or 2.',ch10,&
 954 &   'Action: correct optcell in your input file.'
 955    MSG_ERROR(message)
 956  end if
 957 
 958 
 959 !P
 960  multibinit_dtset%prt_model=0
 961  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prt_model',tread,'INT')
 962  if(tread==1) multibinit_dtset%prt_model=intarr(1)
 963  if(multibinit_dtset%prt_model<0.or.multibinit_dtset%prt_model>4)then
 964    write(message, '(a,i8,a,a,a,a,a)' )&
 965 &   'prt_model is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
 966 &   'are 0, 1 or 2.',ch10,&
 967 &   'Action: correct prt_model in your input file.'
 968    MSG_ERROR(message)
 969  end if
 970 
 971  multibinit_dtset%prt_phfrq=0
 972  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prt_phfrq',tread,'INT')
 973  if(tread==1) multibinit_dtset%prt_phfrq=intarr(1)
 974  if(multibinit_dtset%prt_phfrq<0.or.multibinit_dtset%prt_phfrq>2)then
 975    write(message, '(a,i8,a,a,a,a,a)' )&
 976 &   'prt_phfrq is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
 977 &   'are 0, 1 or 2.',ch10,&
 978 &   'Action: correct prt_phfrq in your input file.'
 979    MSG_ERROR(message)
 980  end if
 981 
 982  multibinit_dtset%fit_initializeData=1
 983  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_initializeData',tread,'INT')
 984  if(tread==1) multibinit_dtset%fit_initializeData=intarr(1)
 985  if(multibinit_dtset%fit_initializeData<0.or.multibinit_dtset%fit_initializeData>1)then
 986    write(message, '(a,i8,a,a,a,a,a)' )&
 987 &   'fit_initializeData is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
 988 &   'are 0, 1 or 2.',ch10,&
 989 &   'Action: correct fit_initializeData in your input file.'
 990    MSG_ERROR(message)
 991  end if
 992 
 993  
 994  multibinit_dtset%fit_generateCoeff=1
 995  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_generateCoeff',tread,'INT')
 996  if(tread==1) multibinit_dtset%fit_generateCoeff=intarr(1)
 997  if(multibinit_dtset%fit_generateCoeff<0.or.multibinit_dtset%fit_generateCoeff>1)then
 998    write(message, '(a,i8,a,a,a,a,a)' )&
 999 &   'fit_generateCoeff is',multibinit_dtset%prtsrlr,', but the only allowed values',ch10,&
1000 &   'are 0, 1 or 2.',ch10,&
1001 &   'Action: correct fit_generateCoeff in your input file.'
1002    MSG_ERROR(message)
1003  end if
1004 
1005 !Default is no output of the real space IFC to file
1006  multibinit_dtset%prt_ifc = 0
1007  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prt_ifc',tread,'INT')
1008  if(tread==1) multibinit_dtset%prt_ifc = intarr(1)
1009  if(multibinit_dtset%prt_ifc < 0 .or. multibinit_dtset%prt_ifc > 1) then
1010    write(message, '(a,i0,a,a,a,a,a)' )&
1011 &   'prtf_ifc is',multibinit_dtset%prt_ifc,'. The only allowed values',ch10,&
1012 &   'are 0 (no output) or 1 (AI2PS format)',ch10,  &
1013 &   'Action: correct prt_ifc in your input file.'
1014    MSG_ERROR(message)
1015  end if
1016 
1017 !Default is no output of the 3rd derivative
1018  multibinit_dtset%strcpling = -1
1019  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'strcpling',tread,'INT')
1020  if(tread==1) multibinit_dtset%strcpling = intarr(1)
1021  if(multibinit_dtset%strcpling < -1 .or. multibinit_dtset%strcpling > 2) then
1022    write(message, '(a,i0,a,a,a,a,a,a,a)' )&
1023 &   'prtf_3rd is ',multibinit_dtset%strcpling,'. The only allowed values',ch10,&
1024 &   'are 0 (no computation), 1 (only computation)',ch10,&
1025 &   'or 2 (computation and print in xml file)',ch10,  &
1026 &   'Action: correct strcpling in your input file.'
1027    MSG_ERROR(message)
1028  end if
1029 
1030 !Q
1031  multibinit_dtset%qrefine=1 ! default is no refinement
1032  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'qrefine',tread,'INT')
1033  if(tread==1) multibinit_dtset%qrefine = intarr(1:3)
1034  do ii=1,3
1035    if(multibinit_dtset%qrefine(ii) < 1) then
1036      write(message, '(a,3i0,a,a,a,a,a)' )&
1037 &     'qrefine is',multibinit_dtset%qrefine,' The only allowed values',ch10,&
1038 &     'are integers >= 1 giving the refinement of the ngqpt grid',ch10,&
1039 &     'Action: correct qrefine in your input file.'
1040      MSG_ERROR(message)
1041    end if
1042  end do
1043 
1044 !R
1045  multibinit_dtset%restartxf=0
1046  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'restartxf',tread,'INT')
1047  if(tread==1) multibinit_dtset%restartxf=intarr(1)
1048  if(multibinit_dtset%restartxf < -3 .or. multibinit_dtset%restartxf > 0)then
1049    write(message, '(a,i8,a,a,a,a,a)' )&
1050 &   'restartxf is',multibinit_dtset%restartxf,', but the only allowed values',ch10,&
1051 &   'is -2 or 0.',ch10,&
1052 &   'Action: correct restartxf in your input file.'
1053    MSG_ERROR(message)
1054  end if
1055 
1056  multibinit_dtset%rfmeth=1
1057  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfmeth',tread,'INT')
1058  if(tread==1) multibinit_dtset%rfmeth=intarr(1)
1059  if(multibinit_dtset%rfmeth<1.or.multibinit_dtset%rfmeth>2)then
1060    write(message, '(a,i0,a,a,a,a,a)' )&
1061 &   'rfmeth is',multibinit_dtset%rfmeth,', but the only allowed values',ch10,&
1062 &   'are 1 or 2 . ',ch10,&
1063 &   'Action: correct rfmeth in your input file.'
1064    MSG_ERROR(message)
1065  end if
1066 
1067  multibinit_dtset%rifcsph=zero
1068  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rifcsph',tread,'DPR')
1069  if(tread==1) multibinit_dtset%rifcsph=dprarr(1)
1070  if(multibinit_dtset%rifcsph<-tol12)then
1071    write(message, '(a,f10.3,a,a,a)' )&
1072 &   'rifcsph is',multibinit_dtset%rifcsph,', which is lower than zero.',ch10,&
1073 &   'Action: correct rifcsph in your input file.'
1074    MSG_ERROR(message)
1075  end if
1076 
1077 !S
1078 
1079  multibinit_dtset%spin_dipdip=0
1080  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_dipdip',tread,'INT')
1081  if(tread==1) multibinit_dtset%spin_dipdip=intarr(1)
1082  if(multibinit_dtset%spin_dipdip>1.or.multibinit_dtset%spin_dipdip<0)then
1083    write(message, '(a,i8,a,a,a,a,a)' )&
1084 &   'spin_dipdip is',multibinit_dtset%spin_dipdip,', but the only allowed values',ch10,&
1085 &   'is 0 or 1.',ch10,&
1086 &   'Action: correct spin_dipdip in your input file.'
1087    MSG_ERROR(message)
1088  end if
1089 
1090 
1091  multibinit_dtset%spin_dt= 1d-16
1092  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_dt',tread,'DPR')
1093  if(tread==1) multibinit_dtset%spin_dt=dprarr(1)
1094  if(multibinit_dtset%spin_dt<0)then
1095     write(message, '(a,es10.2,a,a,a,a,a)' )&
1096          &   'spin_dt is',multibinit_dtset%spin_dt,', but the only allowed values',ch10,&
1097          &   'are superior to 0  .',ch10,&
1098          &   'Action: correct spin_dt in your input file.'
1099     MSG_ERROR(message)
1100  end if
1101 
1102  
1103  multibinit_dtset%spin_dynamics=0
1104  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_dynamics',tread,'INT')
1105  if(tread==1) multibinit_dtset%spin_dynamics=intarr(1)
1106  if(multibinit_dtset%spin_dynamics/=0.and.&
1107       &   multibinit_dtset%spin_dynamics/=1) then
1108     write(message, '(a,i8,a,a,a,a,a)' )&
1109          &   'spin_dynamics is',multibinit_dtset%spin_dynamics,', but the only allowed values',ch10,&
1110          &   'are 0 and 1.',ch10,&
1111          &   'Action: correct spin_dynamics in your input file.'
1112     MSG_ERROR(message)
1113  end if
1114  
1115  multibinit_dtset%spin_mag_field= zero
1116  if(3>marr)then
1117     marr=3
1118     ABI_DEALLOCATE(intarr)
1119     ABI_DEALLOCATE(dprarr)
1120     ABI_ALLOCATE(intarr,(marr))
1121     ABI_ALLOCATE(dprarr,(marr))
1122  end if
1123  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'spin_mag_field',tread,'DPR')
1124  if(tread==1) multibinit_dtset%spin_mag_field(1:3)= dprarr(1:3)
1125 
1126  multibinit_dtset%spin_nctime=100
1127  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_nctime',tread,'INT')
1128  if(tread==1) multibinit_dtset%spin_nctime=intarr(1)
1129  if(multibinit_dtset%spin_nctime<0)then
1130     write(message, '(a,i0,a,a,a)' )&
1131          &   'spin_nctime is',multibinit_dtset%spin_nctime,', which is lower than 0 .',ch10,&
1132          &   'Action: correct spin_nctime in your input file.'
1133     MSG_ERROR(message)
1134  end if
1135 
1136 
1137  multibinit_dtset%spin_ntime=10000
1138  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_ntime',tread,'INT')
1139  if(tread==1) multibinit_dtset%spin_ntime=intarr(1)
1140  if(multibinit_dtset%spin_ntime<0)then
1141     write(message, '(a,i0,a,a,a)' )&
1142          &   'spin_ntime is',multibinit_dtset%spin_ntime,', which is lower than 0 .',ch10,&
1143          &   'Action: correct spin_ntime in your input file.'
1144     MSG_ERROR(message)
1145  end if
1146 
1147 
1148  multibinit_dtset%spin_n1l=1
1149  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_n1l',tread,'INT')
1150  if(tread==1) multibinit_dtset%spin_n1l=intarr(1)
1151  if(multibinit_dtset%spin_n1l<0)then
1152     write(message, '(a,i0,a,a,a)' )&
1153          &   'spin_n1l is',multibinit_dtset%spin_n1l,', which is lower than 0 .',ch10,&
1154          &   'Action: correct spin_n1l in your input file.'
1155     MSG_ERROR(message)
1156  end if
1157 
1158  multibinit_dtset%spin_n2l=0
1159  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_n2l',tread,'INT')
1160  if(tread==1) multibinit_dtset%spin_n2l=intarr(1)
1161  if(multibinit_dtset%spin_n2l<0)then
1162     write(message, '(a,i0,a,a,a)' )&
1163          &   'spin_n2l is',multibinit_dtset%spin_n2l,', which is lower than 0 .',ch10,&
1164          &   'Action: correct spin_n2l in your input file.'
1165     MSG_ERROR(message)
1166  end if
1167  
1168 
1169  multibinit_dtset%spin_qpoint= zero
1170  if(3>marr)then
1171     marr=3
1172     ABI_DEALLOCATE(intarr)
1173     ABI_DEALLOCATE(dprarr)
1174     ABI_ALLOCATE(intarr,(marr))
1175     ABI_ALLOCATE(dprarr,(marr))
1176  end if
1177  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'spin_qpoint',tread,'DPR')
1178  if(tread==1) multibinit_dtset%spin_mag_field(1:3)= dprarr(1:3)
1179 
1180 
1181  multibinit_dtset%spin_temperature=325
1182  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_temperature',tread,'DPR')
1183  if(tread==1) multibinit_dtset%spin_temperature=dprarr(1)
1184  if(multibinit_dtset%spin_temperature<=0)then
1185    write(message, '(a,f10.1,a,a,a,a,a)' )&
1186 &   'spin_temperature is ',multibinit_dtset%spin_temperature,'. The only allowed values',ch10,&
1187 &   'are positives values.',ch10,&
1188 &   'Action: correct spin_semperature in your input file.'
1189    MSG_ERROR(message)
1190  end if
1191 
1192 
1193  multibinit_dtset%spin_tolavg=1d-02
1194  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_tolavg',tread,'DPR')
1195  if(tread==1) multibinit_dtset%spin_tolavg=dprarr(1)
1196  if(multibinit_dtset%spin_tolavg<=0)then
1197     write(message, '(a,f10.1,a,a,a,a,a)' )&
1198          &   'spin_tolavg is ',multibinit_dtset%spin_tolavg,'. The only allowed values',ch10,&
1199          &   'are positives values.',ch10,&
1200          &   'Action: correct spin_tolavg in your input file.'
1201     MSG_ERROR(message)
1202  end if
1203 
1204  multibinit_dtset%spin_tolvar=1d-02
1205  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spin_tolvar',tread,'DPR')
1206  if(tread==1) multibinit_dtset%spin_tolvar=dprarr(1)
1207  if(multibinit_dtset%spin_tolvar<=0)then
1208     write(message, '(a,f10.1,a,a,a,a,a)' )&
1209          &   'spin_tolvar is ',multibinit_dtset%spin_tolvar,'. The only allowed values',ch10,&
1210          &   'are positives values.',ch10,&
1211          &   'Action: correct spin_tolvar in your input file.'
1212     MSG_ERROR(message)
1213  end if
1214 
1215 
1216  multibinit_dtset%symdynmat=1
1217  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'symdynmat',tread,'INT')
1218  if(tread==1) multibinit_dtset%symdynmat=intarr(1)
1219  if(multibinit_dtset%symdynmat/=0.and.multibinit_dtset%symdynmat/=1)then
1220    write(message, '(a,i0,a,a,a,a,a)' )&
1221 &   'symdynmat is',multibinit_dtset%symdynmat,'. The only allowed values',ch10,&
1222 &   'are 0, or 1.',ch10,&
1223 &   'Action: correct symdynmat in your input file.'
1224    MSG_ERROR(message)
1225  end if
1226 
1227 !T
1228 
1229 
1230  multibinit_dtset%temperature=325
1231  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'temperature',tread,'DPR')
1232  if(tread==1) multibinit_dtset%temperature=dprarr(1)
1233  if(multibinit_dtset%temperature<=0)then
1234    write(message, '(a,f10.1,a,a,a,a,a)' )&
1235 &   'Temperature is ',multibinit_dtset%temperature,'. The only allowed values',ch10,&
1236 &   'are positives values.',ch10,&
1237 &   'Action: correct Temperature in your input file.'
1238    MSG_ERROR(message)
1239  end if
1240 
1241 
1242 
1243 !U
1244 
1245 !V
1246 
1247 !W
1248 
1249 !X
1250 
1251 !Y
1252 
1253 !Z
1254 
1255 !=====================================================================
1256 !end non-dependent variables
1257 !=====================================================================
1258 
1259 !=======================================================================
1260 !Read in dependent variables (dependent on dimensions above)
1261 !=======================================================================
1262 
1263 !A
1264  multibinit_dtset%acell= one
1265  if(3>marr)then
1266    marr=3
1267    ABI_DEALLOCATE(intarr)
1268    ABI_DEALLOCATE(dprarr)
1269    ABI_ALLOCATE(intarr,(marr))
1270    ABI_ALLOCATE(dprarr,(marr))
1271  end if
1272  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'acell',tread,'DPR')
1273  if(tread==1) multibinit_dtset%acell(1:3)= dprarr(1:3)
1274  if(any(multibinit_dtset%acell<=tol10))then
1275     write(message, '(3a)' )&
1276 &       'There is negative or zero value for cell ',ch10,&
1277 &       'Action: change acell in your input file.'
1278       MSG_ERROR(message)
1279  end if
1280 
1281  if(6>marr)then
1282    marr=6
1283    ABI_DEALLOCATE(intarr)
1284    ABI_DEALLOCATE(dprarr)
1285    ABI_ALLOCATE(intarr,(marr))
1286    ABI_ALLOCATE(dprarr,(marr))
1287  end if
1288  multibinit_dtset%strtarget(1:6) = zero
1289  call intagm(dprarr,intarr,jdtset,marr,6,string(1:lenstr),'strtarget',tread,'DPR')
1290  if(tread==1) multibinit_dtset%strtarget(1:6)=dprarr(1:6)
1291 
1292 
1293  ABI_ALLOCATE(multibinit_dtset%atifc,(natom))
1294  multibinit_dtset%atifc(:)=0
1295  if(multibinit_dtset%natifc>=1)then
1296    if(multibinit_dtset%natifc>marr)then
1297      marr=multibinit_dtset%natifc
1298      ABI_DEALLOCATE(intarr)
1299      ABI_DEALLOCATE(dprarr)
1300      ABI_ALLOCATE(intarr,(marr))
1301      ABI_ALLOCATE(dprarr,(marr))
1302    end if
1303    call intagm(dprarr,intarr,jdtset,marr,multibinit_dtset%natifc,string(1:lenstr),'atifc',tread,'INT')
1304    if(tread==1) then
1305      multibinit_dtset%atifc(1:multibinit_dtset%natifc)= intarr(1:multibinit_dtset%natifc)
1306    else ! set to the maximum
1307      do iatifc=1,multibinit_dtset%natifc
1308        multibinit_dtset%atifc(iatifc) =  iatifc
1309      end do
1310    end if
1311    ABI_MALLOC(work,(natom))
1312    work(:)=0
1313 
1314    do iatifc=1,multibinit_dtset%natifc
1315      if(multibinit_dtset%atifc(iatifc)<=0.or.multibinit_dtset%atifc(iatifc)>natom)then
1316        write(message, '(a,i0,a,a,a,a,a,i0,a,a,a)' )&
1317 &       'For iatifc=',iatifc,', the number of the atom ifc to be ',ch10,&
1318 &       'analysed is not valid : either negative, ',ch10,&
1319 &       'zero, or larger than natom =',natom,'.',ch10,&
1320 &       'Action: change atifc in your input file.'
1321        MSG_ERROR(message)
1322      end if
1323      work(multibinit_dtset%atifc(iatifc))=1
1324    end do
1325    multibinit_dtset%atifc(1:natom)=int(work(:))
1326    ABI_FREE(work)
1327  end if
1328 
1329 !B
1330 
1331 !C
1332  ABI_ALLOCATE(multibinit_dtset%coefficients,(multibinit_dtset%ncoeff))
1333  if (multibinit_dtset%ncoeff/=0)then
1334    if(multibinit_dtset%ncoeff>marr)then
1335      marr=multibinit_dtset%ncoeff
1336      ABI_DEALLOCATE(intarr)
1337      ABI_DEALLOCATE(dprarr)
1338      ABI_ALLOCATE(intarr,(marr))
1339      ABI_ALLOCATE(dprarr,(marr))
1340    end if
1341    multibinit_dtset%coefficients(:)=zero
1342    call intagm(dprarr,intarr,jdtset,marr,multibinit_dtset%ncoeff,&
1343 &              string(1:lenstr),'coefficients',tread,'DPR')
1344    if(tread==1)then
1345      do ii=1,multibinit_dtset%ncoeff
1346        multibinit_dtset%coefficients(ii)=dprarr(ii)
1347      end do
1348    end if
1349  end if
1350 
1351  ABI_ALLOCATE(multibinit_dtset%conf_cutoff_disp,(multibinit_dtset%natom))
1352  if (multibinit_dtset%natom/=0)then
1353    if(multibinit_dtset%natom>marr)then
1354      marr=multibinit_dtset%natom
1355      ABI_DEALLOCATE(intarr)
1356      ABI_DEALLOCATE(dprarr)
1357      ABI_ALLOCATE(intarr,(marr))
1358      ABI_ALLOCATE(dprarr,(marr))
1359    end if
1360    multibinit_dtset%conf_cutoff_disp(:)=zero
1361    call intagm(dprarr,intarr,jdtset,marr,multibinit_dtset%natom,&
1362 &              string(1:lenstr),'conf_cutoff_disp',tread,'DPR')
1363    if(tread==1)then
1364      do ii=1,multibinit_dtset%natom
1365        multibinit_dtset%conf_cutoff_disp(ii)=dprarr(ii)
1366      end do
1367    end if
1368    if(any(multibinit_dtset%conf_cutoff_disp<zero))then
1369      write(message, '(3a)' )&
1370 &       'There is negative value for conf_cutoff_disp ',ch10,&
1371 &       'Action: change acell in your input file.'
1372      MSG_ERROR(message)
1373    end if
1374  end if
1375 
1376  if(6>marr)then
1377    marr=6
1378    ABI_DEALLOCATE(intarr)
1379    ABI_DEALLOCATE(dprarr)
1380    ABI_ALLOCATE(intarr,(marr))
1381    ABI_ALLOCATE(dprarr,(marr))
1382  end if
1383  multibinit_dtset%conf_cutoff_strain(1:6) = zero
1384  call intagm(dprarr,intarr,jdtset,marr,6,string(1:lenstr),'conf_cutoff_strain',tread,'DPR')
1385  if(tread==1) multibinit_dtset%conf_cutoff_strain(1:6)=dprarr(1:6)
1386  if(any(multibinit_dtset%conf_cutoff_disp<zero))then
1387    write(message, '(3a)' )&
1388 &     'There is negative value for conf_cutoff_strain ',ch10,&
1389 &     'Action: change acell in your input file.'
1390    MSG_ERROR(message)
1391  end if
1392 
1393 !D
1394  multibinit_dtset%dipdip_range(:)= (/0,0,0/)
1395  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'dipdip_range',tread,'INT')
1396  if(tread==1) multibinit_dtset%dipdip_range(1:3)=intarr(1:3)
1397  do ii=1,3
1398    if(multibinit_dtset%dipdip_range(ii)<0.or.multibinit_dtset%dipdip_range(ii)>50)then
1399      write(message, '(a,i0,a,i0,4a,i0,a)' )&
1400 &     'dipdip_range(',ii,') is ',multibinit_dtset%dipdip_range(ii),', which is lower',&
1401 &     ' than 0 of superior than 50.',&
1402 &     ch10,'Action: correct dipdip_range(',ii,') in your input file.'
1403      MSG_ERROR(message)
1404    end if
1405  end do
1406 !E
1407  multibinit_dtset%eivec=0
1408  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'eivec',tread,'INT')
1409  if(tread==1) multibinit_dtset%eivec=intarr(1)
1410  if(multibinit_dtset%eivec<0.or.multibinit_dtset%eivec>4)then
1411    write(message, '(a,i0,a,a,a,a,a)' )&
1412 &   'eivec is',multibinit_dtset%eivec,', but the only allowed values',ch10,&
1413 &   'are 0, 1, 2, 3 or 4.',ch10,&
1414 &   'Action: correct eivec in your input file.'
1415    MSG_ERROR(message)
1416  end if
1417 
1418 !F
1419   multibinit_dtset%fit_anhaStrain=0
1420  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_anhaStrain',tread,'INT')
1421  if(tread==1) multibinit_dtset%fit_anhaStrain=intarr(1)
1422  if(multibinit_dtset%fit_anhaStrain<0.and.multibinit_dtset%fit_anhaStrain>1)then
1423    write(message, '(a,i8,a,a,a,a,a)' )&
1424 &   'fit_anhaStrain is',multibinit_dtset%fit_anhaStrain,', but the only allowed values',ch10,&
1425 &   'are 0 or 1 for multibinit.',ch10,&
1426 &   'Action: correct fit_anhaStrain in your input file.'
1427    MSG_ERROR(message)
1428  end if
1429 
1430  multibinit_dtset%bound_model=0
1431  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bound_model',tread,'INT')
1432  if(tread==1) multibinit_dtset%bound_model=intarr(1)
1433  if(multibinit_dtset%bound_model<0.and.multibinit_dtset%bound_model>1)then
1434    write(message, '(a,i8,a,a,a,a,a)' )&
1435 &   'bound_model is',multibinit_dtset%bound_model,', but the only allowed values',ch10,&
1436 &   'are 0 or 1 for multibinit.',ch10,&
1437 &   'Action: correct bound_model in your input file.'
1438    MSG_ERROR(message)
1439  end if
1440 
1441  multibinit_dtset%bound_anhaStrain=0
1442  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_anhaStrain',tread,'INT')
1443  if(tread==1) multibinit_dtset%bound_anhaStrain=intarr(1)
1444  if(multibinit_dtset%bound_anhaStrain<0.and.multibinit_dtset%bound_anhaStrain>1)then
1445    write(message, '(a,i8,a,a,a,a,a)' )&
1446 &   'fit_anhaStrain is',multibinit_dtset%bound_anhaStrain,', but the only allowed values',ch10,&
1447 &   'are 0 or 1 for multibinit.',ch10,&
1448 &   'Action: correct fit_anhaStrain in your input file.'
1449    MSG_ERROR(message)
1450  end if
1451 
1452  multibinit_dtset%bound_SPCoupling=1
1453  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bound_SPCoupling',tread,'INT')
1454  if(tread==1) multibinit_dtset%bound_SPCoupling=intarr(1)
1455  if(multibinit_dtset%bound_SPCoupling<0.and.multibinit_dtset%bound_SPCoupling>1)then
1456    write(message, '(a,i8,a,a,a,a,a)' )&
1457      &   'fit_SPCoupling is',multibinit_dtset%bound_SPCoupling,&
1458      &   ', but the only allowed values',ch10,&
1459 &   'are 0 or 1 for multibinit.',ch10,&
1460 &   'Action: correct fit_SPCoupling in your input file.'
1461    MSG_ERROR(message)
1462  end if
1463 
1464  multibinit_dtset%fit_SPCoupling=1
1465  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_SPCoupling',tread,'INT')
1466  if(tread==1) multibinit_dtset%fit_SPCoupling=intarr(1)
1467  if(multibinit_dtset%fit_SPCoupling<0.and.multibinit_dtset%fit_SPCoupling>1)then
1468    write(message, '(a,i8,a,a,a,a,a)' )&
1469      &   'fit_SPCoupling is',multibinit_dtset%fit_SPCoupling,&
1470      &   ', but the only allowed values',ch10,&
1471 &   'are 0 or 1 for multibinit.',ch10,&
1472 &   'Action: correct fit_SPCoupling in your input file.'
1473    MSG_ERROR(message)
1474  end if
1475 
1476  multibinit_dtset%bound_cutoff=0
1477  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bound_cutoff',tread,'DPR')
1478  if(tread==1) multibinit_dtset%bound_cutoff=dprarr(1)
1479  if(multibinit_dtset%bound_cutoff<0)then
1480    write(message, '(a,i8,a,a,a,a,a)' )&
1481 &   'bound_cutoff is',multibinit_dtset%bound_cutoff,', but the only allowed values',ch10,&
1482 &   'are positives for multibinit.',ch10,&
1483 &   'Action: correct bound_cutoff in your input file.'
1484    MSG_ERROR(message)
1485  end if
1486 
1487 
1488   multibinit_dtset%bound_maxCoeff=4
1489  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bound_maxCoeff',tread,'INT')
1490  if(tread==1) multibinit_dtset%bound_maxCoeff=intarr(1)
1491  if(multibinit_dtset%bound_maxCoeff<0.and.multibinit_dtset%bound_maxCoeff>1)then
1492    write(message, '(a,i8,a,a,a,a,a)' )&
1493 &   'bound_maxCoeff is',multibinit_dtset%bound_maxCoeff,', but the only allowed values',ch10,&
1494 &   'are 0 or 1 for multibinit.',ch10,&
1495 &   'Action: correct bound_maxCoeff in your input file.'
1496    MSG_ERROR(message)
1497  end if
1498 
1499   multibinit_dtset%bound_temp=325
1500  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bound_temp',tread,'DPR')
1501  if(tread==1) multibinit_dtset%bound_temp=dprarr(1)
1502  if(multibinit_dtset%bound_temp<=0)then
1503    write(message, '(a,f10.1,a,a,a,a,a)' )&
1504 &   'Bound_Temp is ',multibinit_dtset%bound_temp,'. The only allowed values',ch10,&
1505 &   'are positives values.',ch10,&
1506 &   'Action: correct Bound_Temp in your input file.'
1507    MSG_ERROR(message)
1508  end if
1509  multibinit_dtset%bound_step=1000
1510  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bound_step',tread,'INT')
1511  if(tread==1) multibinit_dtset%bound_step=intarr(1)
1512  if(multibinit_dtset%bound_step<0.and.multibinit_dtset%bound_step>1)then
1513    write(message, '(a,i8,a,a,a,a,a)' )&
1514 &   'bound_step is',multibinit_dtset%bound_step,', but the only allowed values',ch10,&
1515 &   'are 0 or 1 for multibinit.',ch10,&
1516 &   'Action: correct bound_step in your input file.'
1517    MSG_ERROR(message)
1518  end if
1519 
1520  multibinit_dtset%fit_coeff=0
1521  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_coeff',tread,'INT')
1522  if(tread==1) multibinit_dtset%fit_coeff=intarr(1)
1523  if(multibinit_dtset%fit_coeff<0.and.multibinit_dtset%fit_coeff>1)then
1524    write(message, '(a,i8,a,a,a,a,a)' )&
1525 &   'fit_coeff is',multibinit_dtset%fit_coeff,', but the only allowed values',ch10,&
1526 &   'are 0 or 1 for multibinit.',ch10,&
1527 &   'Action: correct fit_coeff in your input file.'
1528    MSG_ERROR(message)
1529  end if
1530 
1531  multibinit_dtset%fit_cutoff=0
1532  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_cutoff',tread,'DPR')
1533  if(tread==1) multibinit_dtset%fit_cutoff=dprarr(1)
1534  if(multibinit_dtset%fit_cutoff<0)then
1535    write(message, '(a,i8,a,a,a,a,a)' )&
1536 &   'fit_cutoff is',multibinit_dtset%fit_cutoff,', but the only allowed values',ch10,&
1537 &   'are positives for multibinit.',ch10,&
1538 &   'Action: correct fit_cutoff in your input file.'
1539    MSG_ERROR(message)
1540  end if
1541 
1542  multibinit_dtset%fit_grid(:)= 1
1543  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'fit_grid',tread,'INT')
1544  if(tread==1) multibinit_dtset%fit_grid(1:3)=intarr(1:3)
1545  do ii=1,3
1546    if(multibinit_dtset%fit_grid(ii)<0.or.multibinit_dtset%fit_grid(ii)>20)then
1547      write(message, '(a,i0,a,i0,a,a,a,i0,a)' )&
1548 &     'fit_grid(',ii,') is ',multibinit_dtset%fit_grid(ii),', which is lower',&
1549 &     ' than 0 of superior than 20.',&
1550 &     ch10,'Action: correct fit_grid(',ii,') in your input file.'
1551      MSG_ERROR(message)
1552    end if
1553  end do
1554 
1555  multibinit_dtset%fit_rangePower(:)= (/3,4/)
1556  call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'fit_rangePower',tread,'INT')
1557  if(tread==1) multibinit_dtset%fit_rangePower(1:2)=intarr(1:2)
1558  do ii=1,2
1559    if(multibinit_dtset%fit_rangePower(ii)<0.or.multibinit_dtset%fit_rangePower(ii)>20)then
1560      write(message, '(a,i0,a,i0,a,a,a,i0,a)' )&
1561 &     'fit_rangePower(',ii,') is ',multibinit_dtset%fit_rangePower(ii),', which is lower',&
1562 &     ' than 0 of superior than 20.',&
1563 &     ch10,'Action: correct fit_rangePower(',ii,') in your input file.'
1564      MSG_ERROR(message)
1565    end if
1566  end do
1567 
1568  multibinit_dtset%bound_rangePower(:)= (/6,6/)
1569  call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'bound_rangePower',tread,'INT')
1570  if(tread==1) multibinit_dtset%bound_rangePower(1:2)=intarr(1:2)
1571  do ii=1,2
1572    if(multibinit_dtset%bound_rangePower(ii)<=0.or.multibinit_dtset%bound_rangePower(ii)>20)then
1573      write(message, '(a,i0,a,i0,a,a,a,i0,a)' )&
1574 &     'bound_rangePower(',ii,') is ',multibinit_dtset%bound_rangePower(ii),', which is lower',&
1575 &     ' than 0 of superior than 20.',&
1576 &     ch10,'Action: correct bound_rangePower(',ii,') in your input file.'
1577      MSG_ERROR(message)
1578    end if
1579  end do
1580 
1581   multibinit_dtset%bound_cell(:)= (/6,6,6/)
1582  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'bound_cell',tread,'INT')
1583  if(tread==1) multibinit_dtset%bound_cell(1:3)=intarr(1:3)
1584  do ii=1,3
1585    if(multibinit_dtset%bound_cell(ii)<=0.or.multibinit_dtset%bound_cell(ii)>20)then
1586      write(message, '(a,i0,a,i0,4a,i0,a)' )&
1587 &     'bound_cell(',ii,') is ',multibinit_dtset%bound_cell(ii),', which is lower',&
1588 &     ' than 0 of superior than 20.',&
1589 &     ch10,'Action: correct bound_cell(',ii,') in your input file.'
1590      MSG_ERROR(message)
1591    end if
1592  end do
1593 
1594  multibinit_dtset%fit_tolMSDE=0
1595  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_tolMSDE',tread,'DPR')
1596  if(tread==1) multibinit_dtset%fit_tolMSDE=dprarr(1)
1597  if(multibinit_dtset%fit_tolMSDE<0)then
1598    write(message, '(a,i8,a,a,a,a,a)' )&
1599 &   'fit_tolMSDE is',multibinit_dtset%fit_tolMSDE,', but the only allowed values',ch10,&
1600 &   'are positives for multibinit.',ch10,&
1601 &   'Action: correct fit_tolMSDE in your input file.'
1602    MSG_ERROR(message)
1603  end if
1604 
1605  multibinit_dtset%fit_tolMSDF=0
1606  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_tolMSDF',tread,'DPR')
1607  if(tread==1) multibinit_dtset%fit_tolMSDF=dprarr(1)
1608  if(multibinit_dtset%fit_tolMSDF<0)then
1609    write(message, '(a,i8,a,a,a,a,a)' )&
1610 &   'fit_tolMSDF is',multibinit_dtset%fit_tolMSDF,', but the only allowed values',ch10,&
1611 &   'are positives for multibinit.',ch10,&
1612 &   'Action: correct fit_tolMSDF in your input file.'
1613    MSG_ERROR(message)
1614  end if
1615 
1616  multibinit_dtset%fit_tolMSDS=0
1617  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_tolMSDS',tread,'DPR')
1618  if(tread==1) multibinit_dtset%fit_tolMSDS=dprarr(1)
1619  if(multibinit_dtset%fit_tolMSDS<0)then
1620    write(message, '(a,i8,a,a,a,a,a)' )&
1621 &   'fit_tolMSDS is',multibinit_dtset%fit_tolMSDS,', but the only allowed values',ch10,&
1622 &   'are positives for multibinit.',ch10,&
1623 &   'Action: correct fit_tolMSDS in your input file.'
1624    MSG_ERROR(message)
1625  end if
1626 
1627  multibinit_dtset%fit_tolMSDFS=0
1628  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fit_tolMSDFS',tread,'DPR')
1629  if(tread==1) multibinit_dtset%fit_tolMSDFS=dprarr(1)
1630  if(multibinit_dtset%fit_tolMSDFS<0)then
1631    write(message, '(a,i8,a,a,a,a,a)' )&
1632 &   'fit_tolMSDFS is',multibinit_dtset%fit_tolMSDFS,', but the only allowed values',ch10,&
1633 &   'are positives for multibinit.',ch10,&
1634 &   'Action: correct fit_tolMSDFS in your input file.'
1635    MSG_ERROR(message)
1636  end if
1637 
1638 !G
1639 
1640 !H
1641 
1642 !I
1643 
1644 
1645 !J
1646 
1647 !K
1648 
1649 !L
1650 
1651 !M
1652 !N
1653 
1654  ABI_ALLOCATE(multibinit_dtset%fit_bancoeff,(multibinit_dtset%fit_nbancoeff))
1655  if (multibinit_dtset%fit_nbancoeff >0)then
1656    if(multibinit_dtset%fit_nbancoeff>marr)then
1657      marr=multibinit_dtset%fit_nbancoeff
1658      ABI_DEALLOCATE(intarr)
1659      ABI_ALLOCATE(intarr,(marr))
1660    end if
1661    multibinit_dtset%fit_bancoeff(:)=0
1662    call intagm(dprarr,intarr,jdtset,marr,multibinit_dtset%fit_nbancoeff,&
1663 &              string(1:lenstr),'fit_bancoeff',tread,'INT')
1664    if(tread==1)then
1665      do ii=1,multibinit_dtset%fit_nbancoeff
1666        multibinit_dtset%fit_bancoeff(ii)=intarr(ii)
1667      end do
1668    end if
1669  end if
1670 
1671  ABI_ALLOCATE(multibinit_dtset%fit_fixcoeff,(multibinit_dtset%fit_nfixcoeff))
1672  if (multibinit_dtset%fit_nfixcoeff >0)then
1673    if(multibinit_dtset%fit_nfixcoeff>marr)then
1674      marr=multibinit_dtset%fit_nfixcoeff
1675      ABI_DEALLOCATE(intarr)
1676      ABI_ALLOCATE(intarr,(marr))
1677    end if
1678    multibinit_dtset%fit_fixcoeff(:)=0
1679    call intagm(dprarr,intarr,jdtset,marr,multibinit_dtset%fit_nfixcoeff,&
1680 &              string(1:lenstr),'fit_fixcoeff',tread,'INT')
1681    if(tread==1)then
1682      do ii=1,multibinit_dtset%fit_nfixcoeff
1683        multibinit_dtset%fit_fixcoeff(ii)=intarr(ii)
1684      end do
1685    end if
1686  end if
1687 
1688 !O
1689 
1690 !P
1691 
1692 !Q
1693  ABI_ALLOCATE(multibinit_dtset%qmass,(multibinit_dtset%nnos))
1694  multibinit_dtset%qmass(:)= zero
1695  if(multibinit_dtset%nnos>=1)then
1696    if(multibinit_dtset%nnos>marr)then
1697      marr=multibinit_dtset%nnos
1698      ABI_DEALLOCATE(intarr)
1699      ABI_DEALLOCATE(dprarr)
1700      ABI_ALLOCATE(intarr,(marr))
1701      ABI_ALLOCATE(dprarr,(marr))
1702    end if
1703    call intagm(dprarr,intarr,jdtset,marr,multibinit_dtset%nnos,string(1:lenstr),'qmass',tread,'DPR')
1704    if(tread==1) multibinit_dtset%qmass(:)=dprarr(1:multibinit_dtset%nnos)
1705  end if
1706 
1707  if (multibinit_dtset%nqshft/=0)then
1708    if(3*multibinit_dtset%nqshft>marr)then
1709      marr=3*multibinit_dtset%nqshft
1710      ABI_DEALLOCATE(intarr)
1711      ABI_DEALLOCATE(dprarr)
1712      ABI_ALLOCATE(intarr,(marr))
1713      ABI_ALLOCATE(dprarr,(marr))
1714    end if
1715    ABI_ALLOCATE(multibinit_dtset%q1shft,(3,multibinit_dtset%nqshft))
1716    multibinit_dtset%q1shft(:,:)=zero
1717    call intagm(dprarr,intarr,jdtset,marr,3*multibinit_dtset%nqshft, string(1:lenstr),'q1shft',tread,'DPR')
1718    if(tread==1) multibinit_dtset%q1shft(1:3,1:multibinit_dtset%nqshft)=&
1719 &   reshape(dprarr(1:3*multibinit_dtset%nqshft),(/3,multibinit_dtset%nqshft/))
1720  end if
1721 
1722  ABI_ALLOCATE(multibinit_dtset%qph1l,(3,multibinit_dtset%nph1l))
1723  ABI_ALLOCATE(multibinit_dtset%qnrml1,(multibinit_dtset%nph1l))
1724  if (multibinit_dtset%nph1l/=0)then
1725    if(4*multibinit_dtset%nph1l>marr)then
1726      marr=4*multibinit_dtset%nph1l
1727      ABI_DEALLOCATE(intarr)
1728      ABI_DEALLOCATE(dprarr)
1729      ABI_ALLOCATE(intarr,(marr))
1730      ABI_ALLOCATE(dprarr,(marr))
1731    end if
1732    multibinit_dtset%qph1l(:,:)=zero
1733    multibinit_dtset%qnrml1(:)=zero
1734    call intagm(dprarr,intarr,jdtset,marr,4*multibinit_dtset%nph1l,string(1:lenstr),'qph1l',tread,'DPR')
1735    if(tread==1)then
1736      do iph1=1,multibinit_dtset%nph1l
1737        do ii=1,3
1738          multibinit_dtset%qph1l(ii,iph1)=dprarr(ii+(iph1-1)*4)
1739        end do
1740        multibinit_dtset%qnrml1(iph1)=dprarr(4+(iph1-1)*4)
1741        if(abs(multibinit_dtset%qnrml1(iph1))<DDB_QTOL)then
1742          write(message, '(a,a,a,a,a)' )&
1743 &         'The first list of wavevectors ','should not have non-analytical data.',ch10,&
1744 &         'Action: correct the first list',' of wavevectors in the input file.'
1745          MSG_ERROR(message)
1746        end if
1747      end do
1748    end if
1749  end if
1750 
1751  ABI_ALLOCATE(multibinit_dtset%qph2l,(3,multibinit_dtset%nph2l))
1752  ABI_ALLOCATE(multibinit_dtset%qnrml2,(multibinit_dtset%nph2l))
1753  if (multibinit_dtset%nph2l/=0)then
1754    if(4*multibinit_dtset%nph2l>marr)then
1755      marr=4*multibinit_dtset%nph2l
1756      ABI_DEALLOCATE(intarr)
1757      ABI_DEALLOCATE(dprarr)
1758      ABI_ALLOCATE(intarr,(marr))
1759      ABI_ALLOCATE(dprarr,(marr))
1760    end if
1761    multibinit_dtset%qph2l(:,:)=zero
1762    multibinit_dtset%qnrml2(:)=zero
1763    call intagm(dprarr,intarr,jdtset,marr,4*multibinit_dtset%nph2l,string(1:lenstr),'qph2l',tread,'DPR')
1764    if(tread==1)then
1765      do iph2=1,multibinit_dtset%nph2l
1766        do ii=1,3
1767          multibinit_dtset%qph2l(ii,iph2)=dprarr(ii+(iph2-1)*4)
1768        end do
1769        multibinit_dtset%qnrml2(iph2)=dprarr(4+(iph2-1)*4)
1770        if(abs(multibinit_dtset%qnrml2(iph2))>DDB_QTOL)then
1771          write(message, '(a,a,a,a,a)' )&
1772 &         'The second list of wavevectors',' should have only non-analytical data.',ch10,&
1773 &         'Action: correct the second list','of wavevectors in the input file.'
1774          MSG_ERROR(message)
1775        end if
1776      end do
1777    end if
1778  end if
1779 
1780 !R
1781  if(9>marr)then
1782    marr=9
1783    ABI_DEALLOCATE(intarr)
1784    ABI_DEALLOCATE(dprarr)
1785    ABI_ALLOCATE(intarr,(marr))
1786    ABI_ALLOCATE(dprarr,(marr))
1787  end if
1788  multibinit_dtset%rprim(:,:)= zero
1789  call intagm(dprarr,intarr,jdtset,marr,9,string(1:lenstr),'rprim',tread,'DPR')
1790  if(tread==1) then
1791    multibinit_dtset%rprim(1:3,1:3)= reshape(dprarr(1:9),(/3,3/))
1792 ! check new rprimd
1793    if(all(abs(multibinit_dtset%rprim(1,:))<tol16).or.&
1794 &    all(abs(multibinit_dtset%rprim(2,:))<tol16).or.all(abs(multibinit_dtset%rprim(3,:))<tol16)) then
1795      write(message, '(3a)' )&
1796 &  ' There is a problem with rprim',ch10,&
1797 &   'Action: correct rprim'
1798      MSG_BUG(message)
1799    end if
1800  end if
1801 !S
1802 
1803  if(6>marr)then
1804    marr=6
1805    ABI_DEALLOCATE(intarr)
1806    ABI_DEALLOCATE(dprarr)
1807    ABI_ALLOCATE(intarr,(marr))
1808    ABI_ALLOCATE(dprarr,(marr))
1809  end if
1810  multibinit_dtset%strten_reference(:)= zero
1811  call intagm(dprarr,intarr,jdtset,marr,6,string(1:lenstr),'strten_reference',tread,'DPR')
1812  if(tread==1) multibinit_dtset%strten_reference(1:6)= dprarr(1:6)
1813 
1814 !T
1815 
1816 !U
1817 
1818 !V
1819 
1820 !W
1821 
1822 !X
1823 
1824 !Y
1825 
1826 !Z
1827 
1828 !=======================================================================
1829 !Finished reading in variables - deallocate
1830 !=======================================================================
1831 
1832  ABI_DEALLOCATE(dprarr)
1833  ABI_DEALLOCATE(intarr)
1834 
1835 !=======================================================================
1836 !Check consistency of input variables:
1837 !=======================================================================
1838 
1839  if(multibinit_dtset%prtsrlr/=0 .and. multibinit_dtset%ifcflag/=1) then
1840    write(message, '(3a)' )&
1841 &   'ifcflag must be 1 for the SR/LR decomposition of the phonon frequencies',ch10,&
1842 &   'Action: correct ifcflag in your input file.'
1843    MSG_ERROR(message)
1844  end if
1845 
1846 !FIXME: add check that if freeze_displ /= 0 then you need to be doing ifc and phonon interpolation
1847 
1848  if (multibinit_dtset%ifcflag > 0 .and. sum(abs(multibinit_dtset%ngqpt)) == 0) then
1849    write(message, '(3a)' )&
1850 &   'if you want interatomic force constant output, multibinit needs ngqpt input variable ',ch10,&
1851 &   'Action: set ngqpt in your input file.'
1852    MSG_ERROR(message)
1853  end if
1854 
1855 !check that q-grid refinement is a divisor of ngqpt in each direction
1856  if(any(multibinit_dtset%qrefine(:) > 1) .and. &
1857 &    any(abs(dmod(dble(multibinit_dtset%ngqpt(1:3))/dble(multibinit_dtset%qrefine(:)),one)) > tol10)) then
1858    write(message, '(a,3i0,a,a,a,3i8,a,a)' )&
1859 &   'qrefine is',multibinit_dtset%qrefine,' The only allowed values',ch10,&
1860 &   'are integers which are divisors of the ngqpt grid', multibinit_dtset%ngqpt,ch10,&
1861 &   'Action: correct qrefine in your input file.'
1862    MSG_ERROR(message)
1863  end if
1864 
1865 ! check new rprimd
1866  if(all(multibinit_dtset%acell(:) > one).and.all(abs(multibinit_dtset%rprim(:,:))<tol16))then
1867    write(message, '(3a)' )&
1868 &         ' acell is defined but there is no rprim',ch10,&
1869 &         'Action: add rprim input'
1870    MSG_BUG(message)
1871  end if
1872 
1873 
1874 !check the fit_bancoeff and fit_fixcoeff
1875  do ii=1,multibinit_dtset%fit_nbancoeff
1876    do jj=ii+1,multibinit_dtset%fit_nbancoeff
1877      if (multibinit_dtset%fit_bancoeff(ii) == multibinit_dtset%fit_bancoeff(jj))then
1878        write(message, '(a,I0,a,I0,2a)' )&
1879 &           ' There is two similar numbers for fit_bancoeff: ',multibinit_dtset%fit_bancoeff(ii),&
1880 &           ' and ', multibinit_dtset%fit_bancoeff(jj),ch10,&
1881 &            'Action: change fit_bancoeff'
1882        MSG_BUG(message)
1883      end if
1884    end do
1885  end do
1886 
1887  do ii=1,multibinit_dtset%fit_nfixcoeff
1888    do jj=ii+1,multibinit_dtset%fit_nfixcoeff
1889      if (multibinit_dtset%fit_fixcoeff(ii) == multibinit_dtset%fit_fixcoeff(jj))then
1890        write(message, '(a,I0,a,I0,2a)' )&
1891 &           ' There is two similar numbers for fit_fixcoeff: ',multibinit_dtset%fit_fixcoeff(ii),&
1892 &           ' and ', multibinit_dtset%fit_fixcoeff(jj),ch10,&
1893 &            'Action: change fit_fixcoeff'
1894        MSG_BUG(message)
1895      end if
1896    end do
1897  end do
1898 
1899  do ii=1,3
1900    if(multibinit_dtset%dipdip_range(ii) < multibinit_dtset%ncell(ii)) then
1901      write(message,'(4a,3I3,3a,3I3,6a)') ch10,&
1902 &                 ' --- !WARNING',ch10,&
1903 &                 '     The range of dipdip_range (',multibinit_dtset%dipdip_range(:),')',ch10,&
1904 &                 '     But the range of the cell for the simulation is',&
1905 &                       multibinit_dtset%ncell(:),')',ch10,&
1906 &                 '     dipdip_range is set to ncell.',ch10,&
1907 &                 ' ---',ch10
1908      multibinit_dtset%dipdip_range(:) =  multibinit_dtset%ncell(:)
1909      call wrtout(std_out,message,'COLL')
1910      exit
1911    end if
1912    if(multibinit_dtset%dipdip_range(ii) < multibinit_dtset%bound_cell(ii)) then
1913      write(message,'(4a,3I3,3a,3I3,6a)') ch10,&
1914 &                 ' --- !WARNING',ch10,&
1915 &                 '     The range of dipdip_range (',multibinit_dtset%dipdip_range(:),')',ch10,&
1916 &                 '     But the range of the cell for the simulation is',&
1917 &                       multibinit_dtset%ncell(:),')',ch10,&
1918 &                 '     dipdip_range is set to bound_cell.',ch10,&
1919 &                 ' ---',ch10
1920      multibinit_dtset%dipdip_range(:) =  multibinit_dtset%ncell(:)
1921      call wrtout(std_out,message,'COLL')
1922      exit
1923    end if
1924  end do
1925 
1926 !Check if only one tolerance is specify
1927  if(abs(multibinit_dtset%fit_tolMSDF) >zero .and. abs(multibinit_dtset%fit_tolMSDS) >zero) then
1928    write(message, '(3a)' ) &
1929 &           ' There is two tolerance flags for the fit: fit_tolMSDF and fit_tolMSDS',ch10,&
1930 &            'Action: Put only one tolerance flag'
1931    MSG_BUG(message)
1932  end if
1933  if(abs(multibinit_dtset%fit_tolMSDF) >zero .and. abs(multibinit_dtset%fit_tolMSDE) >zero)then
1934    write(message, '(3a)' ) &
1935 &           ' There is two tolerance flags for the fit: fit_tolMSDF and fit_tolMSDE',ch10,&
1936 &            'Action: Put only one tolerance flag'
1937    MSG_BUG(message)
1938  end if
1939  if(abs(multibinit_dtset%fit_tolMSDF) >zero .and. abs(multibinit_dtset%fit_tolMSDFS) >zero)then
1940    write(message, '(3a)' ) &
1941 &           ' There is two tolerance flags for the fit: fit_tolMSDF and fit_tolMSDFS',ch10,&
1942 &            'Action: Put only one tolerance flag'
1943    MSG_BUG(message)
1944  end if
1945  if(abs(multibinit_dtset%fit_tolMSDS) >zero .and. abs(multibinit_dtset%fit_tolMSDE) >zero)then
1946    write(message, '(3a)' ) &
1947 &           ' There is two tolerance flags for the fit: fit_tolMSDS and fit_tolMSDE',ch10,&
1948 &            'Action: Put only one tolerance flag'
1949    MSG_BUG(message)
1950  end if
1951  if(abs(multibinit_dtset%fit_tolMSDS) >zero .and. abs(multibinit_dtset%fit_tolMSDFS) >zero)then
1952    write(message, '(3a)' ) &
1953 &           ' There is two tolerance flags for the fit: fit_tolMSDS and fit_tolMSDFS',ch10,&
1954 &            'Action: Put only one tolerance flag'
1955    MSG_BUG(message)
1956  end if
1957  if(abs(multibinit_dtset%fit_tolMSDE) >zero .and. abs(multibinit_dtset%fit_tolMSDFS) >zero)then
1958    write(message, '(3a)' ) &
1959 &           ' There is two tolerance flags for the fit: fit_tolMSDE and fit_tolMSDFS',ch10,&
1960 &            'Action: Put only one tolerance flag'
1961    MSG_BUG(message)
1962  end if
1963 
1964 end subroutine invars10

m_multibinit_dataset/multibinit_dtset_free [ Functions ]

[ Top ] [ m_multibinit_dataset ] [ Functions ]

NAME

  multibinit_dtset_free

FUNCTION

  deallocate remaining arrays in the multibinit_dtset datastructure

INPUTS

  multibinit_dtset <type(multibinit_dtset_type)> = multibinit_dataset structure

 OUTPUTS
  multibinit_dtset <type(multibinit_dtset_type)> = multibinit_dataset structure

PARENTS

      multibinit

CHILDREN

SOURCE

398 subroutine multibinit_dtset_free(multibinit_dtset)
399 
400 
401 !This section has been created automatically by the script Abilint (TD).
402 !Do not modify the following lines by hand.
403 #undef ABI_FUNC
404 #define ABI_FUNC 'multibinit_dtset_free'
405 !End of the abilint section
406 
407  implicit none
408 
409 !Arguments ------------------------------------
410 !scalars
411  type(multibinit_dtset_type), intent(inout) :: multibinit_dtset
412 
413 ! *************************************************************************
414 
415  if (allocated(multibinit_dtset%atifc))  then
416    ABI_DEALLOCATE(multibinit_dtset%atifc)
417  end if
418  if (allocated(multibinit_dtset%conf_cutoff_disp))  then
419    ABI_DEALLOCATE(multibinit_dtset%conf_cutoff_disp)
420  end if
421  if (allocated(multibinit_dtset%fit_fixcoeff))  then
422    ABI_DEALLOCATE(multibinit_dtset%fit_fixcoeff)
423  end if
424   if (allocated(multibinit_dtset%fit_bancoeff))  then
425    ABI_DEALLOCATE(multibinit_dtset%fit_bancoeff)
426  end if
427  if (allocated(multibinit_dtset%qmass))  then
428    ABI_DEALLOCATE(multibinit_dtset%qmass)
429  end if
430  if (allocated(multibinit_dtset%coefficients))  then
431    ABI_DEALLOCATE(multibinit_dtset%coefficients)
432  end if
433  if (allocated(multibinit_dtset%qnrml1))  then
434    ABI_DEALLOCATE(multibinit_dtset%qnrml1)
435  end if
436  if (allocated(multibinit_dtset%qnrml2))  then
437    ABI_DEALLOCATE(multibinit_dtset%qnrml2)
438  end if
439  if (allocated(multibinit_dtset%qph1l))  then
440    ABI_DEALLOCATE(multibinit_dtset%qph1l)
441  end if
442  if (allocated(multibinit_dtset%qph2l))  then
443    ABI_DEALLOCATE(multibinit_dtset%qph2l)
444  end if
445  if(allocated(multibinit_dtset%q1shft))then
446    ABI_DEALLOCATE(multibinit_dtset%q1shft)
447  end if
448 
449  !if (allocated(multibinit_dtset%gilbert_damping))  then
450  !  ABI_DEALLOCATE(multibinit_dtset%gilbert_damping)
451  !end if
452 
453  !if (allocated(multibinit_dtset%gyro_ratio))  then
454  !  ABI_DEALLOCATE(multibinit_dtset%gyro_ratio)
455  !end if
456 
457  !if (allocated(multibinit_dtset%qph1l_spin))  then
458  !  ABI_DEALLOCATE(multibinit_dtset%qph1l_spin)
459  !end if
460  !if (allocated(multibinit_dtset%qph2l_spin))  then
461  !  ABI_DEALLOCATE(multibinit_dtset%qph2l_spin)
462  !end if
463 
464 
465 end subroutine multibinit_dtset_free

m_multibinit_dataset/multibinit_dtset_init [ Functions ]

[ Top ] [ m_multibinit_dataset ] [ Functions ]

NAME

 multibinit_dtset_init

FUNCTION

 Init the dtset datatype

INPUTS

 natom=number of atoms, needed for atifc

OUTPUT

 multibinit_dtset <type(multibinit_dtset_type)> = datatype with all the input variables

NOTES

 Should be executed by one processor only.

PARENTS

CHILDREN

SOURCE

240 subroutine multibinit_dtset_init(multibinit_dtset,natom)
241 
242 
243 !This section has been created automatically by the script Abilint (TD).
244 !Do not modify the following lines by hand.
245 #undef ABI_FUNC
246 #define ABI_FUNC 'multibinit_dtset_init'
247 !End of the abilint section
248 
249  implicit none
250 
251 !Arguments -------------------------------
252 !scalars
253  integer,intent(in) :: natom
254  type(multibinit_dtset_type),intent(inout) :: multibinit_dtset
255 !Local variables -------------------------
256 !scalars
257 !arrays
258 
259 !*********************************************************************
260 
261 !copy natom to multibinit_dtset
262  multibinit_dtset%natom=natom
263 
264 !=====================================================================
265 !Scalars
266 !=====================================================================
267  multibinit_dtset%asr=2
268  multibinit_dtset%brav=1
269  multibinit_dtset%bmass=0
270  multibinit_dtset%chneut=0
271  multibinit_dtset%confinement=0
272  multibinit_dtset%conf_power_disp=0
273  multibinit_dtset%conf_power_strain=0
274  multibinit_dtset%conf_power_fact_disp=100
275  multibinit_dtset%conf_power_fact_strain=100
276  multibinit_dtset%delta_df= 1d-02
277  multibinit_dtset%dipdip=1
278  multibinit_dtset%dipdip_prt=0
279  multibinit_dtset%dtion=100
280  multibinit_dtset%dynamics=0
281  multibinit_dtset%eivec=0
282  multibinit_dtset%energy_reference= zero
283  multibinit_dtset%enunit=0
284  multibinit_dtset%fit_anhaStrain=0
285  multibinit_dtset%bound_model=0
286  multibinit_dtset%bound_anhaStrain=0
287  multibinit_dtset%bound_cutoff=0 
288  multibinit_dtset%bound_maxCoeff=4
289  multibinit_dtset%bound_temp=325
290  multibinit_dtset%bound_step=1000
291  multibinit_dtset%bound_SPCoupling=1
292  multibinit_dtset%fit_coeff=0
293  multibinit_dtset%fit_cutoff=0
294  multibinit_dtset%fit_nbancoeff=0
295  multibinit_dtset%fit_ncoeff=0
296  multibinit_dtset%ts_option=0
297  multibinit_dtset%fit_nfixcoeff=0
298  multibinit_dtset%fit_option=0
299  multibinit_dtset%fit_SPCoupling=1
300  multibinit_dtset%fit_generateCoeff=1
301  multibinit_dtset%fit_initializeData=1
302  multibinit_dtset%fit_tolMSDE=zero
303  multibinit_dtset%fit_tolMSDS=zero
304  multibinit_dtset%fit_tolMSDF=zero
305  multibinit_dtset%fit_tolMSDFS=zero
306  multibinit_dtset%ifcana=0
307  multibinit_dtset%ifcflag=1
308  multibinit_dtset%ifcout=-1
309  multibinit_dtset%prtsrlr=0
310  multibinit_dtset%ntime=200
311  multibinit_dtset%nctime=1
312  multibinit_dtset%natifc=natom
313  multibinit_dtset%ncoeff=0
314  multibinit_dtset%nph1l=1
315  multibinit_dtset%nph2l=0
316  multibinit_dtset%nqshft=1
317  multibinit_dtset%nnos=0
318  multibinit_dtset%nsphere=0
319  multibinit_dtset%optcell=0
320  multibinit_dtset%prt_model=0
321  multibinit_dtset%prt_phfrq=0
322  multibinit_dtset%prt_ifc = 0
323  multibinit_dtset%strcpling = -1
324  multibinit_dtset%qrefine=1
325  multibinit_dtset%restartxf=0
326  multibinit_dtset%rfmeth=1
327  multibinit_dtset%rifcsph=zero
328  multibinit_dtset%symdynmat=1
329  multibinit_dtset%temperature=325
330  
331  multibinit_dtset%spin_dipdip=0
332  multibinit_dtset%spin_dynamics=1
333  multibinit_dtset%spin_ntime=10000
334  multibinit_dtset%spin_nctime=100
335  multibinit_dtset%spin_nmatom=0
336  multibinit_dtset%spin_n1l=1
337  multibinit_dtset%spin_n2l=0
338 
339  multibinit_dtset%spin_dt=1d-16
340  !TODO hexu: what is the unit. s. here.  the atomic unit is 2.418884e-17, should we use it ?
341 multibinit_dtset%spin_temperature=325 ! or 0 ?
342 multibinit_dtset%spin_tolavg=1d-2 ! TODO hexu: to be decided. should it be a function of temperature?
343 multibinit_dtset%spin_tolvar=1d-3 ! TODO hexu: as above. 
344 
345 
346 !=======================================================================
347 !Arrays
348 !=======================================================================
349  multibinit_dtset%acell(:) = one
350  multibinit_dtset%conf_cutoff_strain(1:6) = zero
351  multibinit_dtset%dipdip_range(:)= (/0,0,0/)
352  multibinit_dtset%fit_grid(:)= 1
353  multibinit_dtset%fit_rangePower(:)= (/3,4/)
354  multibinit_dtset%bound_rangePower(:)= (/6,6/)
355  multibinit_dtset%bound_cell(:)= (/6,6,6/)
356  multibinit_dtset%ncell(:)= 0
357  multibinit_dtset%ngqpt(:) = 0
358  multibinit_dtset%ng2qpt(:)= 0
359  multibinit_dtset%strtarget(1:6) = zero
360  multibinit_dtset%qmass(:)= zero
361  multibinit_dtset%rprim(:,:)= zero
362  multibinit_dtset%strten_reference(:)= zero
363 
364  multibinit_dtset%spin_mag_field(:)=zero
365  multibinit_dtset%spin_qpoint(:)=zero
366 
367  ABI_ALLOCATE(multibinit_dtset%atifc,(natom))
368  multibinit_dtset%atifc(:)=0
369  ABI_ALLOCATE(multibinit_dtset%conf_cutoff_disp,(multibinit_dtset%natom))
370  multibinit_dtset%conf_cutoff_disp(:)=zero
371  ABI_ALLOCATE(multibinit_dtset%q1shft,(3,multibinit_dtset%nqshft))
372  multibinit_dtset%q1shft(:,:) = zero
373  
374 end subroutine multibinit_dtset_init

m_multibinit_dataset/multibinit_dtset_type [ Types ]

[ Top ] [ m_multibinit_dataset ] [ Types ]

NAME

 multibinit_dtset_type

FUNCTION

 The multibinit_dtset_type structured datatype
 gather all the input variables for the multibinit code.

SOURCE

 59  type multibinit_dtset_type
 60 
 61 ! Integer
 62   integer :: asr
 63   integer :: brav
 64   integer :: chneut
 65   integer :: confinement
 66   integer :: conf_power_disp
 67   integer :: conf_power_strain
 68   integer :: dipdip
 69   integer :: eivec
 70   integer :: elphflag
 71   integer :: enunit
 72   integer :: bound_model
 73   integer :: bound_maxCoeff
 74   integer :: bound_SPCoupling
 75   integer :: bound_AnhaStrain
 76   integer :: bound_step
 77   integer :: fit_anhaStrain
 78   integer :: fit_SPCoupling
 79   integer :: fit_generateCoeff
 80   integer :: fit_initializeData
 81   integer :: fit_coeff
 82   integer :: fit_option
 83   integer :: fit_ncoeff
 84   integer :: fit_nbancoeff
 85   integer :: fit_nfixcoeff
 86   integer :: ts_option
 87   integer :: ifcana
 88   integer :: ifcflag
 89   integer :: ifcout
 90   ! TODO hexu: why integer dtion?
 91   integer :: dtion
 92   integer :: dynamics
 93   ! TODO: use dynamics for spin or add a keyword?
 94   integer :: natifc
 95   integer :: natom
 96   integer :: ncoeff
 97   integer :: nctime
 98   integer :: ntime
 99   integer :: nnos
100   integer :: nph1l
101   integer :: nph2l
102   integer :: nqshft
103   integer :: nsphere
104   integer :: optcell
105   integer :: prt_model
106   integer :: dipdip_prt
107   integer :: prt_phfrq
108   integer :: prt_ifc
109   integer :: strcpling  ! Print the 3rd order in xml file
110   integer :: prtsrlr  ! print the short-range/long-range decomposition of phonon freq.
111   integer :: rfmeth
112   integer :: restartxf
113   integer :: symdynmat
114   integer :: dipdip_range(3)
115   integer :: fit_grid(3)
116   integer :: fit_rangePower(2)
117   integer :: bound_rangePower(2)
118   integer :: bound_cell(3)
119   integer :: ncell(3)
120   integer :: ngqpt(9)             ! ngqpt(9) instead of ngqpt(3) is needed in wght9.f
121   integer :: ng2qpt(3)
122   integer :: kptrlatt(3,3)
123   integer :: kptrlatt_fine(3,3)
124   integer :: qrefine(3)
125 
126 
127   ! TODO hexu: add parameters for spin.
128   integer :: spin_dipdip
129   integer :: spin_dynamics
130   integer :: spin_nctime
131   integer :: spin_ntime
132   integer :: spin_nmatom !TODO hexu: is it needed?
133   integer :: spin_n1l
134   integer :: spin_n2l
135 
136 ! Real(dp)
137   real(dp) :: bmass
138   real(dp) :: conf_power_fact_disp
139   real(dp) :: conf_power_fact_strain
140   real(dp) :: delta_df
141   real(dp) :: energy_reference
142   real(dp) :: bound_cutoff
143   real(dp) :: bound_Temp
144   real(dp) :: fit_cutoff
145   real(dp) :: fit_tolMSDF
146   real(dp) :: fit_tolMSDS
147   real(dp) :: fit_tolMSDE
148   real(dp) :: fit_tolMSDFS
149   real(dp) :: temperature
150   real(dp) :: rifcsph
151   real(dp) :: conf
152   real(dp) :: acell(3)
153   real(dp) :: strten_reference(6)
154   real(dp) :: strtarget(6)
155   real(dp) :: conf_cutoff_strain(6)
156   real(dp) :: rprim(3,3)
157 
158   ! TODO hexu:add parameters for spin
159   real(dp) :: spin_dt
160   real(dp) :: spin_temperature ! TODO hexu: should we consider T(spin)/=T(latt)
161   ! TODO hexu: add spin convergence tol. (or remove it)
162   real(dp) :: spin_tolavg !average
163   real(dp) :: spin_tolvar !covariance
164 
165   real(dp) :: spin_mag_field(3)  ! external magnetic field
166   real(dp) :: spin_qpoint(3) 
167 ! Integer arrays
168   integer, allocatable :: atifc(:)
169   ! atifc(natom)
170   integer, allocatable :: fit_fixcoeff(:)
171   ! fit_fixcoeffs(fit_nfixcoeff)
172 
173   integer, allocatable :: fit_bancoeff(:)
174   ! fit_bancoeffs(fit_nbancoeff)
175 
176   !integer, allocatable :: spin_sublattice(:) ! TODO hexu: difficult to use, better in xml?
177 
178   real(dp), allocatable :: qmass(:)
179   ! qmass(nnos)
180 
181 
182 ! Real arrays
183   real(dp), allocatable :: coefficients(:)
184   ! coefficients(ncoeff)
185 
186   real(dp), allocatable :: conf_cutoff_disp(:)
187   ! conf_cuttoff(natom)
188   
189   real(dp),allocatable  :: q1shft(:,:)
190   !q1shft(3,nqshft)  SHIFT for Q point
191 
192   real(dp), allocatable :: qnrml1(:)
193   ! qnrml1(nph1l)
194 
195   real(dp), allocatable :: qnrml2(:)
196   ! qnrml1(nph1l)
197 
198   real(dp), allocatable :: qph1l(:,:)
199   ! qph1l(3,nph1l)
200 
201   real(dp), allocatable :: qph2l(:,:)
202   ! qph2l(3,nph2l)
203 
204   ! spin part
205   !real(dp), allocatable :: gilbert_damping ! if not provided in xml or override is needed. 
206   !real(dp), allocatable :: gyro_ratio(:) ! if not provided in xml
207 
208   !real(dp), allocatable :: qspin1l(:,:)
209   !real(dp), allocatable :: qspin2l(:,:)
210 
211  end type multibinit_dtset_type

m_multibinit_dataset/outvars_multibinit [ Functions ]

[ Top ] [ m_multibinit_dataset ] [ Functions ]

NAME

 outvars_multibinit

FUNCTION

 Open input file for the multibinit code, then
 echoes the input information.

INPUTS

 multibinit_dtset <type(multibinit_dtset_type)> datatype with all the input variables
 nunit=unit number for input or output

OUTPUT

  (only writing)

NOTES

 Should be executed by one processor only.

PARENTS

      multibinit

CHILDREN

SOURCE

1997 subroutine outvars_multibinit (multibinit_dtset,nunit)
1998 
1999 
2000 !This section has been created automatically by the script Abilint (TD).
2001 !Do not modify the following lines by hand.
2002 #undef ABI_FUNC
2003 #define ABI_FUNC 'outvars_multibinit'
2004 !End of the abilint section
2005 
2006  implicit none
2007 
2008 !Arguments -------------------------------
2009 !scalars
2010  integer,intent(in) :: nunit
2011  type(multibinit_dtset_type),intent(in) :: multibinit_dtset
2012 
2013 !Local variables -------------------------
2014 !Set routine version number here:
2015 !scalars
2016  integer :: ii,iph1,iph2,iqshft
2017 
2018 !*********************************************************************
2019 
2020 !Write the heading
2021  write(nunit,'(a,80a,a)') ch10,('=',ii=1,80),ch10
2022  write(nunit, '(a,a)' )&
2023 & ' -outvars_multibinit: echo values of input variables ----------------------',ch10
2024 
2025 !The flags
2026  if(multibinit_dtset%ifcflag/=0)then
2027    write(nunit,'(a)')' Flags : '
2028    if(multibinit_dtset%ifcflag/=0)write(nunit,'(3x,a9,3i10)')'  ifcflag',multibinit_dtset%ifcflag
2029    if(multibinit_dtset%prt_model/=0)write(nunit,'(3x,a9,3i10)')'prt_model',multibinit_dtset%prt_model
2030    if(multibinit_dtset%prt_phfrq/=0)write(nunit,'(3x,a9,3i10)')'prt_phfrq',multibinit_dtset%prt_phfrq
2031    if(multibinit_dtset%strcpling/=0)write(nunit,'(3x,a9,3i10)')'  strcpling',multibinit_dtset%strcpling
2032    if(multibinit_dtset%strcpling==2)write(nunit,'(3x,a9,3es8.2)')'delta_df',multibinit_dtset%delta_df
2033  end if
2034 
2035  if(multibinit_dtset%dynamics/=0)then
2036    write(nunit,'(a)')' Molecular Dynamics :'
2037    write(nunit,'(3x,a9,3I10.1)')' dynamics',multibinit_dtset%dynamics
2038    write(nunit,'(3x,a9,3F10.1)')'     temp',multibinit_dtset%temperature
2039    write(nunit,'(3x,a9,3I10.1)')'    ntime',multibinit_dtset%ntime
2040    if (multibinit_dtset%nctime /=1)then
2041      write(nunit,'(3x,a9,3I10.1)')'   nctime',multibinit_dtset%nctime
2042    end if
2043    write(nunit,'(3x,a9,3i10)')  '    ncell',multibinit_dtset%ncell
2044    write(nunit,'(3x,a9,3i10)')  '    dtion',multibinit_dtset%dtion
2045    if (multibinit_dtset%restartxf/=0) then
2046      write(nunit,'(3x,a9,3i10)')  'restartxf',multibinit_dtset%restartxf
2047    end if
2048    if(multibinit_dtset%dynamics==13)then
2049      write(nunit,'(3x,a9,3i10)')'  optcell',multibinit_dtset%optcell
2050      write(nunit,'(3x,a9,3F12.1)')'    bmass',multibinit_dtset%bmass
2051      write(nunit,'(3x,a9,3I10)')'     nnos',multibinit_dtset%nnos
2052      write(nunit,'(3x,a12)',advance='no')'    qmass  '
2053      write(nunit,'(3x,15F12.10)') (multibinit_dtset%qmass(ii),ii=1,multibinit_dtset%nnos)
2054    end if
2055  end if
2056 
2057  if(multibinit_dtset%spin_dynamics/=0) then
2058     write(nunit,'(a)')' Spin Dynamics :'
2059     write(nunit,'(2x,a16,3I12.1)')'spin_dynamics',multibinit_dtset%spin_dynamics
2060     write(nunit,'(a18,3es15.5)')'spin_temperature',multibinit_dtset%spin_temperature
2061     write(nunit,'(3x,a15,3I10.1)')'spin_ntime',multibinit_dtset%spin_ntime
2062     write(nunit,'(3x,a15,3I10)')  'ncell',multibinit_dtset%ncell !TODO hexu: duplicate but dynamics can be 0.
2063     write(nunit,'(3x,a15,3es15.5)')  'spin_dt',multibinit_dtset%spin_dt
2064     !write(nunit,'(3x,a14,3es10.5)')  '   spin_tolavg',multibinit_dtset%spin_tolavg
2065     !write(nunit,'(3x,a14,3es10.5)')  '   spin_tolvar',multibinit_dtset%spin_tolvar
2066     write(nunit,'(3x,a15)')   'spin_mag_field'
2067     write(nunit,'(19x,3es12.5)')   (multibinit_dtset%spin_mag_field(ii),ii=1,3)
2068     write(nunit,'(3x,a15)')   'spin_qpoint'
2069     write(nunit,'(19x,3es12.5)')   (multibinit_dtset%spin_qpoint(ii),ii=1,3)
2070 
2071  end if
2072 
2073  if(multibinit_dtset%confinement==1)then
2074    write(nunit,'(a)')' Confinement information :'
2075    write(nunit,'(1x,a22,I5.1)')'       conf_power_disp',multibinit_dtset%conf_power_disp
2076    write(nunit,'(1x,a22,I5.1)')'     conf_power_strain',multibinit_dtset%conf_power_strain
2077    write(nunit,'(1x,a22,3es16.8)')'  conf_power_fact_disp',multibinit_dtset%conf_power_fact_disp
2078    write(nunit,'(1x,a22,3es16.8)')'conf_power_fact_strain',multibinit_dtset%conf_power_fact_strain
2079    write(nunit,'(1x,a22)')'     conf_cutoff_disp'
2080    write(nunit,'(19x,3es16.8)') (multibinit_dtset%conf_cutoff_disp(ii),ii=1,multibinit_dtset%natom)
2081    write(nunit,'(1x,a22)')'    conf_cutoff_strain'
2082    write(nunit,'(19x,3es16.8)') (multibinit_dtset%conf_cutoff_strain(ii),ii=1,6)
2083  end if
2084 
2085  if(multibinit_dtset%fit_coeff/=0)then
2086    write(nunit,'(a)')' Fit the coefficients :'
2087    write(nunit,'(1x,a17,I3.1)')'       fit_coeff',multibinit_dtset%fit_coeff
2088    write(nunit,'(1x,a17,I3.1)')'fit_generateCoeff',multibinit_dtset%fit_generateCoeff
2089    if(multibinit_dtset%fit_initializeData==0)then
2090      write(nunit,'(1x,a17,I3.1)')'fit_initializeData',multibinit_dtset%fit_initializeData
2091    end if
2092    if(multibinit_dtset%fit_tolMSDE > 0)then
2093      write(nunit,'(1x,a17,es16.8)')'    fit_tolMSDE',multibinit_dtset%fit_tolMSDE
2094    end if
2095    if(multibinit_dtset%fit_tolMSDF > 0)then
2096      write(nunit,'(1x,a17,es16.8)')'    fit_tolMSDF',multibinit_dtset%fit_tolMSDF
2097    end if
2098    if(multibinit_dtset%fit_tolMSDS > 0)then
2099      write(nunit,'(1x,a17,es16.8)')'    fit_tolMSDS',multibinit_dtset%fit_tolMSDS
2100    end if
2101    if(multibinit_dtset%fit_tolMSDFS > 0)then
2102      write(nunit,'(1x,a17,es16.8)')'   fit_tolMSDFS',multibinit_dtset%fit_tolMSDFS
2103    end if
2104    write(nunit,'(1x,a17,es16.8)')'      fit_cutoff',multibinit_dtset%fit_cutoff
2105    write(nunit,'(1x,a17,I3.1)')'      fit_option',multibinit_dtset%fit_option
2106    write(nunit,'(1x,a17,2x,I0)')'      fit_ncoeff',multibinit_dtset%fit_ncoeff   
2107    write(nunit,'(1x,a17,3i3)') '        fit_grid',multibinit_dtset%fit_grid
2108    write(nunit,'(1x,a17,I3.1)')'   ts_option',multibinit_dtset%ts_option
2109    write(nunit,'(1x,a17,2i3)') '  fit_rangePower',multibinit_dtset%fit_rangePower
2110    write(nunit,'(1x,a17,I3)')  '  fit_anhaStrain',multibinit_dtset%fit_anhaStrain
2111    write(nunit,'(1x,a17,I3)')  '  fit_SPCoupling',multibinit_dtset%fit_SPCoupling
2112    if(multibinit_dtset%fit_nbancoeff /= 0) then
2113      write(nunit,'(1x,a17,I3)')  '   fit_nbancoeff',multibinit_dtset%fit_nbancoeff
2114      write(nunit,'(1x,a17)',advance='no')'   fit_bancoeff'
2115      write(nunit,'(4x,9i7)') (multibinit_dtset%fit_bancoeff(ii),ii=1,multibinit_dtset%fit_nbancoeff)
2116    end if
2117    if(multibinit_dtset%fit_nfixcoeff /= 0) then    
2118      write(nunit,'(1x,a17,I3)')  '   fit_nfixcoeff',multibinit_dtset%fit_nfixcoeff
2119      write(nunit,'(1x,a17)',advance='no')'   fit_fixcoeff'
2120      write(nunit,'(4x,9i7)') (multibinit_dtset%fit_fixcoeff(ii),ii=1,multibinit_dtset%fit_nfixcoeff)
2121    end if
2122  end if
2123 
2124  if(multibinit_dtset%bound_model /=0)then
2125    write(nunit,'(a)')' Bound the coefficients :'
2126    write(nunit,'(1x,a16,I3)')    'bound_anhaStrain',multibinit_dtset%bound_anhaStrain   
2127    write(nunit,'(1x,a16,I3)')    'bound_SPCoupling',multibinit_dtset%bound_SPCoupling
2128    write(nunit,'(1x,a16,es16.8)')'    bound_cutoff',multibinit_dtset%bound_cutoff
2129    write(nunit,'(1x,a16,1x,3I3)')   '      bound_cell',multibinit_dtset%bound_cell
2130    write(nunit,'(1x,a16,1x,I3)')    '  bound_maxCoeff',multibinit_dtset%bound_maxCoeff
2131    write(nunit,'(1x,a16,es16.8)') '      bound_temp',multibinit_dtset%bound_temp
2132    write(nunit,'(1x,a16,I7)')   '      bound_step',multibinit_dtset%bound_step
2133    write(nunit,'(1x,a16,2I3.1)')'bound_rangePower',multibinit_dtset%bound_rangePower
2134  end if
2135 
2136 !Write the general information
2137  if( multibinit_dtset%rfmeth/=1 .or. &
2138 & multibinit_dtset%enunit/=0 .or. &
2139 & multibinit_dtset%eivec/=0 .or. &
2140 & multibinit_dtset%asr/=0 .or. &
2141 & multibinit_dtset%chneut/=0)then
2142    write(nunit,'(a)')' Miscellaneous information :'
2143    if(multibinit_dtset%rfmeth/=1)write(nunit,'(3x,a9,3i10)')'   rfmeth',multibinit_dtset%rfmeth
2144    if(multibinit_dtset%enunit/=0)write(nunit,'(3x,a9,3i10)')'   enunit',multibinit_dtset%enunit
2145    if(multibinit_dtset%eivec/=0) write(nunit,'(3x,a9,3i10)')'    eivec',multibinit_dtset%eivec
2146    if(multibinit_dtset%asr/=0)   write(nunit,'(3x,a9,3i10)')'      asr',multibinit_dtset%asr
2147    if(multibinit_dtset%chneut/=0)write(nunit,'(3x,a9,3i10)')'   chneut',multibinit_dtset%chneut
2148  end if
2149 
2150 
2151 !For interatomic force constant information
2152  if(multibinit_dtset%ifcflag/=0)then
2153    write(nunit,'(a)')' Interatomic Force Constants Inputs :'
2154    write(nunit,'(3x,a9,3i10)')'   dipdip',multibinit_dtset%dipdip
2155    if(multibinit_dtset%dipdip /= 0)then
2156      write(nunit,'(a12,3i10)') 'dipdip_range',multibinit_dtset%dipdip_range
2157    end if
2158    if(multibinit_dtset%dipdip_prt/=0)then
2159      write(nunit,'(a12,3i10)') 'dipdip_prt',multibinit_dtset%dipdip_prt
2160    end if
2161    if(multibinit_dtset%nsphere/=0)write(nunit,'(3x,a9,3i10)')'  nsphere',multibinit_dtset%nsphere
2162    if(abs(multibinit_dtset%rifcsph)>tol10)write(nunit,'(3x,a9,E16.6)')'  nsphere',multibinit_dtset%rifcsph
2163    write(nunit,'(3x,a9,3i10)')'   ifcana',multibinit_dtset%ifcana
2164    write(nunit,'(3x,a9,3i10)')'   ifcout',multibinit_dtset%ifcout
2165    if(multibinit_dtset%natifc>=1)then
2166      write(nunit,'(3x,a9,3i10)')'   natifc',multibinit_dtset%natifc
2167      write(nunit,'(3x,a12)',advance='no')'    atifc   '
2168      write(nunit,'(3x,15i4)') (multibinit_dtset%atifc(ii)*ii,ii=1,multibinit_dtset%natifc)
2169 
2170    end if
2171    write(nunit,'(a)')' Description of grid 1 :'
2172    write(nunit,'(3x,a9,3i10)')'     brav',multibinit_dtset%brav
2173    write(nunit,'(3x,a9,3i10)')'    ngqpt',multibinit_dtset%ngqpt(1:3)
2174    write(nunit,'(3x,a9,3i10)')'   nqshft',multibinit_dtset%nqshft
2175    if (multibinit_dtset%nqshft/=0)then
2176      write(nunit,'(3x,a9)')'   q1shft'
2177      do iqshft=1,multibinit_dtset%nqshft
2178        write(nunit,'(19x,4es16.8)') (multibinit_dtset%q1shft(ii,iqshft),ii=1,3)
2179      end do
2180    end if
2181    if (any(multibinit_dtset%qrefine(:) > 1)) then
2182      write(nunit,'(3x,a9,3i10)')'  qrefine', multibinit_dtset%qrefine
2183    end if
2184  end if
2185 
2186 
2187 !List of vector 1  (reduced coordinates)
2188  if(multibinit_dtset%nph1l/=0)then
2189    write(nunit,'(a)')' First list of wavevector (reduced coord.) :'
2190    write(nunit,'(3x,a9,3i10)')'    nph1l',multibinit_dtset%nph1l
2191    write(nunit,'(3x,a9)')'    qph1l'
2192    do iph1=1,multibinit_dtset%nph1l
2193      write(nunit,'(19x,3es16.8,2x,es11.3)') &
2194 &     (multibinit_dtset%qph1l(ii,iph1),ii=1,3),multibinit_dtset%qnrml1(iph1)
2195    end do
2196  end if
2197 
2198 !List of vector 2  (cartesian coordinates)
2199  if(multibinit_dtset%nph2l/=0)then
2200    write(nunit,'(a)')' Second list of wavevector (cart. coord.) :'
2201    write(nunit,'(3x,a9,3i10)')'    nph2l',multibinit_dtset%nph2l
2202    write(nunit,'(3x,a9)')'    qph2l'
2203    do iph2=1,multibinit_dtset%nph2l
2204      write(nunit,'(19x,3es16.8,2x,es11.3)') &
2205 &     (multibinit_dtset%qph2l(ii,iph2),ii=1,3),multibinit_dtset%qnrml2(iph2)
2206    end do
2207  end if
2208 
2209  write(nunit,'(a,80a,a)') ch10,('=',ii=1,80),ch10
2210 
2211 end subroutine outvars_multibinit