! This file is part of xtb.
!
! Copyright (C) 2017-2020 Stefan Grimme
!
! xtb is free software: you can redistribute it and/or modify it under
! the terms of the GNU Lesser General Public License as published by
! the Free Software Foundation, either version 3 of the License, or
! (at your option) any later version.
!
! xtb is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
! GNU Lesser General Public License for more details.
!
! You should have received a copy of the GNU Lesser General Public License
! along with xtb.  If not, see <https://www.gnu.org/licenses/>.
module xtb_chargemodel
contains

subroutine new_multieq_model_2019(chrgeq,n,at)
   use xtb_mctc_accuracy, only : wp
   use xtb_type_param
   implicit none
   type(chrg_parameter) :: chrgeq
   integer,intent(in) :: n
   integer,intent(in) :: at(n)
!! ------------------------------------------------------------------------
!  PARAMETRISATION BY S. SPICHER (Mon, 25 Feb 2019 17:06:20 +0100)
!  atomic Hirshfeld dipole moments and partial charges were used in the fit.
!  beta = 3*alpha and dpolc was fixed to the GFN2-xTB values
!! ------------------------------------------------------------------------
   integer, parameter :: max_elem = 86
   real(wp),parameter :: en(max_elem) = [ &
    1.31222670_wp, 1.33049990_wp, 0.80127754_wp, 1.16005622_wp, 1.29157604_wp, &
    1.40450084_wp, 1.46638512_wp, 1.46801039_wp, 1.55348578_wp, 1.37248867_wp, &
    1.37937749_wp, 1.04419231_wp, 1.26263099_wp, 1.11162320_wp, 1.51830508_wp, &
    1.47998887_wp, 1.47394487_wp, 1.34302662_wp, 1.36729589_wp, 1.27794050_wp, &
    1.13199584_wp, 1.16946189_wp, 1.12955539_wp, 1.23794839_wp, 1.03635528_wp, &
    1.21633604_wp, 1.05065523_wp, 1.14699791_wp, 1.10218954_wp, 1.27617283_wp, &
    1.24498441_wp, 1.23286571_wp, 1.41590228_wp, 1.42601999_wp, 1.45841247_wp, &
    1.36946764_wp, 0.99022781_wp, 1.29623400_wp, 1.18855925_wp, 1.24177440_wp, &
    1.08778669_wp, 1.03379538_wp, 1.16318235_wp, 1.15329928_wp, 1.23938470_wp, &
    1.29946952_wp, 1.23925857_wp, 1.29607528_wp, 1.21496789_wp, 1.21020793_wp, &
    1.29904026_wp, 1.41603789_wp, 1.39750387_wp, 1.32362802_wp, 1.31960262_wp, &
    1.30739388_wp, 1.22385600_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, &
    0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, &
    0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, &
    0.56999999_wp, 1.08363456_wp, 1.25089063_wp, 1.28965951_wp, 1.39207680_wp, &
    1.19839481_wp, 1.25938591_wp, 1.34157715_wp, 1.28532469_wp, 1.32452604_wp, &
    1.22277498_wp, 0.92949663_wp, 1.27420937_wp, 1.37079453_wp, 1.37346232_wp, &
    1.38193459_wp ]
   real(wp),parameter :: gam(max_elem) = [ &
   -0.29352787_wp, 1.18228813_wp, 0.81655542_wp, 0.04991939_wp, 0.18077754_wp, &
    0.16236307_wp, 0.15729171_wp, 0.06664232_wp, 0.22420664_wp,-0.62957580_wp, &
    0.40078533_wp, 0.63368871_wp, 0.60184244_wp, 0.15902198_wp, 0.17066813_wp, &
    0.18321894_wp, 0.09316157_wp, 0.31180808_wp, 0.01366019_wp,-0.03419238_wp, &
    0.06025478_wp,-0.06527532_wp, 0.35264616_wp,-0.23266246_wp, 0.24660308_wp, &
    0.30889629_wp, 0.19765403_wp, 0.30423715_wp, 0.89032516_wp,-0.22477300_wp, &
    0.02452835_wp, 0.23273331_wp, 0.15504976_wp, 0.10946405_wp, 0.93073322_wp, &
   -0.40355129_wp, 0.57791901_wp, 0.21101308_wp, 0.11134256_wp, 0.03273833_wp, &
    0.06592765_wp, 0.09853050_wp, 0.02840088_wp,-0.77038368_wp, 0.33861167_wp, &
    0.06232847_wp, 0.17595304_wp,-0.62894354_wp,-0.03859051_wp,-0.14870423_wp, &
   -0.11119651_wp,-0.18203460_wp, 0.05382599_wp,-0.17608856_wp,-0.13200542_wp, &
   -0.23001704_wp, 0.21702733_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, &
    0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, &
    0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, &
    0.11000000_wp,-0.14950188_wp,-0.05301495_wp, 0.32340145_wp, 0.48443586_wp, &
    0.49482288_wp, 0.03181459_wp,-0.08968091_wp, 0.03036240_wp, 0.35118499_wp, &
    0.21196095_wp, 1.00879826_wp, 0.09954288_wp, 0.03013378_wp, 0.01667458_wp, &
    0.07860724_wp ]
   real(wp),parameter :: kappa(max_elem) = [ &
   -0.02597869_wp, 0.66325692_wp,-0.17760442_wp,-0.05425035_wp,-0.01899386_wp, &
    0.03454974_wp, 0.06220207_wp, 0.07518969_wp, 0.17450351_wp, 0.05294182_wp, &
    0.01400123_wp,-0.10199084_wp,-0.04100613_wp,-0.11309997_wp, 0.11651405_wp, &
    0.10094484_wp, 0.10874426_wp, 0.08168330_wp, 0.01390994_wp,-0.02130805_wp, &
   -0.09715052_wp,-0.07461454_wp,-0.08712643_wp,-0.05437354_wp,-0.11329620_wp, &
   -0.06603778_wp,-0.12139278_wp,-0.07741014_wp,-0.08424088_wp,-0.02682389_wp, &
   -0.04497087_wp,-0.05460499_wp, 0.04342269_wp, 0.04099114_wp, 0.09000060_wp, &
    0.00071681_wp,-0.15290585_wp, 0.06720034_wp,-0.05999122_wp,-0.03190377_wp, &
   -0.10955588_wp,-0.11437713_wp,-0.06969508_wp,-0.08803104_wp,-0.03470070_wp, &
   -0.00198100_wp,-0.02170821_wp,-0.01661924_wp,-0.05852491_wp,-0.06697080_wp, &
   -0.05895723_wp, 0.04217060_wp, 0.05124334_wp,-0.01482324_wp, 0.01491499_wp, &
   -0.01601696_wp,-0.00766699_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp, &
   -0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp, &
   -0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp, &
   -0.11000000_wp,-0.09612654_wp,-0.05424292_wp,-0.03255023_wp, 0.01293739_wp, &
   -0.04887012_wp,-0.03931052_wp, 0.00161909_wp,-0.03376011_wp,-0.00976174_wp, &
   -0.05304306_wp, 0.01236187_wp,-0.02297437_wp, 0.03918072_wp, 0.05434160_wp, &
    0.00830592_wp ]
   real(wp),parameter :: alpha(max_elem) = [ &
    0.83051332_wp, 0.27934189_wp, 0.57717327_wp, 5.41127387_wp, 2.09405910_wp, &
    3.77537318_wp, 3.48318023_wp, 2.54214604_wp, 2.73404785_wp, 0.16750746_wp, &
    1.15468837_wp, 3.42988344_wp, 1.01513110_wp, 3.11624566_wp, 2.62048544_wp, &
    1.09480440_wp, 3.46065237_wp, 4.48704043_wp, 3.88654194_wp, 1.75864066_wp, &
    1.95337522_wp, 1.49408444_wp, 1.38950107_wp, 1.31196889_wp, 3.04235459_wp, &
    0.71788381_wp, 1.15658896_wp, 2.27724965_wp, 0.86379032_wp, 0.81398266_wp, &
    0.66709843_wp, 0.87908763_wp, 3.24925489_wp, 1.58803868_wp, 1.06623475_wp, &
    0.98174111_wp, 1.75011798_wp, 5.42115002_wp, 3.91931444_wp, 3.11270895_wp, &
    4.01156044_wp, 4.09746782_wp, 2.87234516_wp, 0.43018109_wp, 2.19964334_wp, &
    2.69978488_wp, 5.23629932_wp, 0.72482212_wp, 1.59437558_wp, 1.75692514_wp, &
    0.94058335_wp, 1.06309445_wp, 3.51583993_wp, 1.60898068_wp, 0.72873185_wp, &
    1.47627148_wp, 3.99569972_wp, 1.64000000_wp, 1.64000000_wp, 1.64000000_wp, &
    1.64000000_wp, 1.64000000_wp, 1.64000000_wp, 1.64000000_wp, 1.64000000_wp, &
    1.64000000_wp, 1.64000000_wp, 1.64000000_wp, 1.64000000_wp, 1.64000000_wp, &
    1.64000000_wp, 1.19273782_wp, 1.17975130_wp, 0.47635729_wp, 3.33070704_wp, &
    3.26261926_wp, 2.81175921_wp, 1.67471893_wp, 3.42329532_wp, 1.26308728_wp, &
    1.71990410_wp, 2.78739037_wp, 4.45639288_wp, 2.68108583_wp, 2.31545094_wp, &
    3.72300764_wp ]
   real(wp),parameter :: dpolc(max_elem) = [ &
    0.5563889E-01_wp,-0.1000000E-01_wp,-0.5000000E-02_wp,-0.6133410E-02_wp, &
   -0.4811860E-02_wp,-0.4116740E-02_wp, 0.3521273E-01_wp,-0.4935670E-01_wp, &
   -0.8339183E-01_wp, 0.1000000E+00_wp, 0.0000000E+00_wp,-0.8200500E-03_wp, &
    0.2633341E-01_wp,-0.2575000E-03_wp, 0.2110225E-01_wp,-0.1511170E-02_wp, &
   -0.2536958E-01_wp,-0.2077329E-01_wp,-0.1033830E-02_wp,-0.2366750E-02_wp, &
   -0.5151770E-02_wp,-0.4345060E-02_wp,-0.3500000E-02_wp, 0.1496690E-02_wp, &
   -0.7591680E-02_wp, 0.4129290E-02_wp,-0.2479380E-02_wp,-0.1261887E-01_wp, &
   -0.7000000E-02_wp,-0.1000000E-02_wp, 0.2672190E-02_wp, 0.1084600E-02_wp, &
   -0.2012940E-02_wp,-0.2886480E-02_wp,-0.1088586E-01_wp,-0.8893570E-02_wp, &
   -0.9332800E-03_wp,-0.4599250E-02_wp,-0.6372910E-02_wp,-0.5996150E-02_wp, &
   -0.2887290E-02_wp, 0.3463270E-02_wp,-0.4584160E-02_wp,-0.8192200E-03_wp, &
    0.7016000E-04_wp,-0.3103610E-02_wp,-0.8003140E-02_wp,-0.1053640E-02_wp, &
    0.9510790E-02_wp, 0.8502900E-03_wp,-0.1551900E-03_wp,-0.2634140E-02_wp, &
   -0.6036480E-02_wp,-0.2144470E-02_wp,-0.8000000E-03_wp,-0.2600000E-02_wp, &
   -0.3951980E-02_wp,-0.7238060E-02_wp,-0.7048190E-02_wp,-0.6858320E-02_wp, &
   -0.6668450E-02_wp,-0.6478580E-02_wp,-0.6288710E-02_wp,-0.6098840E-02_wp, &
   -0.5908970E-02_wp,-0.5719100E-02_wp,-0.5529230E-02_wp,-0.5339360E-02_wp, &
   -0.5149490E-02_wp,-0.4959610E-02_wp,-0.4769740E-02_wp,-0.5376850E-02_wp, &
   -0.2003430E-02_wp, 0.6588600E-03_wp,-0.5876360E-02_wp,-0.5100900E-02_wp, &
   -0.6738220E-02_wp,-0.4236840E-02_wp, 0.3934180E-02_wp,-0.2500000E-02_wp, &
    0.3740180E-02_wp, 0.1007016E-01_wp,-0.7372520E-02_wp,-0.1344854E-01_wp, &
   -0.3481230E-02_wp,-0.1675970E-02_wp ]
   integer :: i

   call chrgeq%allocate(n,extended=.true.)

   do i = 1, n
      chrgeq%en   (i) = en   (at(i))
      chrgeq%gam  (i) = gam  (at(i))
      chrgeq%kappa(i) = kappa(at(i))
      chrgeq%alpha(i) = alpha(at(i))
      chrgeq%beta (i) = alpha(at(i))*3
      chrgeq%dpol (i) = dpolc(at(i))
   enddo

