read_fort.F90

Path: ILDGClass/OLDS/read_fort.F90
Last Update: Thu May 26 17:40:55 +0900 2011
dot/f_3.png

Required files

Methods

plaq   read_fort  

Included Modules

lattice_class file_tools_class field_gauge_class

Public Instance methods

Subroutine :
u(COL,COL,NDIM,NTX,NTY,NTZ,NTT) :complex(DP),intent(in)
plq :real(DP), intent(out)

[Source]

subroutine plaq(u,plq)
  use lattice_class
  implicit none
  !
  ! u_(a,b)_mu(x,y,z,t) = u(b,a,mu,x,y,z,t) ILDG format
  !
  complex(DP),intent(in)  :: u(COL,COL,NDIM,NTX,NTY,NTZ,NTT)
  real(DP),   intent(out) :: plq
  complex(DP) :: u12(COL,COL),u41(COL,COL),u1241(COL,COL)
  integer :: ic,jc,mu,nu
  integer :: it,iz,iy,ix
  integer :: it2,iz2,iy2,ix2
  integer :: it4,iz4,iy4,ix4

  plq = 0.0_DP
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,it,mu,nu,jc,ic, &
!$OMP&                    ix2,iy2,iz2,it2,   &
!$OMP&                    ix4,iy4,iz4,it4,   &
!$OMP&                    u12,u41,u1241) REDUCTION(+:plq)
  do it=1,NTT
  do iz=1,NTZ
  do iy=1,NTY
  do ix=1,NTX
    do mu=1,NDIM
    do nu=1,NDIM
      if (mu==nu) cycle

      ix2 = mod(ix-1 + isx(mu),NTX)+1
      iy2 = mod(iy-1 + isy(mu),NTY)+1
      iz2 = mod(iz-1 + isz(mu),NTZ)+1
      it2 = mod(it-1 + ist(mu),NTT)+1
      ix4 = mod(ix-1 + isx(nu),NTX)+1
      iy4 = mod(iy-1 + isy(nu),NTY)+1
      iz4 = mod(iz-1 + isz(nu),NTZ)+1
      it4 = mod(it-1 + ist(nu),NTT)+1

      !
      ! u12 = u(x1) * u(x2)
      !
      do jc=1,COL
      do ic=1,COL
        u12(ic,jc) = u(1,ic,mu,ix,iy,iz,it) * u(jc,1,nu,ix2,iy2,iz2,it2) + u(2,ic,mu,ix,iy,iz,it) * u(jc,2,nu,ix2,iy2,iz2,it2) + u(3,ic,mu,ix,iy,iz,it) * u(jc,3,nu,ix2,iy2,iz2,it2)
      enddo
      enddo

      !
      ! u41 = u(x4)^dag * u(x1)^dag
      !
      do jc=1,COL
      do ic=1,COL
        u41(ic,jc) = u(ic,1,mu,ix4,iy4,iz4,it4) * u(1,jc,nu,ix,iy,iz,it) + u(ic,2,mu,ix4,iy4,iz4,it4) * u(2,jc,nu,ix,iy,iz,it) + u(ic,3,mu,ix4,iy4,iz4,it4) * u(3,jc,nu,ix,iy,iz,it)
        u41(ic,jc) = conjg(u41(ic,jc))
      enddo
      enddo

      do jc=1,COL
      do ic=1,COL
        u1241(ic,jc) = u12(ic,1) * u41(1,jc) + u12(ic,2) * u41(2,jc) + u12(ic,3) * u41(3,jc)
      enddo
      enddo

      plq = plq + REAL(u1241(1,1) + u1241(2,2) + u1241(3,3))

    enddo
    enddo
  enddo
  enddo
  enddo
  enddo

  plq = plq/(NTX*NTY*NTZ*NTT*NDIM*(NDIM-1)*COL)

  return
end subroutine
Main Program :

[Source]

