From 2afa9a35eb29715fd1af52657d379e4af773f77a Mon Sep 17 00:00:00 2001 From: QuanSheng Wu Date: Tue, 7 May 2024 15:49:07 +0800 Subject: [PATCH] fixed a bug --- examples/Cu/wt.in-OHE-evolve | 2 +- src/Makefile.intel-mpi-with-ARPACK | 2 +- src/module.f90 | 18 +++++++++--------- src/readinput.f90 | 3 +++ src/sigma_OHE.f90 | 8 ++++---- 5 files changed, 18 insertions(+), 15 deletions(-) diff --git a/examples/Cu/wt.in-OHE-evolve b/examples/Cu/wt.in-OHE-evolve index 550fbba1..e2261be1 100644 --- a/examples/Cu/wt.in-OHE-evolve +++ b/examples/Cu/wt.in-OHE-evolve @@ -17,7 +17,7 @@ NumOccupied = 6 ! set it anyway even don't use it. Nk1 = 51 ! Kmesh(1) for KCUBE_BULK BTauNum= 100 ! Number of B*tau we calculate BTauMax = 40.0 ! The maximum B*tau, starting from Btau=0. -Nslice_BTau_Max = 20000 ! increase this number if negative magnetoresistance occurs, default =5000 +Nslice_BTau_Max = 5000 ! increase this number if negative magnetoresistance occurs, default =5000 / LATTICE diff --git a/src/Makefile.intel-mpi-with-ARPACK b/src/Makefile.intel-mpi-with-ARPACK index 1def7572..ecabdec9 100644 --- a/src/Makefile.intel-mpi-with-ARPACK +++ b/src/Makefile.intel-mpi-with-ARPACK @@ -15,7 +15,7 @@ OBJ = module.o sparse.o wt_aux.o math_lib.o symmetry.o \ # compiler F90 = mpif90 -fpp -DMPI -fpe3 -O3 -DARPACK -DINTELMKL #F90 = ifort -fpp -DINTELMKL -fpe3 -check all -traceback -g -CC = mpicc -fpp -O3 -DINTELMKL +CC = mpicc -cpp -O3 -DINTELMKL INCLUDE = -I${MKLROOT}/include WFLAG = -nogen-interface diff --git a/src/module.f90 b/src/module.f90 index 9d3a0358..28597d3e 100644 --- a/src/module.f90 +++ b/src/module.f90 @@ -450,11 +450,7 @@ module para logical :: valley_projection_calc ! for valley projection, only for sparse hamiltonina, you need to provide the valley operator logical :: w3d_nested_calc = .false. - - !> an integer to print the messages - !> iprint_level=3 : print all the debug messages - integer :: iprint_level = 1 - + namelist / Control / BulkBand_calc, BulkFS_calc, BulkFS_Plane_calc, & BulkFS_plane_stack_calc, BulkGap_plane_calc, & QPI_unfold_plane_calc, & @@ -486,7 +482,7 @@ module para LandauLevel_B_dos_calc,LanczosBand_calc,LanczosDos_calc, & LandauLevel_B_calc, LandauLevel_kplane_calc,landau_chern_calc, & FermiLevel_calc,ANE_calc, export_newhr,export_maghr,w3d_nested_calc, & - valley_projection_calc, iprint_level + valley_projection_calc integer :: Nslab ! Number of slabs for 2d Slab system integer :: Nslab1 ! Number of slabs for 1D wire system @@ -596,14 +592,18 @@ module para !> a real number to control when it's a cycle in subroutine RKF45_pack !> by default RKF45_PERIODIC_LEVEL= 1 - real(dp) :: RKF45_PERIODIC_LEVE + real(dp) :: RKF45_PERIODIC_LEVEL + + !> an integer to print the messages + !> iprint_level=3 : print all the debug messages + integer :: iprint_level = 1 !> namelist parameters namelist /PARAMETERS/ Eta_Arc,EF_broadening, OmegaNum, OmegaNum_unfold, OmegaMin, OmegaMax, & E_arc, Nk1, Nk2, Nk3, NP, Gap_threshold, Tmin, Tmax, NumT, & NBTau, BTauNum, BTauMax, Rcut, Magp, Magq, Magp_min, Magp_max, Nslice_BTau_Max, & - wcc_neighbour_tol, wcc_calc_tol, Beta,NumLCZVecs, & - Relaxation_Time_Tau, symprec, arpack_solver, RKF45_PERIODIC_LEVE, & + wcc_neighbour_tol, wcc_calc_tol, Beta,NumLCZVecs, iprint_level, & + Relaxation_Time_Tau, symprec, arpack_solver, RKF45_PERIODIC_LEVEL, & NumRandomConfs, NumSelectedEigenVals, projection_weight_mode, topsurface_atom_index real(Dp) :: E_fermi ! Fermi energy, search E-fermi in OUTCAR for VASP, set to zero for Wien2k diff --git a/src/readinput.f90 b/src/readinput.f90 index f8d25055..d40162c3 100644 --- a/src/readinput.f90 +++ b/src/readinput.f90 @@ -570,6 +570,7 @@ subroutine readinput topsurface_atom_index= 0 arpack_solver= 'zndrv1' RKF45_PERIODIC_LEVEL= 1 + iprint_level = 1 !> by default, we only project on atoms for a given wave function @@ -624,6 +625,8 @@ subroutine readinput write(stdout, '(1x, a, f16.5)')'Relaxation_Time_Tau (ps)', Relaxation_Time_Tau write(stdout, '(1x, a, f16.5)')'Rcut', Rcut write(stdout, '(1x, a, i16 )')'Magp', Magp + write(stdout, '(1x, a, i16 )')'iprint_level', iprint_level + write(stdout, '(1x, a, f16.2)')'RKF45_PERIODIC_LEVEL', RKF45_PERIODIC_LEVEL write(stdout, '(1x, a, i16 )')'Magp_min', Magp_min write(stdout, '(1x, a, i16 )')'Magp_max', Magp_max write(stdout, '(1x, a, f16.5)')'wcc_calc_tol', wcc_calc_tol diff --git a/src/sigma_OHE.f90 b/src/sigma_OHE.f90 index 0fd5c2b5..ebf020fd 100644 --- a/src/sigma_OHE.f90 +++ b/src/sigma_OHE.f90 @@ -1526,7 +1526,7 @@ end subroutine velocity_calc_iband !> a pack for rkf45 subroutine RKF45_pack(magnetic_field, iband, NSlice_Btau, k_start, Btau_start, Btau_final, kout, icycle, fail) use wmpi - use para, only : dp, eps9, eps8, cpuid, stdout, eps6, Bdirection, eps3 + use para, only : dp, eps9, eps8, cpuid, stdout, eps6, Bdirection, eps3, RKF45_PERIODIC_LEVEL implicit none !> inout variable @@ -1810,7 +1810,7 @@ subroutine evolve_k_ohe do ib=1, Nband_Fermi_Level if (cpuid.eq.0) then write(bandname, '(i10)')bands_fermi_level(ib) - write(evolvefilename, '(3a)')'evolve_band_', trim(adjustl(bandname)), '.dat' + write(evolvefilename, '(3a)')'evolve_band_', trim(adjustl(bandname)), '.txt' open(unit=myfileindex(ib), file=evolvefilename) write(myfileindex(ib), '(a10, i5, a16, f16.8)')'# evolve k ', bands_fermi_level(ib), 'energy level', Ek(ib) write(myfileindex(ib), '("#", a13, 24a16)')'BTau (T.ps)', & @@ -1889,7 +1889,7 @@ subroutine evolve_k_ohe write(myfileindex(ib), '(1X, a, 3f12.5, a, i8)')"# Starting k point (fractional coordinates) :", k_start, ' ith-band ', bands_fermi_level(ib) write(myfileindex(ib), '(a, f16.8, a)')'# At energy level', Ek(ib)/eV2Hartree, ' eV' write(myfileindex(ib), '("#", a13, 20a16)')'BTau (T.ps)', & - 'Omega*Tau', 'kx', 'ky', 'kz', 'vx', 'vy', 'vz', 'k1', 'k2','k3' + 'kx', 'ky', 'kz', 'vx', 'vy', 'vz', 'k1', 'k2','k3', "vx'", "vy'", "vz'" do it=1, NSlice_Btau_all(ik) k= kout_all(:, it, ik) !call velocity_calc_iband(ib+bands_fermi_level(1)-1, k, v_t) @@ -1898,7 +1898,7 @@ subroutine evolve_k_ohe call direct_cart_rec(kout_all(:, it, ik), k) call project_k3_to_kplane_defined_by_direction(v_t, Bdirection, v_t2) write(myfileindex(ib), '(100f16.8)')Btau*Magneticfluxdensity_atomic/Relaxation_Time_Tau, & - k, v_t, kout_all(:, it, ik), E_iband, v_t2 + k, v_t, kout_all(:, it, ik), v_t2 enddo endif endif