end subroutine new_multieq_model_2019

subroutine new_charge_model_2019(chrgeq,n,at)
   use xtb_mctc_accuracy, only : wp
   use xtb_type_param
   implicit none
   type(chrg_parameter) :: chrgeq
   integer,intent(in) :: n
   integer,intent(in) :: at(n)
!! ------------------------------------------------------------------------
!  PARAMETRISATION BY S. SPICHER (Fri, 14 Dec 2018 16:13:08 +0100)
!! ------------------------------------------------------------------------
   integer, parameter :: max_elem = 86
   real(wp),parameter :: en(max_elem) = (/ &
    1.23695041_wp, 1.26590957_wp, 0.54341808_wp, 0.99666991_wp, 1.26691604_wp,&
    1.40028282_wp, 1.55819364_wp, 1.56866440_wp, 1.57540015_wp, 1.15056627_wp,&
    0.55936220_wp, 0.72373742_wp, 1.12910844_wp, 1.12306840_wp, 1.52672442_wp,&
    1.40768172_wp, 1.48154584_wp, 1.31062963_wp, 0.40374140_wp, 0.75442607_wp,&
    0.76482096_wp, 0.98457281_wp, 0.96702598_wp, 1.05266584_wp, 0.93274875_wp,&
    1.04025281_wp, 0.92738624_wp, 1.07419210_wp, 1.07900668_wp, 1.04712861_wp,&
    1.15018618_wp, 1.15388455_wp, 1.36313743_wp, 1.36485106_wp, 1.39801837_wp,&
    1.18695346_wp, 0.36273870_wp, 0.58797255_wp, 0.71961946_wp, 0.96158233_wp,&
    0.89585296_wp, 0.81360499_wp, 1.00794665_wp, 0.92613682_wp, 1.09152285_wp,&
    1.14907070_wp, 1.13508911_wp, 1.08853785_wp, 1.11005982_wp, 1.12452195_wp,&
    1.21642129_wp, 1.36507125_wp, 1.40340000_wp, 1.16653482_wp, 0.34125098_wp,&
    0.58884173_wp, 0.68441115_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp,&
    0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp,&
    0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp, 0.56999999_wp,&
    0.56999999_wp, 0.87936784_wp, 1.02761808_wp, 0.93297476_wp, 1.10172128_wp,&
    0.97350071_wp, 1.16695666_wp, 1.23997927_wp, 1.18464453_wp, 1.14191734_wp,&
    1.12334192_wp, 1.01485321_wp, 1.12950808_wp, 1.30804834_wp, 1.33689961_wp,&
    1.27465977_wp /)
   real(wp),parameter :: gam(max_elem) = (/ &
   -0.35015861_wp, 1.04121227_wp, 0.09281243_wp, 0.09412380_wp, 0.26629137_wp,&
    0.19408787_wp, 0.05317918_wp, 0.03151644_wp, 0.32275132_wp, 1.30996037_wp,&
    0.24206510_wp, 0.04147733_wp, 0.11634126_wp, 0.13155266_wp, 0.15350650_wp,&
    0.15250997_wp, 0.17523529_wp, 0.28774450_wp, 0.42937314_wp, 0.01896455_wp,&
    0.07179178_wp,-0.01121381_wp,-0.03093370_wp, 0.02716319_wp,-0.01843812_wp,&
   -0.15270393_wp,-0.09192645_wp,-0.13418723_wp,-0.09861139_wp, 0.18338109_wp,&
    0.08299615_wp, 0.11370033_wp, 0.19005278_wp, 0.10980677_wp, 0.12327841_wp,&
    0.25345554_wp, 0.58615231_wp, 0.16093861_wp, 0.04548530_wp,-0.02478645_wp,&
    0.01909943_wp, 0.01402541_wp,-0.03595279_wp, 0.01137752_wp,-0.03697213_wp,&
    0.08009416_wp, 0.02274892_wp, 0.12801822_wp,-0.02078702_wp, 0.05284319_wp,&
    0.07581190_wp, 0.09663758_wp, 0.09547417_wp, 0.07803344_wp, 0.64913257_wp,&
    0.15348654_wp, 0.05054344_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp,&
    0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp,&
    0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp, 0.11000000_wp,&
    0.11000000_wp,-0.02786741_wp, 0.01057858_wp,-0.03892226_wp,-0.04574364_wp,&
   -0.03874080_wp,-0.03782372_wp,-0.07046855_wp, 0.09546597_wp, 0.21953269_wp,&
    0.02522348_wp, 0.15263050_wp, 0.08042611_wp, 0.01878626_wp, 0.08715453_wp,&
    0.10500484_wp /)
   real(wp),parameter :: kappa(max_elem) = (/ &
    0.04916110_wp, 0.10937243_wp,-0.12349591_wp,-0.02665108_wp,-0.02631658_wp,&
    0.06005196_wp, 0.09279548_wp, 0.11689703_wp, 0.15704746_wp, 0.07987901_wp,&
   -0.10002962_wp,-0.07712863_wp,-0.02170561_wp,-0.04964052_wp, 0.14250599_wp,&
    0.07126660_wp, 0.13682750_wp, 0.14877121_wp,-0.10219289_wp,-0.08979338_wp,&
   -0.08273597_wp,-0.01754829_wp,-0.02765460_wp,-0.02558926_wp,-0.08010286_wp,&
   -0.04163215_wp,-0.09369631_wp,-0.03774117_wp,-0.05759708_wp, 0.02431998_wp,&
   -0.01056270_wp,-0.02692862_wp, 0.07657769_wp, 0.06561608_wp, 0.08006749_wp,&
    0.14139200_wp,-0.05351029_wp,-0.06701705_wp,-0.07377246_wp,-0.02927768_wp,&
   -0.03867291_wp,-0.06929825_wp,-0.04485293_wp,-0.04800824_wp,-0.01484022_wp,&
    0.07917502_wp, 0.06619243_wp, 0.02434095_wp,-0.01505548_wp,-0.03030768_wp,&
    0.01418235_wp, 0.08953411_wp, 0.08967527_wp, 0.07277771_wp,-0.02129476_wp,&
   -0.06188828_wp,-0.06568203_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,&
   -0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,&
   -0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,&
   -0.11000000_wp,-0.03585873_wp,-0.03132400_wp,-0.05902379_wp,-0.02827592_wp,&
   -0.07606260_wp,-0.02123839_wp, 0.03814822_wp, 0.02146834_wp, 0.01580538_wp,&
   -0.00894298_wp,-0.05864876_wp,-0.01817842_wp, 0.07721851_wp, 0.07936083_wp,&
    0.05849285_wp /)
   real(wp),parameter :: alpha(max_elem) = (/ &
    0.55159092_wp, 0.66205886_wp, 0.90529132_wp, 1.51710827_wp, 2.86070364_wp,&
    1.88862966_wp, 1.32250290_wp, 1.23166285_wp, 1.77503721_wp, 1.11955204_wp,&
    1.28263182_wp, 1.22344336_wp, 1.70936266_wp, 1.54075036_wp, 1.38200579_wp,&
    2.18849322_wp, 1.36779065_wp, 1.27039703_wp, 1.64466502_wp, 1.58859404_wp,&
    1.65357953_wp, 1.50021521_wp, 1.30104175_wp, 1.46301827_wp, 1.32928147_wp,&
    1.02766713_wp, 1.02291377_wp, 0.94343886_wp, 1.14881311_wp, 1.47080755_wp,&
    1.76901636_wp, 1.98724061_wp, 2.41244711_wp, 2.26739524_wp, 2.95378999_wp,&
    1.20807752_wp, 1.65941046_wp, 1.62733880_wp, 1.61344972_wp, 1.63220728_wp,&
    1.60899928_wp, 1.43501286_wp, 1.54559205_wp, 1.32663678_wp, 1.37644152_wp,&
    1.36051851_wp, 1.23395526_wp, 1.65734544_wp, 1.53895240_wp, 1.97542736_wp,&
    1.97636542_wp, 2.05432381_wp, 3.80138135_wp, 1.43893803_wp, 1.75505957_wp,&
    1.59815118_wp, 1.76401732_wp, 1.63999999_wp, 1.63999999_wp, 1.63999999_wp,&
    1.63999999_wp, 1.63999999_wp, 1.63999999_wp, 1.63999999_wp, 1.63999999_wp,&
    1.63999999_wp, 1.63999999_wp, 1.63999999_wp, 1.63999999_wp, 1.63999999_wp,&
    1.63999999_wp, 1.47055223_wp, 1.81127084_wp, 1.40189963_wp, 1.54015481_wp,&
    1.33721475_wp, 1.57165422_wp, 1.04815857_wp, 1.78342098_wp, 2.79106396_wp,&
    1.78160840_wp, 2.47588882_wp, 2.37670734_wp, 1.76613217_wp, 2.66172302_wp,&
    2.82773085_wp /)
   integer :: i

   call chrgeq%allocate(n)

   do i = 1, n
      chrgeq%en   (i) = en   (at(i))
      chrgeq%gam  (i) = gam  (at(i))
      chrgeq%kappa(i) = kappa(at(i))
      chrgeq%alpha(i) = alpha(at(i))
   enddo

