Commit a62bce59 authored by Jonas Schwab's avatar Jonas Schwab

Update Tutorial/Solutions/Exercise_2/ALF/Analysis/

parent a4b80ffb
.PHONY : all tidy clean
BINS= cov_scal.out cov_tau.out cov_tau_ph.out Max_SAC.out cov_eq.out
OBJS= cov_scal.o cov_tau.o cov_tau_ph.o Max_SAC.o cov_eq.o
OBJS1 = Predefined_Latt_mod.o
MODS=predefined_lattices.mod
BINS=cov_eq.out cov_scal.out cov_tau.out cov_tau_ph.out Max_SAC.out ana.out #repack_latt.out
OBJS=cov_eq.o cov_scal.o cov_tau.o cov_tau_ph.o Max_SAC.o ana.o
OBJS1=Predefined_Latt_mod.o ana_mod.o
MODS=predefined_lattices.mod ana_mod.mod
all: $(OBJS1) $(BINS)
cov_tau_ph.F90 : cov_tau.F90
cpp -w -P -DPartHole cov_tau.F90 -o cov_tau_ph.F90
all: $(BINS)
Predefined_Latt_mod.o: ../Prog/Predefined_Latt_mod.F90
$(ALF_FC) -c -o $@ $(ALF_FLAGS_ANA) $<
cov_tau_ph.o: cov_tau.F90
$(ALF_FC) -c -o $@ $(ALF_FLAGS_ANA) -DPartHole $<
%.o: %.F90
$(ALF_FC) -c $*.F90 $(ALF_FLAGS_ANA)
$(ALF_FC) -c -o $@ $(ALF_FLAGS_ANA) $<
%.out: %.F90
$(ALF_FC) -c $(ALF_FLAGS_ANA) $<
$(ALF_FC) -o $@ $*.o $(OBJS1) $(ALF_LIB)
%.out: $(OBJS1) %.o
$(ALF_FC) -o $@ $^ $(ALF_LIB)
tidy:
rm -f $(OBJS) $(OBJS1) $(MODS) cov_tau_ph.F90
......
! Copyright (C) 2016-2019 The ALF project
! Copyright (C) 2016-2020 The ALF project
!
! The ALF project is free software: you can redistribute it and/or modify
! it under the terms of the GNU General Public License as published by
......@@ -42,9 +42,8 @@
!--------------------------------------------------------------------
Use MaxEnt_stoch_mod
use iso_fortran_env, only: output_unit, error_unit
Implicit None
!Implicit Real (Kind=Kind(0.d0)) (A-G,O-Z)
!Implicit Integer (H-N)
Interface
Subroutine Rescale ( XCOV, XQMC,XTAU, Ntau_st, Ntau_en, Tolerance, Ntau)
......@@ -67,7 +66,7 @@
Integer :: Ngamma, Ndis, NBins, NSweeps, Nwarm, N_alpha, N_cov
Integer :: N_skip, N_rebin, Norb
Integer :: N_skip, N_rebin, N_Back, N_auto, Norb
Real (Kind=Kind(0.d0)) :: OM_st, OM_en, alpha_st, R, Tolerance
Logical :: Checkpoint
Character (Len=2) :: Channel
......@@ -78,9 +77,9 @@
Real (Kind=Kind(0.d0)) :: Zero
NAMELIST /VAR_Max_Stoch/ Ngamma, Ndis, NBins, NSweeps, Nwarm, N_alpha, &
& OM_st, OM_en, alpha_st, R, Checkpoint, Channel, Tolerance
& OM_st, OM_en, alpha_st, R, Checkpoint, Tolerance
NAMELIST /VAR_errors/ N_skip, N_rebin, N_cov
NAMELIST /VAR_errors/ N_skip, N_rebin, N_cov, N_Back, N_auto
open(unit=30,file='parameters',status='old',action='read', iostat=io_error)
if (io_error.eq.0) then
......@@ -92,26 +91,9 @@
endif
close(30)
Open(unit=50,File='Info_MaxEnt',Status="unknown")
write(50,*) 'Channel :: ', Channel
If (Channel == "ph" ) then
Write(50,*) 'Om_start is set to zero. PH channel corresponds to symmetric data '
Om_st = 0.d0
endif
Write(50, "('Covariance :: ',I2)") N_cov
Write(50, "('Checkpoint :: ',L )") Checkpoint
Write(50, "('Om_st, Om_en :: ',2x,F12.6,2x,F12.6)") Om_st, Om_en
Write(50, "('Delta Om :: ',2x,F12.6)") (Om_en - Om_st)/real(Ndis,kind(0.d0))
Write(50, "('Bins, Sweeps, Warm :: ',2x,I4,2x,I4,2x,I4)") NBins, NSweeps, Nwarm
If (N_alpha <= 10 ) then
Write(error_unit,*) 'Not enough temperatures: N_alpha has to be bigger than 10'
error stop 1
Endif
Write(50, "('N_Alpha, Alpha_st,R:: ',2x,I4,F12.6,2x,F12.6)") N_alpha, alpha_st, R
open (unit=10,File="g_dat", status="unknown")
read(10,*) ntau, nbin_qmc, Beta, Norb
read(10,*) ntau, nbin_qmc, Beta, Norb, Channel
Allocate ( XCOV(NTAU,NTAU), XQMC(NTAU),XTAU(NTAU) )
XCOV = 0.d0
Do nt = 1,NTAU
......@@ -128,6 +110,23 @@
close(10)
dtau = Xtau(2) - Xtau(1)
Open(unit=50,File='Info_MaxEnt',Status="unknown")
write(50,*) 'Channel :: ', Channel
If (Channel == "ph" ) then
Write(50,*) 'Om_start is set to zero. PH channel corresponds to symmetric data '
Om_st = 0.d0
endif
Write(50, "('Covariance :: ',I2)") N_cov
Write(50, "('Checkpoint :: ',L1)") Checkpoint
Write(50, "('Om_st, Om_en :: ',2x,F12.6,2x,F12.6)") Om_st, Om_en
Write(50, "('Delta Om :: ',2x,F12.6)") (Om_en - Om_st)/real(Ndis,kind(0.d0))
Write(50, "('Bins, Sweeps, Warm :: ',2x,I4,2x,I4,2x,I4)") NBins, NSweeps, Nwarm
If (N_alpha <= 10 ) then
Write(error_unit,*) 'Not enough temperatures: N_alpha has to be bigger than 10'
error stop 1
Endif
Write(50, "('N_Alpha, Alpha_st,R:: ',2x,I4,F12.6,2x,F12.6)") N_alpha, alpha_st, R
Zero = 1.D-10
pi = acos(-1.d0)
Ntau_st = 1
......@@ -138,7 +137,7 @@
Case ("PP")
xmom1 = 2.d0* pi * xqmc(1)
Case ("P")
xmom1 = pi * real(norb,Kind(0.d0))
xmom1 = pi * ( xqmc(1) + xqmc(ntau) )
! Remove the tau = beta point from the data since it is correlated
! due to the sum rule with the tau=0 data point. Also if the tau = 0
! data point has no fluctations (due to particle-hole symmetry for instance)
......@@ -222,9 +221,9 @@
If ( .not. Checkpoint ) then
Command = "rm dump*"
Call System (Command)
Call EXECUTE_COMMAND_LINE(Command)
Command = "ls"
Call System (Command)
Call EXECUTE_COMMAND_LINE(Command)
endif
Open (Unit=10,File="energies",status="unknown")
......@@ -295,7 +294,7 @@
Z = Z + A(nwp)/cmplx( om - omp, delta, kind(0.d0))
enddo
Z = Z * dom
write(43,"('X'2x,F14.7,2x,F16.8,2x,F16.8)" ) xom(nw), dble(Z), -Aimag(Z)/pi
write(43,"('X',2x,F14.7,2x,F16.8,2x,F16.8)") xom(nw), dble(Z), -Aimag(Z)/pi
enddo
close(43)
......
../Prog/Predefined_Latt_mod.F90
\ No newline at end of file
! Copyright (C) 2019 The ALF project
!
! The ALF project is free software: you can redistribute it and/or modify
! it under the terms of the GNU General Public License as published by
! the Free Software Foundation, either version 3 of the License, or
! (at your option) any later version.
!
! The ALF project is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
! GNU General Public License for more details.
!
! You should have received a copy of the GNU General Public License
! along with Foobar. If not, see http://www.gnu.org/licenses/.
!
! Under Section 7 of GPL version 3 we require you to fulfill the following additional terms:
!
! - It is our hope that this program makes a contribution to the scientific community. Being
! part of that community we feel that it is reasonable to require you to give an attribution
! back to the original authors if you have benefitted from this program.
! Guidelines for a proper citation can be found on the project's homepage
! http://alf.physik.uni-wuerzburg.de .
!
! - We require the preservation of the above copyright notice and this license in all original files.
!
! - We prohibit the misrepresentation of the origin of the original source files. To obtain
! the original source files please visit the homepage http://alf.physik.uni-wuerzburg.de .
!
! - If you make substantial changes to the program we require you to either consider contributing
! to the ALF project or to mark your material in a reasonable way as different from the original version.
Program ana
!--------------------------------------------------------------------
!> @author
!> ALF-project
!
!> @brief
!> Wrapper program for analyzing Monte Carlo bins.
!>
!> @details
!> Analyzes scalar observables, equal-time and timedisplaced correlation functions.
!>
!> Hand over any number of file names as command line argumens to this program.
!> If will analyze names ending in "_scal" as scalar observables, "_eq" asequal-time
!> correlation functions and "_tau" as timedisplaced correlation functions.
!>
!> For timedisplaced observables that are listen in "names_PH" in namelist "VAR_PH_names"
!> in file "parameters", it is assumed that they are symmetric in imaginary time around tau = beta/2
!>
!
!--------------------------------------------------------------------
use ana_mod
implicit none
Integer :: i, n, nargs
Character (len=64) :: name
Character (len=64), allocatable :: names(:)
nargs = COMMAND_ARGUMENT_COUNT()
allocate( names(nargs) )
do i = 1, nargs
CALL GET_COMMAND_ARGUMENT(i, name)
names(i) = name
enddo
do n=1, size(names)
name = names(n)
i = len(trim(name)) -4
if ( i > 1 ) then
if ( name(i:) == '_scal' ) then
print *, ''
print '(A,A)', "analyzing scalar observable ", name
call Cov_vec(name)
endif
endif
enddo
do n=1, size(names)
name = names(n)
i = len(trim(name)) -2
if ( i > 1 ) then
if ( name(i:) == '_eq' ) then
print *, ''
print '(A,A)', "analyzing equal time correlation ", name
call Cov_eq(name)
endif
endif
enddo
do n=1, size(names)
name = names(n)
i = len(trim(name)) -3
if ( i > 1 ) then
if ( name(i:) == '_tau' ) then
print *, ''
print '(A,A)', "analyzing time displaced correlation ", name
call Cov_tau(name)
endif
endif
enddo
deallocate( names )
end Program ana
This diff is collapsed.
......@@ -44,7 +44,7 @@
Use MyMats
Use Matrix
Use Lattices_v3
Use Predefined_Lattices
Use Predefined_Lattices
use iso_fortran_env, only: output_unit, error_unit
Implicit none
......@@ -76,6 +76,8 @@
NAMELIST /VAR_errors/ n_skip, N_rebin, N_Cov, N_Back, N_auto
N_skip = 1
N_rebin = 1
N_Back = 1
N_auto = 0
OPEN(UNIT=5,FILE='parameters',STATUS='old',ACTION='read',IOSTAT=ierr)
......@@ -114,7 +116,7 @@
Close(10)
Write(6,*) "# of bins: ", Nbins
nbins = Nbins - n_skip
Write(6,*) "Effective # of bins: ", Nbins
Write(6,*) "Effective # of bins: ", Nbins/N_rebin
N_auto=min(N_auto,Nbins/3)
if(Nbins <= 1) then
write (error_unit,*) "Effective # of bins smaller than 2. Analysis impossible!"
......
......@@ -61,6 +61,8 @@
NAMELIST /VAR_errors/ n_skip, N_rebin, N_Cov, N_Back, N_auto
N_skip = 1
N_rebin = 1
N_auto=0
OPEN(UNIT=5,FILE='parameters',STATUS='old',ACTION='read',IOSTAT=ierr)
IF (ierr /= 0) THEN
......@@ -109,7 +111,7 @@
endif
OPEN (UNIT=21, FILE='Var_scalJ', STATUS='unknown')
WRITE(21,*) 'Effective number of bins, and bins: ', Nbins_eff, Nbins
WRITE(21,*) 'Effective number of bins, and bins: ', Nbins_eff/N_rebin, Nbins
ALLOCATE (EN(Nbins_eff), SIGN(Nbins_eff))
DO IOBS = 1,NOBS
WRITE(21,*)
......
......@@ -63,13 +63,14 @@
Real (Kind=Kind(0.d0)), allocatable :: Xk_p(:,:)
Character (len=64) :: File_out
NAMELIST /VAR_errors/ n_skip, N_rebin, N_Cov, N_Back, N_auto
NAMELIST /VAR_errors/ N_skip, N_rebin, N_Cov, N_Back, N_auto
N_skip = 1
N_rebin = 1
N_Back = 1
N_auto = 0
N_Cov = 0
OPEN(UNIT=5,FILE='parameters',STATUS='old',ACTION='read',IOSTAT=ierr)
IF (ierr /= 0) THEN
WRITE(error_unit,*) 'unable to open <parameters>',ierr
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment