!to cal energy potential LJ for xx,xy,xz

module energypotmodule
contains
subroutine energypot(N,rx,ry,rz,t,energy)

implicit none

integer, intent(IN) :: N,t  	!global parameters
double precision, dimension(:), intent(IN) :: rx,ry,rz   !global parameters
double precision, dimension(:), intent(INOUT) :: energy
double precision rsd,b,dx,dy,dz
integer p1,p2,e1,e2,n1,n2,i,j
double precision mp,mne,V,KE,rs,P


mp=1836.1526    !proton/electron mass ratio
mne=1838.684    !neutron/electron mass ratio

p1=7
p2=12
e1=1
e2=6
n1=13
n2=18

!energy(t)=0.0
P=0.0 

!!!!! Calculating Potential Energy  !!!!!!!!

!Electron-Proton interaction

do i=e1,e2        !i particle and calculating for this only
    
    do j=p1,p2        !j particle pairing with i for pairs (i,j)

!!! pair (p-e) !!!
!i=electron
!j=proton

            dx=rx(i)-rx(j)
            dy=ry(i)-ry(j)
            dz=rz(i)-rz(j)



            rs=(dx*dx)+(dy*dy)+(dz*dz)     !new rs will consider MIC rs's also that are qualifed for cutoff filtrations
            rs=sqrt(rs)
            !p=rs

            V=4.0*((1.0/(1.5*(rs**1.5)))+(8.0*rs))
      
            !energy(t)=energy(t)+V
            P=P+V   !storing total energy at each step in P

    end do
end do


!Electron-Electron interaction

do i=e1,e2-1        !i particle and calculating for this only
    do j=i+1,e2        !j particle pairing with i for pairs (i,j)

!!! pair (e-e) !!!
!i=electron
!j=electron

        !if (i/=j) then

            dx=rx(i)-rx(j)
            dy=ry(i)-ry(j)
            dz=rz(i)-rz(j)



            rs=(dx*dx)+(dy*dy)+(dz*dz)     !new rs will consider MIC rs's also that are qualifed for cutoff filtrations
            rs=sqrt(rs)
            !p=rs

            V=4.0*((-1*(rs**0.8)/0.8)+(1.5*rs))
      
            !energy(t)=energy(t)+V
            P=P+V

        !end if

    end do
end do



!proton-proton, neutron-neutron, proton-neutron interactions

do i=p1,n2-1        !i particle and calculating for this only
    do j=i+1,n2        !j particle pairing with i for pairs (i,j)

!!! pair (e-e) !!!
!i=electron
!j=electron

        !if (i/=j) then

            dx=rx(i)-rx(j)
            dy=ry(i)-ry(j)
            dz=rz(i)-rz(j)



            rs=(dx*dx)+(dy*dy)+(dz*dz)     !new rs will consider MIC rs's also that are qualifed for cutoff filtrations
            rs=sqrt(rs)
            !p=rs

            V=4.0*((-2.75*(rs**0.6))+(5.0*rs))
      
            !energy(t)=energy(t)+V
            P=P+V

        !end if

    end do
end do

energy(t)=P

end subroutine energypot
end module energypotmodule
	  
