From 05265f77a35120fba76d3dcf4d0222dd2c689d19 Mon Sep 17 00:00:00 2001 From: jean paul nshuti Date: Mon, 6 Oct 2025 13:48:56 +0200 Subject: [PATCH] Update the keys --- src/model/keys.f90 | 68 +++++++++--------- src/model/model.f90 | 164 +++++++++++++++++++------------------------- 2 files changed, 102 insertions(+), 130 deletions(-) diff --git a/src/model/keys.f90 b/src/model/keys.f90 index 128512d..ea9d2d1 100644 --- a/src/model/keys.f90 +++ b/src/model/keys.f90 @@ -10,7 +10,7 @@ contains character(len=1) prefix(4) parameter (prefix=['N','P','A','S']) !character (len=20) key(4,25) - integer,parameter:: np=34 + integer,parameter:: np=32 character(len=16) parname(np) integer i,j ! Defining keys for potential @@ -24,59 +24,57 @@ contains parname(4)='LXYVA2O2' parname(5)='LXYVE1O2' parname(6)='LXYVE2O2' - parname(7)='LXYVA2O3' - parname(8)='LXYVE1O3' - parname(9)='LXYVE2O3' + ! W & Z of E1 - parname(10)='LXYWZE1O0' - parname(11)='LXYWZE1O1' - parname(12)='LXYWZE1O2' - parname(13)='LXYWZE1O3' - parname(14)='LXYWZE2O0' - parname(15)='LXYWZE2O1' - parname(16)='LXYWZE2O2' - parname(17)='LXYWZE2O3' + parname(7)='LXYWZE1O0' + parname(8)='LXYWZE1O1' + parname(9)='LXYWZE1O2' + parname(10)='LXYWZE2O0' + parname(11)='LXYWZE2O1' + parname(12)='LXYWZE2O2' + ! WW and Z Pseudo between E1 and E2 ! p STANDS FOR PSEUDO JAHN-TELLER - parname(18)='LXYPE1E2O0' - parname(19)='LXYPE1E2O1' - parname(20)='LXYPE1E2O2' + parname(13)='LXYPE1E2O0' + parname(14)='LXYPE1E2O1' + parname(15)='LXYPE1E2O2' ! no order 3 ! PSEUDO A2 & E1 - parname(21)='LXYPA2E1O0' - parname(22)='LXYPA2E1O1' - parname(23)='LXYPA2E1O2' + parname(16)='LXYPA2E1O0' + parname(17)='LXYPA2E1O1' + parname(18)='LXYPA2E1O2' ! Pseudo JAHN-TELLER BETWEEN A2 AND E1 - parname(24)='LXYPA2E2O0' - parname(25)='LXYPA2E2O1' - parname(26)='LXYPA2E2O2' - parname(27)='LXYPA2E2O3' + parname(19)='LXYPA2E2O0' + parname(20)='LXYPA2E2O1' + parname(21)='LXYPA2E2O2' + ! keys for lz - parname(28)='LZWZE1O1' - parname(29)='LZWZE1O2' - parname(30)='LZWZE2O1' - parname(31)='LZWZE2O2' - parname(33)='LZPE1E2O0' - parname(34)='LZPE1E2O1' - parname(35)='LZPE1E2O2' - parname(36)='LZPA2E1O1' - parname(37)='LZPA2E2O2' - parname(39)='LZPA2E2O1' + parname(22)='LZWZE1O1' + parname(23)='LZWZE1O2' + parname(24)='LZWZE2O1' + parname(25)='LZWZE2O2' + parname(26)='LZPE1E2O0' + parname(27)='LZPE1E2O1' + parname(28)='LZPE1E2O2' + parname(29)='LZPA2E1O1' + parname(30)='LZPA2E2O2' + parname(31)='LZPA2E2O1' - parname(34)='TYPE_CAL'! TYPE OF THE CALCULATION WHETHER IT IS THE TRACE OR SOMETHING ELSE + parname(32)='TYPE_CAL'! TYPE OF THE CALCULATION WHETHER IT IS THE TRACE OR SOMETHING ELSE do i=1,np do j=1,4 - key(j, i)=prefix(j)//trim(parname(i))//':' ! first 86 keys are the potential keys - + key(j, i)=prefix(j)//trim(parname(i))//':' + write(8,*) key(j,i) enddo + write(8,*) '' enddo diff --git a/src/model/model.f90 b/src/model/model.f90 index 3b831da..5aeae8a 100644 --- a/src/model/model.f90 +++ b/src/model/model.f90 @@ -6,6 +6,8 @@ module diab_mod: use accuracy_constants, only: dp, idp use dim_parameter, only: ndiab, nstat, ntot,npar,qn,pst + private + public :: Lx_diab, Ly_diab, Lz_diab implicit none contains @@ -15,9 +17,11 @@ module diab_mod: real(dp),dimension(:),intent(in):: q,t real(dp),dimension(npar),intent(in):: p real(dp):: xs,ys,xb,yb,a,b - real(dp):: v3_vec(8) + real(dp):: v3_vec(8), v2(6) integer(idp):: i,j,id + + ! check the dimension of the matrix if (size(E,1) .ne. ndiab) then write(*,*) " Error in Lx_diab: wrong dimension of L matrix ", size(E,1) @@ -26,7 +30,24 @@ module diab_mod: ! rewrite the coordinate array q into symmetry adapted coordinates call rewrite_coord(q,a,xs,ys,xb,yb,b,1) - + v2(1)=xs**2-ys**2 + v2(2)=xb**2-yb**2 + v2(3)=xs*xb-ys*yb + v2(4)=2*xs*ys + v2(5)=2*xb*yb + v2(6)=xs*yb+xb*ys + + v3( 1) = xs*(xs**2-3*ys**2) + v3( 2) = xb*(xb**2-3*yb**2) + v3( 3) = xb*(xs**2-ys**2) - 2*yb*xs*ys + v3( 4) = xs*(xb**2-yb**2) - 2*ys*xb*yb + v3( 5) = ys*(3*xs**2-ys**2) + v3( 6) = yb*(3*xb**2-yb**2) + v3( 7) = yb*(xs**2-ys**2)+2*xb*xs*ys + v3( 8) = ys*(xb**2-yb**2)+2*xs*xb*yb + + e = 0.0_dp + id = 1 e(1,1)=e(1,1)+p(pst(1,id))*xs+p(pst(1,id)+1)*xb ! 2 param id=id+1 ! 2 e(2,2)=e(2,2)+p(pst(1,id))*xs+p(pst(1,id)+1)*xb ! 2 p @@ -50,43 +71,27 @@ module diab_mod: +p(pst(1,id)+2)*(xs*xb-ys*yb) e(5,5)=e(5,5)+p(pst(1,id))*(xs**2-ys**2)+p(pst(1,id)+1)*(xb**2-yb**2) & +p(pst(1,id)+2)*(xs*xb-ys*yb) - ! order 3 - id=id+1 ! 7 - e(1,1)=e(1,1)+p(pst(1,id))*xs*ss+p(pst(1,id)+1)*xb *sb ! 2 param - id=id+1 ! 8 - e(2,2)=e(2,2)+p(pst(1,id))*xs*ss+p(pst(1,id)+1)*xb *sb ! 2 p - e(3,3)=e(3,3)+p(pst(1,id))*xs*ss+p(pst(1,id)+1)*xb*sb - id =id+1 ! 9 - e(4,4)=e(4,4)+p(pst(1,id))*xs*ss+p(pst(1,id)+1)*xb*sb ! 2 p - e(5,5)=e(5,5)+p(pst(1,id))*xs*ss+p(pst(1,id)+1)*xb*sb - - ! W and Z term of E1 ! order 0 - id=id+1 ! 10 + id=id+1 ! 7 e(2,2)=e(2,2)+p(pst(1,id)) e(3,3)=e(3,3)-p(pst(1,id)) !e(2,3)=e(2,3) ! order 1 - id=id+1 ! 11 ! 2 param + id=id+1 ! 8 ! 2 param e(2,2)=e(2,2)+ p(pst(1,id))*xs+p(pst(1,id)+1)*xb e(3,3)=e(3,3)- (p(pst(1,id))*xs+p(pst(1,id)+1)*xb) e(2,3)=e(2,3)- p(pst(1,id))*ys -p(pst(1,id)+1)*yb ! order 2 - id=id+1 ! 12 ! 3p + id=id+1 ! 9 ! 3p do i=1,3 e(2,2)=e(2,2)+p(pst(1,id)+(i-1))*v2(i) e(3,3)=e(3,3)-p(pst(1,id)+(i-1))*v2(i) e(2,3)=e(2,3)+ p(pst(1,id)+(i-1))*v2(i+3) enddo ! order 3 - id=id+1 ! 13 ! 8 param - do i=1,4 - e(2,2)=e(2,2)+(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - e(3,3)=e(3,3)-(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - e(2,3)=e(2,3)+(-p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i+4) - enddo + ! try the testing of higher order terms !e(2,3)=e(2,3)- p(pst(1,id))*ys*ss +p(pst(1,id)+1)*ss*2*xs*ys @@ -97,32 +102,24 @@ module diab_mod: ! W and Z for E2 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - id=id+1 ! 14 + id=id+1 ! 10 e(4,4)=e(4,4)+p(pst(1,id)) e(5,5)=e(5,5)-p(pst(1,id)) e(4,5)=e(4,5) ! order 1 - id=id+1 ! 2 param 15 + id=id+1 !112 param 15 e(4,4)=e(4,4)+ p(pst(1,id))*xs+p(pst(1,id)+1)*xb e(5,5)=e(5,5)- (p(pst(1,id))*xs+p(pst(1,id)+1)*xb) e(4,5)=e(4,5)- p(pst(1,id))*ys-p(pst(1,id)+1)*yb ! order 2 - id=id+1 ! 16 ! 3p + id=id+1 ! 12 ! 3p do i=1,3 e(4,4)=e(4,4)+p(pst(1,id)+(i-1))*v2(i) e(5,5)=e(5,5)-p(pst(1,id)+(i-1))*v2(i) e(4,5)=e(4,5)+ p(pst(1,id)+(i-1))*v2(i+3) enddo - ! order 3 - id=id+1 ! 17 ! 8 param - do i=1,4 - e(4,4)=e(4,4)+(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - e(5,5)=e(5,5)-(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - e(4,5)=e(4,5)+(-p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i+4) - enddo - ! make the dipole E = b* E e = b * e @@ -131,20 +128,20 @@ module diab_mod: ! WW and ZZ - id =id+1 ! 18 + id =id+1 ! 13 e(2,4)=e(2,4)+p(pst(1,id))*b e(3,5)=e(3,5)-p(pst(1,id))*b ! ORDER 1 - id=id+1 ! 19 ! 6 parama + id=id+1 ! 14 ! 6 parama e(2,4)=e(2,4)+b*((p(pst(1,id))+p(pst(1,id)+1)+p(pst(1,id)+2))*xs+(p(pst(1,id)+3)+p(pst(1,id)+4)+p(pst(1,id)+5))*xb) e(3,5)=e(3,5)+b*((p(pst(1,id))+p(pst(1,id)+1)-p(pst(1,id)+2))*xs+(p(pst(1,id)+3)+p(pst(1,id)+4)-p(pst(1,id)+5))*xb) e(2,5)=e(2,5)+b*((p(pst(1,id))-p(pst(1,id)+1)-p(pst(1,id)+2))*ys+(p(pst(1,id)+3)-p(pst(1,id)+4)-p(pst(1,id)+5))*yb) e(3,4)=e(3,4)+b*((-p(pst(1,id))+p(pst(1,id)+1)-p(pst(1,id)+2))*ys+(-p(pst(1,id)+3)+p(pst(1,id)+4)-p(pst(1,id)+5))*yb) ! order 2 - id=id+1 ! 20 + id=id+1 ! 15 - do i=1,3 ! 9 param + do i=1,3 ! param e(2,4)=e(2,4)+b*(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+2))+p(pst(1,id)+(i+5)))*v2(i) e(3,5)=e(3,5)+b*(-p(pst(1,id)+(i-1))+p(pst(1,id)+(i+2))+p(pst(1,id)+(i+5)))*v2(i) e(2,5)=e(2,5)+b*(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+2))-p(pst(1,id)+(i+5)))*v2(i+3) @@ -155,17 +152,17 @@ module diab_mod: ! ################################################## !################################################### ! order 0 - id=id+1 ! 1 param ! 21 + id=id+1 ! 1 param ! 16 e(1,3)=e(1,3)+b*(p(pst(1,id))) ! order 1 - id = id +1 ! 22 + id = id +1 ! 17 e(1,2)=e(1,2)-b*(p(pst(1,id))*ys + p(pst(1,id)+1)*yb) e(1,3)=e(1,3)+b*(p(pst(1,id))*xs + p(pst(1,id)+1)*xb) ! order 2 - id=id+1 ! 23 + id=id+1 ! 18 e(1,2)=e(1,2)+b*(p(pst(1,id))*(2*xs*ys)+p(pst(1,id)+1)*(2*xb*yb)+p(pst(1,id)+2)*(xs*yb+xb*ys)) e(1,3)=e(1,3)+b*(p(pst(1,id))*(xs**2-ys**2) + p(pst(1,id)+1)*(xb**2-yb**2) & @@ -176,28 +173,21 @@ module diab_mod: !########################################################################################################## ! order 0 - id =id+1 !24 + id =id+1 !19 e(1,5)=e(1,5)+p(pst(1,id)) ! order 1 - id = id +1 ! 25 + id = id +1 ! 20 e(1,4)=e(1,4)-(p(pst(1,id))*ys + p(pst(1,id)+1)*yb) e(1,5)=e(1,5)+(p(pst(1,id))*xs + p(pst(1,id)+1)*xb) ! order 2 - id=id+1 ! 26 + id=id+1 ! 21 e(1,4)=e(1,4)+p(pst(1,id))*(2*xs*ys)+p(pst(1,id)+1)*(2*xb*yb)+p(pst(1,id)+2)*(xs*yb+xb*ys) e(1,5)=e(1,5)+p(pst(1,id))*(xs**2-ys**2) + p(pst(1,id)+1)*(xb**2-yb**2) & +p(pst(1,id)+2)*(xs*xb-ys*yb) - ! order 3 - id=id+1 ! 27 ! 8 param - do i=1,4 - e(1,4)=e(1,4)+(p(pst(1,id)+(i-1))-p(pst(1,id)+(i+3)))*v3(i+4) - e(1,5)=e(1,5)+(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - - enddo e(1,4:5) = b* e(1,4:5) @@ -222,9 +212,12 @@ module diab_mod: ! rewrite the coordinate array q into symmetry adapted coordinates call rewrite_coord(q,a,xs,ys,xb,yb,b,1) - e=0.0d0 - ss=xs**2+ys**2 ! totaly symmetric term - sb=xb**2+yb**2 + v2(1)=xs**2-ys**2 + v2(2)=xb**2-yb**2 + v2(3)=xs*xb-ys*yb + v2(4)=2*xs*ys + v2(5)=2*xb*yb + v2(6)=xs*yb+xb*ys v3( 1) = xs*(xs**2-3*ys**2) v3( 2) = xb*(xb**2-3*yb**2) @@ -233,9 +226,10 @@ module diab_mod: v3( 5) = ys*(3*xs**2-ys**2) v3( 6) = yb*(3*xb**2-yb**2) v3( 7) = yb*(xs**2-ys**2)+2*xb*xs*ys - v3( 8) = ys*(xb**2-yb**2)+2*xs*xb*yb + v3( 8) = ys*(xb**2-yb**2)+2*xs*xb*yb - e=0.0d0 + e = 0.0_dp + ! V-term id=1 ! 1 @@ -256,87 +250,65 @@ module diab_mod: id=id+1 ! 6 e(4,4)=e(4,4)-(p(pst(1,id))*(2*xs*ys)+p(pst(1,id)+1)*(2*xb*yb)+p(pst(1,id)+2)*(xs*yb+xb*ys)) e(5,5)=e(5,5)-(p(pst(1,id))*(2*xs*ys)+p(pst(1,id)+1)*(2*xb*yb)+p(pst(1,id)+2)*(xs*yb+xb*ys)) - ! order 3 - id=id+1 ! 7 - - e(1,1)=e(1,1)+p(pst(1,id))*ys*ss + p(pst(1,id)+1)*yb*sb - id=id+1 ! 2 - e(2,2)=e(2,2)+p(pst(1,id))*ys*ss + p(pst(1,id)+1)*yb*sb - e(3,3)=e(3,3)+p(pst(1,id))*ys*ss + p(pst(1,id)+1)*yb*sb - id =id+1 ! 3 - e(4,4)=e(4,4)+p(pst(1,id))*ys*ss + p(pst(1,id)+1)*yb*sb - e(5,5)=e(5,5)+p(pst(1,id))*ys*ss + p(pst(1,id)+1)*yb*sb + ! W and Z of E1 ! order 0 - id=id+1 ! 10 + id=id+1 ! 7 e(2,3)=e(2,3)+p(pst(1,id)) ! order 1 - id=id+1 ! + id=id+1 ! 8 e(2,2)=e(2,2)-p(pst(1,id))*ys -p(pst(1,id)+1)*yb e(3,3)=e(3,3)+p(pst(1,id))*ys+ p(pst(1,id)+1)*yb e(2,3)=e(2,3)-p(pst(1,id))*xs -p(pst(1,id)+1)*xb ! order 2 - id=id+1 ! 12 + id=id+1 ! 9 do i=1,3 e(2,2)=e(2,2)+p(pst(1,id)+(i-1))*v2(i+3) e(3,3)=e(3,3)-p(pst(1,id)+(i-1))*v2(i+3) e(2,3)=e(2,3)-p(pst(1,id)+(i-1))*v2(i) enddo - id=id+1 ! 8 - do i=1,4 - e(2,2)=e(2,2)+(p(pst(1,id)+(i-1))-p(pst(1,id)+(i+3)))*v3(i+4) - e(3,3)=e(3,3)-(p(pst(1,id)+(i-1))-p(pst(1,id)+(i+3)))*v3(i+4) - e(2,3)=e(2,3)+(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - enddo - !! W and Z of E2 !!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!! ! order 0 - id=id+1 ! 14 + id=id+1 ! 10 e(4,5)=e(4,5)+p(pst(1,id)) ! order 1 - id=id+1 ! 15 + id=id+1 ! 11 e(4,4)=e(4,4)-p(pst(1,id))*ys -p(pst(1,id)+1)*yb e(5,5)=e(5,5)+p(pst(1,id))*ys+ p(pst(1,id)+1)*yb e(4,5)=e(4,5)-p(pst(1,id))*xs -p(pst(1,id)+1)*xb ! order 2 - id=id+1 ! 16 + id=id+1 ! 12 do i=1,3 e(4,4)=e(4,4)+p(pst(1,id)+(i-1))*v2(i+3) e(5,5)=e(5,5)-p(pst(1,id)+(i-1))*v2(i+3) e(4,5)=e(4,5)-p(pst(1,id)+(i-1))*v2(i) enddo - id=id+1 ! 17 - do i=1,4 - e(4,4)=e(4,4)+(p(pst(1,id)+(i-1))-p(pst(1,id)+(i+3)))*v3(i+4) - e(5,5)=e(5,5)-(p(pst(1,id)+(i-1))-p(pst(1,id)+(i+3)))*v3(i+4) - e(4,5)=e(4,5)+(p(pst(1,id)+(i-1))+p(pst(1,id)+(i+3)))*v3(i) - enddo ! PSEUDO JAHN-TELLER E1 AND E2 e = b* e !ORDER 0 - id=id+1 ! 18 + id=id+1 ! 13 e(2,5)=e(2,5)+p(pst(1,id)) e(3,4)=e(3,4)+p(pst(1,id)) ! order 1 - id=id+1 + id=id+1 ! 14 e(2,4)=e(2,4)+((p(pst(1,id))+p(pst(1,id)+1)-p(pst(1,id)+2))*ys+(p(pst(1,id)+3)+p(pst(1,id)+4)-p(pst(1,id)+5))*yb) e(3,5)=e(3,5)+((p(pst(1,id))+p(pst(1,id)+1)+p(pst(1,id)+2))*ys+(p(pst(1,id)+3)+p(pst(1,id)+4)+p(pst(1,id)+5))*yb) e(2,5)=e(2,5)+((-p(pst(1,id))+p(pst(1,id)+1)-p(pst(1,id)+2))*xs+(-p(pst(1,id)+3)+p(pst(1,id)+4)-p(pst(1,id)+5))*xb) e(3,4)=e(3,4)+((p(pst(1,id))-p(pst(1,id)+1)-p(pst(1,id)+2))*xs+(+p(pst(1,id)+3)-p(pst(1,id)+4)-p(pst(1,id)+5))*xb) ! order 2 - id=id+1 + id=id+1 ! 15 e(2,4)=e(2,4)+(p(pst(1,id)+(i-1))-p(pst(1,id)+(i+2))-p(pst(1,id)+(i+5)))*v2(i+3) e(3,5)=e(3,5)+(-p(pst(1,id)+(i-1))-p(pst(1,id)+(i+2))-p(pst(1,id)+(i+5)))*v2(i+3) @@ -350,16 +322,16 @@ module diab_mod: ! ##################### ! order 0 - id=id+1 + id=id+1 ! 16 e(1,2)=e(1,2)+(p(pst(1,id))) ! order 1 - id=id+1 + id=id+1 ! 17 e(1,2)=e(1,2)-(p(pst(1,id))*xs + p(pst(1,id)+1)*xb) e(1,3)=e(1,3)-(p(pst(1,id))*ys + p(pst(1,id)+1)*yb) ! order 2 - id=id+1 + id=id+1 !18 e(1,2)=e(1,2)-(p(pst(1,id))*(xs**2-ys**2)+p(pst(1,id)+1)*(xb**2-yb**2) & +p(pst(1,id)+2)*(xs*xb-ys*yb)) e(1,3)=e(1,3)+(p(pst(1,id))*(2*xs*ys)+p(pst(1,id)+1)*(2*xb*yb)+ & @@ -370,16 +342,16 @@ module diab_mod: !############################################################################### ! order 0 - id = id+1 + id = id+1 !19 e(1,4)=e(1,4)+p(pst(1,id)) ! order 1 - id=id+1 + id=id+1 ! 20 e(1,4)=e(1,4)-(p(pst(1,id))*xs + p(pst(1,id)+1)*xb) e(1,5)=e(1,5)-(p(pst(1,id))*ys + p(pst(1,id)+1)*yb) ! order 2 - id=id+1 + id=id+1 ! 21 e(1,4)=e(1,4)-(p(pst(1,id))*(xs**2-ys**2)+p(pst(1,id)+1)*(xb**2-yb**2) & +p(pst(1,id)+2)*(xs*xb-ys*yb)) e(1,5)=e(1,5)+(p(pst(1,id))*(2*xs*ys)+p(pst(1,id)+1)*(2*xb*yb)+ & @@ -408,9 +380,11 @@ module diab_mod: stop endif call rewrite_coord(q,a,xs,xb,yb,b,1) + + e = 0.0_dp ! id for lz - id = ! has to be + id = 22 ! has to be ! the diagonal terms ! the v-term is 0th order and 3rd order. @@ -418,7 +392,7 @@ module diab_mod: ! w and z of E'' ! order 1 - id =id +1 + id = id e(2,2) = e(2,2) + p(pst(1,id))*ys + p(pst(1,id)+1)*yb e(3,3) = e(3,3) - p(pst(1,id))*ys - p(pst(1,id)+1)*yb e(2,3) = e(2,3) - p(pst(1,id))*xs -p(pst(1,id)+1)*xb