end subroutine new_charge_model_2019

subroutine new_charge_model_2018(chrgeq,n,at)
   use xtb_mctc_accuracy, only : wp
   use xtb_type_param
   implicit none
   type(chrg_parameter) :: chrgeq
   integer,intent(in) :: n
   integer,intent(in) :: at(n)
!! ------------------------------------------------------------------------
!  PARAMETRISATION BY S. SPICHER          (Fri, 19 Oct 2018 12:15:49 +0200)
!! ------------------------------------------------------------------------
   integer, parameter :: max_elem = 86
   real(wp),parameter :: en(max_elem) = (/ &
    1.23246797_wp, 1.23678380_wp, 0.54354054_wp, 1.07316582_wp, 1.13569914_wp,&
    1.37937538_wp, 1.53150709_wp, 1.54032780_wp, 1.46882833_wp, 1.13302128_wp,&
    0.63327514_wp, 0.68730302_wp, 0.99823617_wp, 1.07950372_wp, 1.50449652_wp,&
    1.39951019_wp, 1.52598012_wp, 1.27907542_wp, 0.50724547_wp, 0.85250848_wp,&
    0.63084942_wp, 0.97122001_wp, 0.89213568_wp, 1.06848779_wp, 0.88675761_wp,&
    0.92643940_wp, 1.07658822_wp, 0.79667533_wp, 1.02060727_wp, 0.99094458_wp,&
    1.13024862_wp, 1.14478047_wp, 1.34552412_wp, 1.36515367_wp, 1.44560491_wp,&
    1.15603238_wp, 0.30776713_wp, 0.55488955_wp, 0.72061253_wp, 0.82208038_wp,&
    0.91871553_wp, 0.74958769_wp, 0.98634094_wp, 0.93563427_wp, 1.00922031_wp,&
    1.08092750_wp, 1.09077840_wp, 1.10155149_wp, 0.97377281_wp, 1.15043898_wp,&
    1.20792358_wp, 1.36527027_wp, 1.38690869_wp, 1.12021226_wp, 0.37885185_wp,&
    0.81197551_wp, 0.60282207_wp, 0.60300000_wp, 0.60300000_wp, 0.60300000_wp,&
    0.60300000_wp, 0.60300000_wp, 0.60300000_wp, 0.60300000_wp, 0.60300000_wp,&
    0.60300000_wp, 0.60300000_wp, 0.60300000_wp, 0.60300000_wp, 0.60300000_wp,&
    0.60300000_wp, 0.60457044_wp, 1.03860983_wp, 1.05078223_wp, 1.16035684_wp,&
    1.21782330_wp, 1.28662970_wp, 1.25584923_wp, 1.17736525_wp, 1.15718537_wp,&
    1.05728283_wp, 1.01478684_wp, 1.10230954_wp, 1.23575846_wp, 1.29741267_wp,&
    1.26158938_wp/)
   real(wp),parameter :: gam(max_elem) = (/ &
   -0.27489752_wp, 1.00704831_wp, 0.34071460_wp, 0.08105579_wp, 0.23471517_wp,&
    0.12761720_wp, 0.05932592_wp, 0.08277108_wp, 0.34239102_wp, 0.78945378_wp,&
    0.28268286_wp, 0.15198052_wp, 0.15391686_wp, 0.14203269_wp, 0.27190586_wp,&
    0.11133624_wp, 0.46712500_wp, 0.40340963_wp, 0.43192407_wp, 0.13848799_wp,&
    0.18409419_wp, 0.08022090_wp, 0.12834219_wp, 0.11194514_wp, 0.11124648_wp,&
    0.09114324_wp, 0.01245244_wp, 0.03607465_wp, 0.06443861_wp, 0.23075980_wp,&
    0.10908422_wp, 0.08345121_wp, 0.17433437_wp, 0.12874012_wp, 0.19799687_wp,&
    0.35525360_wp, 0.43550143_wp, 0.24189880_wp, 0.09134220_wp, 0.11710135_wp,&
    0.06107508_wp, 0.08816239_wp, 0.05368855_wp, 0.08338086_wp, 0.05683249_wp,&
    0.13991799_wp, 0.12938431_wp, 0.09081899_wp, 0.08946107_wp, 0.04486041_wp,&
    0.13635071_wp, 0.15548135_wp, 0.09717308_wp, 0.17235987_wp, 0.63306207_wp,&
    0.21376153_wp, 0.12277197_wp, 0.12700000_wp, 0.12700000_wp, 0.12700000_wp,&
    0.12700000_wp, 0.12700000_wp, 0.12700000_wp, 0.12700000_wp, 0.12700000_wp,&
    0.12700000_wp, 0.12700000_wp, 0.12700000_wp, 0.12700000_wp, 0.12700000_wp,&
    0.12700000_wp, 0.13318693_wp, 0.05965057_wp, 0.01671840_wp, 0.02491941_wp,&
    0.01351197_wp, 0.02341562_wp, 0.04219035_wp, 0.07839352_wp, 0.15950374_wp,&
    0.17359289_wp, 0.16380548_wp, 0.11025180_wp, 0.05546045_wp, 0.06672905_wp,&
    0.11297489_wp/)
   real(wp),parameter :: kappa(max_elem) = (/ &
    0.01422129_wp, 0.36987075_wp,-0.12738757_wp,-0.03351538_wp,-0.04216033_wp,&
    0.06676096_wp, 0.14070447_wp, 0.14684002_wp, 0.11360012_wp, 0.24745874_wp,&
   -0.09997174_wp,-0.10191510_wp,-0.05658005_wp,-0.03042332_wp, 0.17934618_wp,&
    0.09229482_wp, 0.19499461_wp, 0.19494756_wp,-0.05405802_wp,-0.07626752_wp,&
   -0.12341715_wp,-0.03639558_wp,-0.09048218_wp,-0.06284416_wp,-0.13515558_wp,&
   -0.10699213_wp,-0.06715865_wp,-0.20918026_wp,-0.07667280_wp,-0.00551953_wp,&
    0.00027244_wp,-0.01922636_wp, 0.08732735_wp, 0.08486002_wp, 0.13964063_wp,&
    0.16391237_wp,-0.11757842_wp,-0.11185355_wp,-0.10742621_wp,-0.07585379_wp,&
   -0.06280485_wp,-0.12142065_wp,-0.06814013_wp,-0.07531689_wp,-0.06075124_wp,&
    0.01955447_wp, 0.04077176_wp, 0.01340060_wp,-0.05917210_wp,-0.00740372_wp,&
    0.02571361_wp, 0.09533328_wp, 0.08639724_wp, 0.07409935_wp,-0.08851343_wp,&
   -0.02042462_wp,-0.10384451_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,&
   -0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,&
   -0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,-0.11000000_wp,&
   -0.11000000_wp,-0.12184554_wp,-0.05144041_wp,-0.06006418_wp,-0.02817541_wp,&
   -0.00792045_wp, 0.01627009_wp, 0.00884773_wp,-0.00049183_wp, 0.02652209_wp,&
   -0.01258868_wp,-0.04974423_wp,-0.02081982_wp, 0.02927793_wp, 0.05314950_wp,&
    0.06732131_wp/)
   real(wp),parameter :: alpha(max_elem) = (/ &
    0.57707019_wp, 0.51575989_wp, 0.99515576_wp, 1.80197223_wp, 2.31419253_wp,&
    1.76206788_wp, 1.30582582_wp, 1.25023901_wp, 1.81619582_wp, 0.55372930_wp,&
    1.36990551_wp, 1.26851926_wp, 2.22945756_wp, 1.72383310_wp, 1.71034451_wp,&
    1.78237180_wp, 2.06376648_wp, 1.40024022_wp, 1.25534997_wp, 2.63503044_wp,&
    1.38566080_wp, 1.57025576_wp, 1.37262637_wp, 1.79160917_wp, 1.59635110_wp,&
    1.01726540_wp, 1.36035014_wp, 1.08587893_wp, 1.17481116_wp, 1.47067441_wp,&
    2.05436095_wp, 2.22351909_wp, 2.34304426_wp, 2.25772299_wp, 3.05020818_wp,&
    1.29437190_wp, 1.06290969_wp, 1.84660949_wp, 1.58520176_wp, 1.60093239_wp,&
    1.68794490_wp, 1.36082945_wp, 1.87731510_wp, 1.42281352_wp, 1.48513495_wp,&
    1.38964955_wp, 1.35988990_wp, 1.85268772_wp, 1.83356928_wp, 2.01115973_wp,&
    2.31017191_wp, 2.39609492_wp, 3.42888191_wp, 1.54613874_wp, 1.41710487_wp,&
    1.97377465_wp, 1.61360236_wp, 1.40000000_wp, 1.40000000_wp, 1.40000000_wp,&
    1.40000000_wp, 1.40000000_wp, 1.40000000_wp, 1.40000000_wp, 1.40000000_wp,&
    1.40000000_wp, 1.40000000_wp, 1.40000000_wp, 1.40000000_wp, 1.40000000_wp,&
    1.40000000_wp, 1.23539671_wp, 1.95821108_wp, 1.70664862_wp, 1.89862769_wp,&
    1.85150533_wp, 1.96977793_wp, 1.44086039_wp, 2.12023827_wp, 2.32373422_wp,&
    2.25901981_wp, 2.45035320_wp, 2.57320085_wp, 2.15313970_wp, 2.53046115_wp,&
    2.57377591_wp/)
   integer :: i

   call chrgeq%allocate(n)

   do i = 1, n
      chrgeq%en   (i) = en   (at(i))
      chrgeq%gam  (i) = gam  (at(i))
      chrgeq%kappa(i) = kappa(at(i))
      chrgeq%alpha(i) = alpha(at(i))
   enddo

