******************************************************************************** * * * 2DHF version 1-2003 * * Copyright (C) 1996 Jacek Kobus, Leif Laaksonen, Dage Sundholm * * * * This software may be used and distributed according to the terms * * of the GNU General Public License, see README and COPYING. * * * ******************************************************************************** c ### iniaddr ### c c Determines dimensions of arrays and addresses of particular c orbitals and coulomb potentials in cw_orb and cw_coul arrays c addresses of exchange potentials v(i,j) in the cw_exch array. c c With each cwork array are associated its address arrays: c c cw_orb - iaddr1(iorb,ip) c cw_coul - iaddr2(icoul,ip) c cw_exch - iaddr3(iexch,ip) c cw_suppl - iaddr4(isuppl) c cw_sctch - iaddr5(iscratch) c c In the case od iaddr1 the first index refers to orbital number c and the second one to the division point, i.e. c iaddr1(iorb,1)=i1b(iorb) c iaddr1(iorb,2)=i1e(iorb) c iaddr1(iorb,3)=i1si(iorb) c iaddr1(iorb,4)=i1ng(iorb) c c i1addr1(iorb,1) - address of the first value of orbital c iorb within cw_orb array c (address of the first value of orbital c iorb defined on grid #1) c (starting address of orbital iorb) c iaddr1(iorb,2) - address of the last value of orbital c iorb within cw_orb array c iaddr1(iorb,3) - number of values defining orbital iorb c iaddr1(iorb,4) - number of grids used by the orbital iorb c (up to 9 grids can be used) c iaddr1(iorb,ip) - ip=5-10 reserved for further use c c Similar definitions hold for iaddr2-3. c We assume no more than 60 orbitals and Coulomb potentials and c consequently no more than 1830 (iorb*(iorb+1)/2) exchange potentials. c subroutine iniaddr implicit integer*4 (i-n) implicit real*8 (a-h,o-z) include 'commons8.inc' c warning!!!! istrt(1)=0 and i1b(1)=1, i2b(1)=1 i1b(1)=1 ngrid=nni*iemu(i1ng(1)) i1e(1)=ngrid i1si(1)=ngrid i1mu(1)=iemu(i1ng(1)) c current value of cw_orb length is stored in l1cur l1cur=ngrid do iorb=2,norb i1b(iorb)=i1e(iorb-1)+1 ngrid=nni*iemu(i1ng(iorb)) i1e(iorb)=i1b(iorb)+ngrid-1 i1si(iorb)=ngrid i1mu(iorb)=iemu(i1ng(iorb)) l1cur=l1cur+ngrid enddo c it is assumed for a time being that potential functions are defined c on the same number of grids as the orbitals i2b(1)=1 ngrid=nni*iemu(i2ng(1)) i2e(1)=ngrid i2si(1)=ngrid i2mu(1)=iemu(i1ng(1)) c current value of cw_coul length is stored in l2cur l2cur=ngrid do iorb=2,norb i2b(iorb)=i2e(iorb-1)+1 ngrid=nni*iemu(i2ng(iorb)) i2e(iorb)=i2b(iorb)+ngrid-1 i2si(iorb)=ngrid i2mu(iorb)=iemu(i2ng(iorb)) l2cur=l2cur+ngrid enddo c warning! the value of k=i+i*(i-1)/2 for sigma type orbitals is undefined c dimension of i3??? is norb*(norb+1)/2 to allow for two c exchange potentials between non-sigma orbitals c exchange potential corresponding to orbitals iorb1 and iorb2 c defined on 1ing(iorb1) and i1ng(iorb2) grids respectively, is c defined on smaller of these grids, i.e i3ng=min(i1ng(iorb1),i1ng(iorb2)) iend=0 l=0 l3cur=0 i3boffset=0 do iorb1=1,norb do iorb2=iorb1,norb k=iorb1+iorb2*(iorb2-1)/2 ilc(k)=0 i3xk(iorb1,iorb2)=0 i3xk(iorb2,iorb1)=0 if (iorb1.eq.iorb2.and.ll(iorb1).eq.0) goto 50 ilc(k)=1 l=l+1 ngds=min(i1ng(iorb1),i1ng(iorb2)) i3ng(k)=ngds ngrid=nni*iemu(ngds) iend=iend+ngrid i3e(k)=iend i3b(k)=i3e(k)-ngrid+1 i3si(k)=ngrid i3mu(k)=iemu(i3ng(k)) i3bfin=i3e(k) i3xk(iorb1,iorb2)=k i3xk(iorb2,iorb1)=k i3orb1(k)=iorb1 i3orb2(k)=iorb2 l3cur=l3cur+ngrid if (iorb1.eq.iorb2) goto 50 if (ll(iorb1).eq.0.or.ll(iorb2).eq.0) goto 50 ilc(k)=2 l=l+1 iend=iend+ngrid l3cur=l3cur+ngrid if (k.gt.1830) then write(*,*) 'iniaddr:' write(*,*) 'address array for exchange potentials ', & 'is too short' stop 'iniaddr' endif 00050 continue enddo i3boffset=i3bfin enddo nexch=l c When exchange potential functions are stored as separate files each c containing a single potential then their labels are are generated c in the fashon defined below. c c prepare a two-dim array containing numerical labels of c exchange potential functions being needed for every orbital c c all exchange potentials have to be defined on the same grid c c let k=i3xk(iorb1,iorb2) be an index associated uniquely with a given c pair of orbitals c if k=0 - no ex.pot (no record) c if k.ne.0 and ilc(k)=1 - ex.pot record label is stored in i3xrec1 c if k.ne.0 and ilc(k)=2 - ex.pot record label is stored in i3xrec2 ngrid=nni*mxnmu l=0 do iorb1=1,norb ibeg=1 do iorb2=iorb1,norb k=i3xk(iorb1,iorb2) if (k.ne.0) then if (ilc(k).ne.0) then l=l+1 i3xrec1(k)=l if (ilc(k).eq.2) then l=l+1 i3xrec2(k)=l endif endif endif enddo enddo c determine exchange potentials associated with a particular orbital c and store the record lables in i3btv array; the number of ex. pot. for c each orbital is kept in i3nexp array do iorb1=1,norb l=0 ibeg=1 do iorb2=1,maxorb i3breck(iorb1,iorb2)=0 enddo do iorb2=1,norb k=i3xk(iorb1,iorb2) if (k.ne.0) then if (ilc(k).ne.0) then l=l+1 if (l.gt.maxorb) then print *,'iniaddr: the second dimension of arrays ', & 'i3xpair,i3brec can not exceed',maxorb stop endif i3xpair(iorb1,l)=i3xrec1(k) i3brec (iorb1,l)=ibeg i3breck(iorb1,l)=k i3btv(iorb1,iorb2)=ibeg i3btv(iorb2,iorb1)=ibeg c before a new portion of exc. pot.is being read in exc.pot c address array i3b has to be modified accordingly, i.e. c i3b(k) (where k is determined by i3breck and must be non zero) i3nexcp(iorb1)=l ibeg=ibeg+ngrid if (ilc(k).eq.2) then l=l+1 i3xpair(iorb1,l)=i3xrec2(k) i3brec(iorb1,l)=ibeg i3breck(iorb1,l)=k i3nexcp(iorb1)=l ibeg=ibeg+ngrid endif endif endif enddo enddo c checking i3btv entries iprint=0 if (iprint.ne.0) then print *,'iniaddr: ngrid ',nni*mxnmu print *,'iniaddr: nexch ',nexch print *, 'iniaddr: iorb1,iorb2,k,i3xk,ilc(k),i3b' do iorb1=1,norb do iorb2=iorb1,norb k=iorb1+iorb2*(iorb2-1)/2 if (iorb1.eq.iorb2) print *,'===' write(*,1200) iorb1,iorb2,k,i3xk(iorb1,iorb2),ilc(k), & i3b(k) enddo print *, ' ' enddo print *,'iaddr: iorb,i3nexcp' do iorb1=1,norb print *,iorb1,i3nexcp(iorb1) enddo print *,' ' print *,'iaddr: iorb // l,i3xpair,i3brec,i3breck' do iorb1=1,norb print *,'iorb1: ',iorb1 do l=1,i3nexcp(iorb1) k=i3breck(iorb1,l) write(*,1210) l, i3xpair(iorb1,l), & i3brec(iorb1,l),k enddo enddo 01200 format(2i4,i6,3i12) 01210 format(i4,5i12) print *,'iaddr: iorb1,iorb2,k,i3xrec1,i3xrec2' do iorb1=1,norb do iorb2=iorb1,norb k=i3xk(iorb1,iorb2) if (ilc(k).ne.0) then print *,iorb1,iorb2,k,i3xrec1(k),i3xrec2(k) endif enddo enddo print *, ' ' print *,'iorb1 // iorb2,k,irec,i3beg' do iorb1=1,norb print *,' ' print *, 'iorb1: ',iorb1 do iorb2=1,i3nexcp(iorb1) k=i3breck(iorb1,iorb2) i3beg=i3brec(iorb1,iorb2) irec=i3xpair(iorb1,iorb2) ngrid=i3si(k) print *,iorb2,k,irec,i3beg enddo enddo endif c is working array cw_suppl large enough? if (14.gt.length4) then write(iout6,*) 'cw_suppl is too short!' write(iout6,*) 'declared length is:',length4 write(iout6,*) 'needed length is 14' stop 'prepvar' endif do i=1,14 i4b(i)=(i-1)*mxsize+1 i4e(i)=iorb*mxsize i4si(i)=mxsize i4ng(i)=ngrids enddo c current value of cw_sctch length is stored in l5cur mxsize8=(nni+8)*(mxnmu+8) do i=1,10 i5b(i)=(i-1)*mxsize8+1 i5e(i)=i*mxsize8 i5si(i)=mxsize8 i5ng(i)=ngrids enddo l5cur=i5e(20) c are working arrays large enough? if (l1cur.gt.length1) then write(iout6,*) 'cw_orb is too short!' write(iout6,*) 'declared length is:',length1 write(iout6,*) 'needed length is:',l1cur if (idbg(550).eq.0) stop 'iniaddr' endif if (l2cur.gt.length2) then write(iout6,*) 'cw_coul is too short!' write(iout6,*) 'declared length is:',length2 write(iout6,*) 'needed length is:',l2cur if (idbg(550).eq.0) stop 'iniaddr' endif c l3cur is calculated as if all exchange potential functions were c to be kept in core. This is the case if iform=1 or iform=3. if (imethod.eq.1) then if (iform.eq.1.or.iform.eq.3) then if (l3cur.gt.length3) then write(iout6,*) 'cw_exch is too short!' write(iout6,*) 'declared length is:',length3 write(iout6,*) 'needed length is:',l3cur if (idbg(550).eq.0) stop 'iniaddr' endif else nons=0 do iorb=1,norb if (ll(iorb).ne.0) nons=nons+1 enddo c l3cur=nni*maxmu*(norb+nons) l3cur=mxsize*(norb+nons) if (l3cur.gt.length3) then write(iout6,*) 'cw_exch is too short!' write(iout6,*) 'declared length is:',length3 write(iout6,*) 'needed length is:',l3cur stop 'iniaddr' endif endif endif if (l5cur.gt.length5) then write(iout6,*) 'cw_sctch is too short!' write(iout6,*) 'declared length is:',length5 write(iout6,*) 'needed length is',10*mxsize8 c stop 'iniaddr' endif return end