program read_fort
  use lattice_class
  use file_tools_class
  use field_gauge_class
  implicit none
  integer, parameter :: NARGS = 10

  !
  ! u_(a,b)_mu(x,y,z,t) = u(b,a,mu,x,y,z,t) ILDG format
  !
  complex(DP) :: u(COL,COL,NDIM,NTX,NTY,NTZ,NTT)

  type(vfield_gluon_wg) :: uu
  real(DP)    :: plq,Kud,Ks,CSW,beta,cRG0
  character(len=CHARLEN) :: fname,ftraj,fhead,fpath,fnode,cdate,ctime,czone,comment,arg(10)
  integer :: traj,i
  integer :: iout
  integer :: ic,jc,mu,iarg
  integer :: ipx,ipy,ipz
  integer :: ix,iy,iz,it
  integer :: itx,ity,itz,itt
  integer :: ith,ieo
  integer :: nodeid0

  iarg = COMMAND_ARGUMENT_COUNT()
  if (iarg /= NARGS) then
    write(*,'("Usage:")')
    write(*,'("% read_foort <input.bin> <comment> <beta> <crg0> <Kud> <Ks> <Csw> <traj> <output dir> <output head>")')
    stop
  endif
  do i=1,NARGS
    call GET_COMMAND_ARGUMENT(i,arg(i))
  enddo
  read(arg(1),'(A)')fname
  read(arg(2),'(A)')comment
  read(arg(3),*)beta
  read(arg(4),*)cRG0
  read(arg(5),*)Kud
  read(arg(6),*)Ks
  read(arg(7),*)Csw
  read(arg(8),*)traj
  read(arg(9), '(A)')fpath
  read(arg(10),'(A)')fhead

  write(*,'(A)')TRIM(fname)
  write(*,'(A)')TRIM(comment)
  write(*,'("Beta=",F16.7)')beta
  write(*,'("cRG0=",F16.7)')cRG0
  write(*,'(" Kud=",F16.7)')Kud
  write(*,'(" Ks_=",F16.7)')Ks
  write(*,'(" Csw=",F16.7)')Csw
  write(*,'("Traj=",I6)')traj
  write(*,'(A)')TRIM(fpath)
  write(*,'(A)')TRIM(fhead)

  call DATE_AND_TIME(cdate,ctime,czone)

  write(*,'("SIZEOF(U)=",I12)')SIZEOF(u)

  iout = search_free_file_unit()
  !
  ! Stream read from C-binary file
  !
  open(iout,file=fname,access="stream",form="unformatted")
  read(iout)u
  close(iout)

  call plaq(u,plq)

  write(*,'("PLQ=",F20.15)')plq

  call new(uu)
  call clear(uu)

  write(ftraj,'(I6.6)')traj
  fpath = TRIM(fpath)//"/"//TRIM(ftraj)

  do ipx=0,NDIMX-1
  do ipy=0,NDIMY-1
  do ipz=0,NDIMZ-1
    nodeid0 = ipx + ipy *NDIMX + ipz*NDIMX*NDIMY
    do mu=1,NDIM
!$OMP PARALLEL DO PRIVATE(ix,iy,iz,it,itx,ity,itz,itt,ith,ieo,jc,ic)
    do ix=1,NX
    do iy=1,NY
    do iz=1,NZ
    do it=1,NT
      itx = ix + ipx*NX
      ity = iy + ipy*NY
      itz = iz + ipz*NZ
      itt = it
      ith = it/2
      ieo = mod(itt+itz+ity+itx,2)
      do jc=1,COL
      do ic=1,COL
        uu%eo(ieo)%mu(mu)%s(ith,iz,iy,ix)%u(ic,jc) = u(jc,ic,mu,itx,ity,itz,itt)  ! ILDG ordering
      enddo
      enddo
    enddo
    enddo
    enddo
    enddo
    enddo

    write(fnode,'(".x",z1,"y",z1,"z",z1)')ipx,ipy,ipz
    fname=TRIM(fpath)//"/"//TRIM(fhead)//TRIM(fnode)
    write(*,'("Open:",A)')TRIM(fname)
    iout = search_free_file_unit()

    open(iout,file=fname,form="unformatted")
    write(iout) uu
    write(iout) NTX,NTY,NTZ,NTT
    write(iout) NDIMX,NDIMY,NDIMZ
    write(iout) plq
    close(iout)

  enddo
  enddo
  enddo

  stop
end program