end subroutine new_charge_model_2018

subroutine gfn0_charge_model(chrgeq,n,at,jData)
   use xtb_mctc_accuracy, only : wp
   use xtb_type_param
   use xtb_xtb_data
   implicit none
   type(TCoulombData), intent(in) :: jData
   type(chrg_parameter) :: chrgeq
   integer, intent(in)  :: n
   integer, intent(in)  :: at(n)
   integer :: i

   call chrgeq%allocate(n)

   do i = 1, n
      chrgeq%en   (i) = jData%electronegativity(at(i))
      chrgeq%gam  (i) = jData%chemicalHardness(at(i))
      chrgeq%kappa(i) = jData%kCN(at(i))
      chrgeq%alpha(i) = jData%chargeWidth(at(i))
   enddo

end subroutine gfn0_charge_model

subroutine read_charge_model(chrgeq,n,at,fname)
   use xtb_mctc_accuracy, only : wp
   use xtb_type_param
   implicit none
   type(chrg_parameter) :: chrgeq
   integer, intent(in)  :: n
   integer, intent(in)  :: at(n)
   integer, parameter   :: max_elem = 86
   character(len=*),intent(in) :: fname
   integer :: i,ifile,err
   logical :: exist

   real(wp) :: en(max_elem)    = 0.0_wp
   real(wp) :: gam(max_elem)   = 0.0_wp
   real(wp) :: kappa(max_elem) = 0.0_wp
   real(wp) :: alpha(max_elem) = 0.0_wp
   real(wp) :: beta(max_elem)  = 0.0_wp
   real(wp) :: dpolc(max_elem) = 0.0_wp

   call chrgeq%allocate(n,extended=.true.)

   inquire(file=fname,exist=exist)
   open(newunit=ifile,file=fname)
   do i = 1, max_elem
      read(ifile,*,iostat=err) en(i),gam(i),kappa(i),alpha(i),beta(i),dpolc(i)
      if (err.ne.0) exit
   enddo
   close(ifile)

   do i = 1, n
      chrgeq%en   (i) = en   (at(i))
      chrgeq%gam  (i) = gam  (at(i))
      chrgeq%kappa(i) = kappa(at(i))
      chrgeq%alpha(i) = alpha(at(i))
      chrgeq%beta (i) = beta (at(i))
      chrgeq%dpol (i) = dpolc(at(i))
   enddo

end subroutine read_charge_model

end module xtb_chargemodel
