Class | ildg_class |
In: |
ILDGClass/ildg_class.F90
|
$Id: ildg_class.F90,v 1.6 2011/05/27 04:35:47 ishikawa Exp $
Subroutine : | |
fpath : | character(len=*), intent(in) |
fhead : | character(len=*), intent(in) |
u : | type(gauge_link_ildg), intent(in) |
plq : | real(DP), intent(in) |
Convert ILDG config to desired format
Splitting for parallel
Reordering for even/odd
fpath : directory path fhead : config file name u : ILDG ordering config plq : plaquette value
When
fpath = "hogehoge/ugauga/config" fhead = "myconfig"
are given, config files are stored as
hogehoge/ugauga/config/myconfig.x?y?z?
subroutine convert_and_save_config(fpath,fhead,u,plq) ! ! Convert ILDG config to desired format ! ! Splitting for parallel ! Reordering for even/odd ! ! fpath : directory path ! fhead : config file name ! u : ILDG ordering config ! plq : plaquette value ! ! When ! fpath = "hogehoge/ugauga/config" ! fhead = "myconfig" ! are given, config files are stored as ! ! hogehoge/ugauga/config/myconfig.x?y?z? ! use file_tools_class use field_gauge_class implicit none character(len=*), intent(in) :: fpath,fhead type(gauge_link_ildg), intent(in) :: u real(DP), intent(in) :: plq character(len=CHARLEN) :: fnode,fname type(vfield_gluon_wg), asynchronous, allocatable :: uu integer :: io_list(0:NDIMX-1,0:NDIMY-1,0:NDIMZ-1) integer :: iout integer :: ipx,ipy,ipz integer :: mu,ix,iy,iz,it,ith,ic,jc integer :: itx,ity,itz,itt,ieo type(timer) :: etime !----------------------------- ! open parallel config files !----------------------------- do ipz=0,NDIMZ-1 do ipy=0,NDIMY-1 do ipx=0,NDIMX-1 io_list(ipx,ipy,ipz) = search_free_file_unit() iout = io_list(ipx,ipy,ipz) write(fnode,'(".x",z1,"y",z1,"z",z1)')ipx,ipy,ipz fname=TRIM(fpath)//"/"//TRIM(fhead)//TRIM(fnode) write(*,'("Open:",A)')TRIM(fname) open(unit=iout,file=fname,form='unformatted',status='unknown',asynchronous='yes') enddo enddo enddo write(*,'("Converting ILDG ordering to our ordering....")',advance='no') call new(etime,master_only=.true.) call tic(etime) allocate(uu) call new(uu) call clear(uu) !----------------------------- ! splitting, reordering ! and storing !----------------------------- do ipz=0,NDIMZ-1 do ipy=0,NDIMY-1 do ipx=0,NDIMX-1 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%u(jc,ic,mu,itx,ity,itz,itt) enddo enddo enddo enddo enddo enddo enddo iout = io_list(ipx,ipy,ipz) write(unit=iout,asynchronous='yes') uu write(unit=iout,asynchronous='yes') NTX,NTY,NTZ,NTT write(unit=iout,asynchronous='yes') NDIMX,NDIMY,NDIMZ write(unit=iout,asynchronous='yes') plq enddo enddo enddo deallocate(uu) call toc(etime) write(*,'(" done. ETIME:",F10.4)')get_elapse(etime) call delete(etime) write(*,'("Waiting for Storing data....")',advance='no') call new(etime,master_only=.true.) call tic(etime) do ipz=0,NDIMZ-1 do ipy=0,NDIMY-1 do ipx=0,NDIMX-1 iout = io_list(ipx,ipy,ipz) close(unit=iout) enddo enddo enddo call toc(etime) write(*,'(" done. ETIME:",F10.4)')get_elapse(etime) call delete(etime) return end subroutine
Subroutine : | |
u : | type(gauge_link_ildg), intent(inout) |
subroutine delete_ildg_config(u) implicit none type(gauge_link_ildg), intent(inout) :: u deallocate(u%u) u%sizeof = COL*COL*NDIM*NTX*NTY*NTZ*NTT*SIZE_OF_COMPLEX return end subroutine
Derived Type : | |
u(:,:,:,:,:,:,:) : | complex(DP), allocatable |
sizeof = COL*COL*NDIM*NTX*NTY*NTZ*NTT*SIZE_OF_COMPLEX : | integer |
ILDG format ordering
u_(a,b)_mu(x,y,z,t) = u(b,a,mu,x,y,z,t) ILDG format
Subroutine : | |
u : | type(gauge_link_ildg), intent(inout) |
subroutine new_ildg_config(u) implicit none type(gauge_link_ildg), intent(inout) :: u allocate(u%u(COL,COL,NDIM,NTX,NTY,NTZ,NTT)) u%sizeof = COL*COL*NDIM*NTX*NTY*NTZ*NTT*SIZE_OF_COMPLEX return end subroutine
Function : | |
plq : | real(DP) |
u : | type(gauge_link_ildg), intent(in) |
Compute plaquette value
plq = Sum_{mu,nu,n} ReTr[U_mu,nu(n)]/(COL*NTT*NTZ*NTY*NTX*(NDIM-1)*NDIM)
function plaquette(u) result(plq) ! ! Compute plaquette value ! ! plq = Sum_{mu,nu,n} ReTr[U_mu,nu(n)]/(COL*NTT*NTZ*NTY*NTX*(NDIM-1)*NDIM) ! ! implicit none type(gauge_link_ildg), intent(in) :: u real(DP) :: 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%u(1,ic,mu,ix,iy,iz,it) * u%u(jc,1,nu,ix2,iy2,iz2,it2) + u%u(2,ic,mu,ix,iy,iz,it) * u%u(jc,2,nu,ix2,iy2,iz2,it2) + u%u(3,ic,mu,ix,iy,iz,it) * u%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%u(ic,1,mu,ix4,iy4,iz4,it4) * u%u(1,jc,nu,ix,iy,iz,it) + u%u(ic,2,mu,ix4,iy4,iz4,it4) * u%u(2,jc,nu,ix,iy,iz,it) + u%u(ic,3,mu,ix4,iy4,iz4,it4) * u%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)) + REAL(u1241(2,2)) + REAL(u1241(3,3)) enddo enddo enddo enddo enddo enddo plq = plq/(NTX*NTY*NTZ*NTT*NDIM*(NDIM-1)*COL) return end function
Subroutine : | |
fname : | character(len=*), intent(in) |
u : | type(gauge_link_ildg), intent(inout) |
endian_conversion : | logical, optional, intent(in) |
Read ILDG config binary data
This uses stream read from C-binary file
subroutine read_ildg_config(fname,u,endian_conversion) ! ! Read ILDG config binary data ! ! This uses stream read from C-binary file ! use file_tools_class implicit none character(len=*), intent(in) :: fname type(gauge_link_ildg), intent(inout) :: u logical, optional, intent(in) :: endian_conversion integer :: iout integer :: itx,ity,itz,itt,mu,ic,jc real(DP) :: ure,uim complex(DP) :: uu type(timer) :: etime character(len=1) :: data(1:8) write(*,'("Reading ILDG binary config data....")',advance='no') call new(etime,master_only=.true.) call tic(etime) iout = search_free_file_unit() open(iout,file=TRIM(ADJUSTL(fname)),access="stream",form="unformatted",status="old") read(iout)u%u close(iout) call toc(etime) write(*,'(" done. ETIME:",F10.4)')get_elapse(etime) call delete(etime) if (present(endian_conversion)) then if (endian_conversion) then write(*,'("Endian converting....")',advance='no') call new(etime,master_only=.true.) call tic(etime) ! ! Do endian conversion ! !$OMP PARALLEL DO PRIVATE(itt,itz,ity,itx,mu,ic,jc,uu,ure,uim,data) do itt=1,NTT do itz=1,NTZ do ity=1,NTY do itx=1,NTX do mu=1,NDIM do ic=1,COL do jc=1,COL uu = u%u(jc,ic,mu,itx,ity,itz,itt) ure = REAL(uu) uim = AIMAG(uu) data(:) = TRANSFER(ure,data) ! convert to byte data data(1:8) = data(8:1:-1) ! endian conversion ure = TRANSFER(data,ure) ! convert to double data data(:) = TRANSFER(uim,data) ! convert to byte data data(1:8) = data(8:1:-1) ! endian conversion uim = TRANSFER(data,uim) ! convert to double data u%u(jc,ic,mu,itx,ity,itz,itt) = CMPLX(ure,uim,kind=DP) enddo enddo enddo enddo enddo enddo enddo call toc(etime) write(*,'(" done. ETIME:",F10.4)')get_elapse(etime) call delete(etime) endif endif return end subroutine