TABLE OF CONTENTS
- ABINIT/m_chebfiwf
- m_chebfiwf/chebfiwf2
- m_chebfiwf/getBm1X
- m_chebfiwf/getghc_gsc1
- m_chebfiwf/precond1
ABINIT/m_chebfiwf [ Functions ]
NAME
m_chebfiwf
FUNCTION
This module contains a routine updating the whole wave functions at a given k-point, using the Chebyshev filtering method (2021 implementation using xG abstraction layer) for a given spin-polarization, from a fixed hamiltonian but might also simply compute eigenvectors and eigenvalues at this k point. it will also update the matrix elements of the hamiltonian.
COPYRIGHT
Copyright (C) 2018-2022 ABINIT group (BS) This file is distributed under the terms of the gnu general public license, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . for the initials of contributors, see ~abinit/doc/developers/contributors.txt .
SOURCE
21 #if defined HAVE_CONFIG_H 22 #include "config.h" 23 #endif 24 25 #include "abi_common.h" 26 27 ! nvtx related macro definition 28 #include "nvtx_macros.h" 29 30 module m_chebfiwf 31 32 use defs_abitypes 33 use defs_basis 34 use m_abicore 35 use m_errors 36 use m_fstrings 37 use m_time 38 39 use m_chebfi 40 use m_chebfi2 41 use m_invovl 42 43 use m_cgtools, only : dotprod_g 44 use m_dtset, only : dataset_type 45 46 use m_hamiltonian, only : gs_hamiltonian_type 47 use m_pawcprj, only : pawcprj_type 48 use m_nonlop, only : nonlop 49 use m_prep_kgb, only : prep_getghc, prep_nonlop 50 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free 51 use m_getghc, only : multithreaded_getghc 52 53 use m_xg 54 use m_xgTransposer 55 56 #if defined(HAVE_GPU_CUDA) && defined(HAVE_GPU_NVTX_V3) 57 use m_nvtx_data 58 #endif 59 60 use, intrinsic :: iso_c_binding, only: c_associated,c_loc,c_ptr,c_f_pointer 61 62 use m_xmpi 63 use m_xomp 64 #ifdef HAVE_OPENMP 65 use omp_lib 66 #endif 67 68 implicit none 69 70 private 71 72 integer, parameter :: l_tim_getghc=7 73 real(dp), parameter :: inv_sqrt2 = 1/sqrt2 74 75 ! For use in getghc_gsc1 76 integer, save :: l_cpopt 77 integer, save :: l_icplx 78 integer, save :: l_istwf 79 integer, save :: l_npw 80 integer, save :: l_nband_filter 81 integer, save :: l_nspinor 82 logical, save :: l_paw 83 integer, save :: l_prtvol 84 integer, save :: l_sij_opt 85 integer, save :: l_paral_kgb 86 integer, save :: l_useria 87 integer, save :: l_block_sliced 88 real(dp), allocatable,save :: l_pcon(:) 89 type(mpi_type),pointer,save :: l_mpi_enreg 90 type(gs_hamiltonian_type),pointer,save :: l_gs_hamk 91 92 integer, parameter :: DEBUG_ROWS = 5 93 integer, parameter :: DEBUG_COLUMNS = 5 94 95 public :: chebfiwf2 96 97 CONTAINS !========================================================================================
m_chebfiwf/chebfiwf2 [ Functions ]
[ Top ] [ m_chebfiwf ] [ Functions ]
NAME
chebfiwf2
FUNCTION
This routine updates the whole wave functions set at a given k-point, using the Chebfi method (2021 version using xG abstraction layer)
INPUTS
dtset= input variables for this dataset kinpw(npw)= kinetic energy for each plane wave (Hartree) mpi_enreg= MPI-parallelisation information nband= number of bands at this k point npw= number of plane waves at this k point nspinor= number of spinorial components of the wavefunctions prtvol= control print volume and debugging
OUTPUT
eig(nband)= eigenvalues (hartree) for all bands enl_out(nband)= contribution of each band to the nl part of energy resid(nband)= residuals for each band
SIDE EFFECTS
cg(2,npw*nspinor*nband)= planewave coefficients of wavefunctions gs_hamk <type(gs_hamiltonian_type)>=all data for the hamiltonian at k
SOURCE
128 subroutine chebfiwf2(cg,dtset,eig,enl_out,gs_hamk,kinpw,mpi_enreg,& 129 & nband,npw,nspinor,prtvol,resid) 130 131 implicit none 132 133 !Arguments ------------------------------------ 134 integer,intent(in) :: nband,npw,prtvol,nspinor 135 type(mpi_type),target,intent(inout) :: mpi_enreg 136 real(dp),target,intent(inout) :: cg(2,npw*nspinor*nband) 137 real(dp),intent(in) :: kinpw(npw) 138 real(dp),target,intent(out) :: resid(nband) 139 real(dp),intent(out) :: enl_out(nband) 140 real(dp),target,intent(out) :: eig(nband) 141 type(dataset_type),intent(in) :: dtset 142 type(gs_hamiltonian_type),target,intent(inout) :: gs_hamk 143 144 !Local variables------------------------------- 145 !scalars 146 integer, parameter :: tim_chebfiwf2 = 1750 147 integer :: ipw,space,blockdim,nline,total_spacedim,ierr,nthreads 148 real(dp) :: cputime,walltime,localmem 149 type(c_ptr) :: cptr 150 type(chebfi_t) :: chebfi 151 type(xgBlock_t) :: xgx0,xgeigen,xgresidu 152 !arrays 153 real(dp) :: tsec(2),chebfiMem(2) 154 real(dp),pointer :: eig_ptr(:,:) => NULL() 155 real(dp),pointer :: resid_ptr(:,:) => NULL() 156 real(dp), allocatable :: l_gvnlxc(:,:) 157 158 !Stupid things for NC 159 integer,parameter :: choice=1, paw_opt=0, signs=1 160 real(dp) :: gsc_dummy(0,0) 161 type(pawcprj_type) :: cprj_dum(gs_hamk%natom,0) 162 163 ! ********************************************************************* 164 165 !################ INITIALIZATION ##################################### 166 !###################################################################### 167 168 call timab(tim_chebfiwf2,1,tsec) 169 cputime = abi_cpu_time() 170 walltime = abi_wtime() 171 172 !Set module variables 173 l_paw = (gs_hamk%usepaw==1) 174 l_cpopt=-1;l_sij_opt=0;if (l_paw) l_sij_opt=1 175 l_istwf=gs_hamk%istwf_k 176 l_npw = npw 177 l_nspinor = nspinor 178 l_prtvol = prtvol 179 l_mpi_enreg => mpi_enreg 180 l_gs_hamk => gs_hamk 181 l_nband_filter = nband 182 l_paral_kgb = dtset%paral_kgb 183 l_block_sliced = dtset%diago_apply_block_sliced 184 185 !Variables 186 nline=dtset%nline 187 blockdim=l_mpi_enreg%nproc_band*l_mpi_enreg%bandpp 188 !for debug 189 l_useria=dtset%useria 190 191 !Depends on istwfk 192 if ( l_istwf == 2 ) then ! Real only 193 ! SPACE_CR mean that we have complex numbers but no re*im terms only re*re 194 ! and im*im so that a vector of complex is consider as a long vector of real 195 ! therefore the number of data is (2*npw*nspinor)*nband 196 ! This space is completely equivalent to SPACE_R but will correctly set and 197 ! get the array data into the xgBlock 198 space = SPACE_CR 199 l_icplx = 2 200 else ! complex 201 space = SPACE_C 202 l_icplx = 1 203 end if 204 205 !Memory info 206 if ( prtvol >= 3 ) then 207 if (l_mpi_enreg%paral_kgb == 1) then 208 total_spacedim = l_icplx*l_npw*l_nspinor 209 call xmpi_sum(total_spacedim,l_mpi_enreg%comm_bandspinorfft,ierr) 210 else 211 total_spacedim = 0 212 end if 213 chebfiMem = chebfi_memInfo(nband,l_icplx*l_npw*l_nspinor,space,l_mpi_enreg%paral_kgb, & 214 & total_spacedim,l_mpi_enreg%bandpp) !blockdim 215 localMem = (l_npw+2*l_npw*l_nspinor+2*nband)*kind(1.d0) !blockdim 216 write(std_out,'(1x,A,F10.6,1x,A)') "Each MPI process calling chebfi should need around ", & 217 (localMem+sum(chebfiMem))/1e9,"GB of peak memory as follows :" 218 write(std_out,'(4x,A,F10.6,1x,A)') "Permanent memory in chebfiwf : ",(localMem)/1e9,"GB" 219 write(std_out,'(4x,A,F10.6,1x,A)') "Permanent memory in m_chebfi : ",(chebfiMem(1))/1e9,"GB" 220 write(std_out,'(4x,A,F10.6,1x,A)') "Temporary memory in m_chebfi : ",(chebfiMem(2))/1e9,"GB" 221 end if 222 223 !For preconditionning 224 ABI_MALLOC(l_pcon,(1:l_icplx*npw)) 225 226 !$omp parallel do schedule(static), shared(l_pcon,kinpw) 227 do ipw=1-1,l_icplx*npw-1 228 if(kinpw(ipw/l_icplx+1)>huge(0.0_dp)*1.d-11) then 229 l_pcon(ipw+1)=0.d0 230 else 231 l_pcon(ipw+1) = (27+kinpw(ipw/l_icplx+1)*(18+kinpw(ipw/l_icplx+1)*(12+8*kinpw(ipw/l_icplx+1)))) & 232 & / (27+kinpw(ipw/l_icplx+1)*(18+kinpw(ipw/l_icplx+1)*(12+8*kinpw(ipw/l_icplx+1))) + 16*kinpw(ipw/l_icplx+1)**4) 233 end if 234 end do 235 236 call xgBlock_map(xgx0,cg,space,l_icplx*l_npw*l_nspinor,nband,l_mpi_enreg%comm_bandspinorfft) 237 238 if ( l_istwf == 2 ) then ! Real only 239 ! Scale cg 240 call xgBlock_scale(xgx0,sqrt2,1) !ALL MPI processes do this 241 242 ! This is possible since the memory in cg and xgx0 is the same 243 ! Don't know yet how to deal with this with xgBlock 244 !MPI HANDLES THIS AUTOMATICALLY (only proc 0 is me_g0) 245 if(l_mpi_enreg%me_g0 == 1) cg(:, 1:npw*nspinor*nband:npw) = cg(:, 1:npw*nspinor*nband:npw) * inv_sqrt2 246 end if 247 248 !Trick with C is to change rank of arrays (:) to (:,:) 249 cptr = c_loc(eig) 250 call c_f_pointer(cptr,eig_ptr,(/ nband,1 /)) 251 call xgBlock_map(xgeigen,eig_ptr,SPACE_R,nband,1,l_mpi_enreg%comm_bandspinorfft) 252 !Trick the with C to change rank of arrays (:) to (:,:) 253 cptr = c_loc(resid) 254 call c_f_pointer(cptr,resid_ptr,(/ nband,1 /)) 255 call xgBlock_map(xgresidu,resid_ptr,SPACE_R,nband,1,l_mpi_enreg%comm_bandspinorfft) 256 257 ! ABI_MALLOC(l_gvnlxc,(2,l_npw*l_nspinor*l_nband_filter)) 258 call timab(tim_chebfiwf2,2,tsec) 259 260 cputime = abi_cpu_time() - cputime 261 walltime = abi_wtime() - walltime 262 263 nthreads = xomp_get_num_threads(open_parallel = .true.) 264 265 if ( cputime/walltime/dble(nthreads) < 0.75 .and. (int(cputime/walltime)+1) /= nthreads) then 266 if ( prtvol >= 3 ) then 267 write(std_out,'(a)',advance='no') sjoin(" Chebfi took", sec2str(cputime), "of cpu time") 268 write(std_out,*) sjoin("for a wall time of", sec2str(walltime)) 269 write(std_out,'(a,f6.2)') " -> Ratio of ", cputime/walltime 270 end if 271 ABI_COMMENT(sjoin("You should set the number of threads to something close to",itoa(int(cputime/walltime)+1))) 272 end if 273 274 call chebfi_init(chebfi,nband,l_icplx*l_npw*l_nspinor,dtset%tolwfr,dtset%ecut, & 275 & dtset%paral_kgb,l_mpi_enreg%nproc_band,l_mpi_enreg%bandpp, & 276 & l_mpi_enreg%nproc_fft,nline, space,1,l_gs_hamk%istwf_k, & 277 & l_mpi_enreg%comm_bandspinorfft,l_mpi_enreg%me_g0,l_paw) 278 279 !################ RUUUUUUUN ##################################### 280 !###################################################################### 281 282 call chebfi_run(chebfi,xgx0,getghc_gsc1,getBm1X,precond1,xgeigen,xgresidu,l_mpi_enreg) 283 284 !Free preconditionning since not needed anymore 285 ABI_FREE(l_pcon) 286 287 !Compute enlout (nonlocal energy for each band if necessary) This is the best 288 ! quick and dirty trick to compute this part in NC. gvnlc cannot be part of 289 ! chebfi algorithm 290 if ( .not. l_paw ) then 291 !Check l_gvnlc size 292 !if ( size(l_gvnlxc) < 2*nband*l_npw*l_nspinor ) then 293 !if ( size(l_gvnlxc) /= 0 ) then 294 ! ABI_FREE(l_gvnlxc) 295 ABI_MALLOC(l_gvnlxc,(0,0)) 296 !end if 297 298 ABI_NVTX_START_RANGE(NVTX_CHEBFI2_NONLOP) 299 !Call nonlop 300 call nonlop(choice,l_cpopt,cprj_dum,enl_out,l_gs_hamk,0,eig,mpi_enreg,nband,1,paw_opt,& 301 & signs,gsc_dummy,l_tim_getghc,cg,l_gvnlxc) 302 ABI_NVTX_END_RANGE() 303 ABI_FREE(l_gvnlxc) 304 end if 305 306 !Free chebfi 307 call chebfi_free(chebfi) 308 309 !################ SORRY IT'S ALREADY FINISHED : ) ################# 310 !###################################################################### 311 312 call timab(tim_chebfiwf2,2,tsec) 313 314 cputime = abi_cpu_time() - cputime 315 walltime = abi_wtime() - walltime 316 317 nthreads = xomp_get_num_threads(open_parallel = .true.) 318 319 if ( cputime/walltime/dble(nthreads) < 0.75 .and. (int(cputime/walltime)+1) /= nthreads) then 320 if ( prtvol >= 3 ) then 321 write(std_out,'(a)',advance='no') sjoin(" Chebfi took", sec2str(cputime), "of cpu time") 322 write(std_out,*) sjoin("for a wall time of", sec2str(walltime)) 323 write(std_out,'(a,f6.2)') " -> Ratio of ", cputime/walltime 324 end if 325 ABI_COMMENT(sjoin("You should set the number of threads to something close to",itoa(int(cputime/walltime)+1))) 326 end if 327 328 DBG_EXIT("COLL") 329 330 end subroutine chebfiwf2
m_chebfiwf/getBm1X [ Functions ]
[ Top ] [ m_chebfiwf ] [ Functions ]
NAME
getBm1X
FUNCTION
This routine computes S^-1|C> for a given wave function C. It acts as a driver for apply_invovl.
SIDE EFFECTS
X <type(xgBlock_t)>= memory block containing |C> Bm1X <type(xgBlock_t)>= memory block containing S^-1|C> transposer <type(xgTransposer_t)>= data used for array transpositions
SOURCE
463 subroutine getBm1X(X,Bm1X,transposer) 464 465 implicit none 466 467 !Arguments ------------------------------------ 468 type(xgBlock_t), intent(inout) :: X 469 type(xgBlock_t), intent(inout) :: Bm1X 470 type(xgTransposer_t), optional, intent(inout) :: transposer 471 472 !Local variables------------------------------- 473 !scalars 474 integer :: blockdim 475 integer :: spacedim 476 integer :: cpuRow 477 !arrays 478 real(dp), pointer :: ghc_filter(:,:) 479 real(dp), pointer :: gsm1hc_filter(:,:) 480 type(pawcprj_type), allocatable :: cwaveprj_next(:,:) !dummy 481 482 ! ********************************************************************* 483 484 call xgBlock_getSize(X,spacedim,blockdim) 485 486 spacedim = spacedim/l_icplx 487 488 call xgBlock_reverseMap(X,ghc_filter,l_icplx,spacedim*blockdim) 489 490 call xgBlock_reverseMap(Bm1X,gsm1hc_filter,l_icplx,spacedim*blockdim) 491 492 !scale back cg 493 if(l_istwf == 2) then 494 call xgBlock_scale(X,inv_sqrt2,1) 495 if (l_paral_kgb == 0) then 496 if(l_mpi_enreg%me_g0 == 1) ghc_filter(:, 1:spacedim*blockdim:l_npw) = ghc_filter(:, 1:spacedim*blockdim:l_npw) * sqrt2 497 else 498 cpuRow = xgTransposer_getRank(transposer, 2) 499 if (cpuRow == 0) then 500 ghc_filter(:, 1:spacedim*blockdim:spacedim) = ghc_filter(:, 1:spacedim*blockdim:spacedim) * sqrt2 501 end if 502 end if 503 if(l_paw) then 504 call xgBlock_scale(Bm1X,inv_sqrt2,1) 505 if (l_paral_kgb == 0) then 506 if(l_mpi_enreg%me_g0 == 1) & 507 & gsm1hc_filter(:, 1:spacedim*blockdim:l_npw) = gsm1hc_filter(:, 1:spacedim*blockdim:l_npw) * sqrt2 508 else 509 if (cpuRow == 0) then 510 gsm1hc_filter(:, 1:spacedim*blockdim:spacedim) = gsm1hc_filter(:, 1:spacedim*blockdim:spacedim) * sqrt2 511 end if 512 end if 513 end if 514 end if 515 516 if(l_paw) then 517 ABI_MALLOC(cwaveprj_next, (l_gs_hamk%natom,l_nspinor*blockdim)) 518 call pawcprj_alloc(cwaveprj_next,0,l_gs_hamk%dimcprj) 519 call apply_invovl(l_gs_hamk, ghc_filter(:,:), gsm1hc_filter(:,:), cwaveprj_next(:,:), & 520 & spacedim, blockdim, l_mpi_enreg, l_nspinor, l_block_sliced) 521 else 522 gsm1hc_filter(:,:) = ghc_filter(:,:) 523 end if 524 525 !Scale cg, ghc, gsc 526 if ( l_istwf == 2 ) then 527 call xgBlock_scale(X,sqrt2,1) 528 if (l_paral_kgb == 0) then 529 if(l_mpi_enreg%me_g0 == 1) then 530 ghc_filter(:, 1:spacedim*blockdim:l_npw) = ghc_filter(:, 1:spacedim*blockdim:l_npw) * inv_sqrt2 531 endif 532 else 533 if (cpuRow == 0) then 534 ghc_filter(:, 1:spacedim*blockdim:spacedim) = ghc_filter(:, 1:spacedim*blockdim:spacedim) * inv_sqrt2 535 end if 536 end if 537 if(l_paw) then 538 call xgBlock_scale(Bm1X,sqrt2,1) 539 if (l_paral_kgb == 0) then 540 if(l_mpi_enreg%me_g0 == 1) & 541 & gsm1hc_filter(:, 1:spacedim*blockdim:l_npw) = gsm1hc_filter(:, 1:spacedim*blockdim:l_npw) * inv_sqrt2 542 else 543 if (cpuRow == 0) then 544 gsm1hc_filter(:, 1:spacedim*blockdim:spacedim) = gsm1hc_filter(:, 1:spacedim*blockdim:spacedim) * inv_sqrt2 545 end if 546 end if 547 end if 548 end if 549 550 if (l_paw) then 551 if (l_useria /= 121212) then 552 call pawcprj_free(cwaveprj_next) 553 ABI_FREE(cwaveprj_next) 554 end if 555 end if 556 557 end subroutine getBm1X
m_chebfiwf/getghc_gsc1 [ Functions ]
[ Top ] [ m_chebfiwf ] [ Functions ]
NAME
getghc_gsc1
FUNCTION
This routine computes H|C> and possibly S|C> for a given wave function C. It acts as a driver for getghc, taken into account parallelism, multithreading, etc.
SIDE EFFECTS
X <type(xgBlock_t)>= memory block containing |C> AX <type(xgBlock_t)>= memory block containing H|C> BX <type(xgBlock_t)>= memory block containing S|C> transposer <type(xgTransposer_t)>= data used for array transpositions
SOURCE
351 subroutine getghc_gsc1(X,AX,BX,transposer) 352 353 implicit none 354 355 !Arguments ------------------------------------ 356 type(xgBlock_t), intent(inout) :: X 357 type(xgBlock_t), intent(inout) :: AX 358 type(xgBlock_t), intent(inout) :: BX 359 type(xgTransposer_t), optional, intent(inout) :: transposer 360 integer :: blockdim 361 integer :: spacedim 362 type(pawcprj_type) :: cprj_dum(l_gs_hamk%natom,0) 363 364 !Local variables------------------------------- 365 !scalars 366 integer :: cpuRow 367 real(dp) :: eval 368 !arrays 369 real(dp), pointer :: cg(:,:) 370 real(dp), pointer :: ghc(:,:) 371 real(dp), pointer :: gsc(:,:) 372 real(dp), allocatable :: l_gvnlxc(:,:) 373 374 ! ********************************************************************* 375 376 ABI_NVTX_START_RANGE(NVTX_GETGHC) 377 378 call xgBlock_getSize(X,spacedim,blockdim) 379 380 spacedim = spacedim/l_icplx 381 382 call xgBlock_reverseMap(X,cg,l_icplx,spacedim*blockdim) 383 call xgBlock_reverseMap(AX,ghc,l_icplx,spacedim*blockdim) 384 call xgBlock_reverseMap(BX,gsc,l_icplx,spacedim*blockdim) 385 386 !Scale back cg 387 if(l_istwf == 2) then 388 call xgBlock_scale(X,inv_sqrt2,1) 389 390 if (l_paral_kgb == 0) then 391 if(l_mpi_enreg%me_g0 == 1) cg(:, 1:spacedim*blockdim:l_npw) = cg(:, 1:spacedim*blockdim:l_npw) * sqrt2 392 else 393 cpuRow = xgTransposer_getRank(transposer, 2) 394 if (cpuRow == 0) then 395 cg(:, 1:spacedim*blockdim:spacedim) = cg(:, 1:spacedim*blockdim:spacedim) * sqrt2 396 end if 397 end if 398 end if 399 400 !if ( size(l_gvnlxc) < 2*blockdim*spacedim ) then 401 ! ABI_FREE(l_gvnlxc) 402 ! ABI_MALLOC(l_gvnlxc,(2,blockdim*spacedim)) 403 !end if 404 ABI_MALLOC(l_gvnlxc,(0,0)) 405 406 call multithreaded_getghc(l_cpopt,cg,cprj_dum,ghc,gsc,& 407 l_gs_hamk,l_gvnlxc,eval,l_mpi_enreg,blockdim,l_prtvol,l_sij_opt,l_tim_getghc,0) 408 409 ABI_FREE(l_gvnlxc) 410 411 !Scale cg, ghc, gsc 412 if ( l_istwf == 2 ) then 413 call xgBlock_scale(X,sqrt2,1) 414 call xgBlock_scale(AX,sqrt2,1) 415 416 if (l_paral_kgb == 0) then 417 if(l_mpi_enreg%me_g0 == 1) then 418 cg(:, 1:spacedim*blockdim:l_npw) = cg(:, 1:spacedim*blockdim:l_npw) * inv_sqrt2 419 ghc(:, 1:spacedim*blockdim:l_npw) = ghc(:, 1:spacedim*blockdim:l_npw) * inv_sqrt2 420 endif 421 else 422 if (cpuRow == 0) then 423 cg(:, 1:spacedim*blockdim:spacedim) = cg(:, 1:spacedim*blockdim:spacedim) * inv_sqrt2 424 ghc(:, 1:spacedim*blockdim:spacedim) = ghc(:, 1:spacedim*blockdim:spacedim) * inv_sqrt2 425 end if 426 end if 427 if(l_paw) then 428 call xgBlock_scale(BX,sqrt2,1) 429 if (l_paral_kgb == 0) then 430 if(l_mpi_enreg%me_g0 == 1) gsc(:, 1:spacedim*blockdim:l_npw) = gsc(:, 1:spacedim*blockdim:l_npw) * inv_sqrt2 431 else 432 if (cpuRow == 0) then 433 gsc(:, 1:spacedim*blockdim:spacedim) = gsc(:, 1:spacedim*blockdim:spacedim) * inv_sqrt2 434 end if 435 end if 436 end if 437 end if 438 439 if ( .not. l_paw ) call xgBlock_copy(X,BX) 440 441 ABI_NVTX_END_RANGE() 442 443 end subroutine getghc_gsc1
m_chebfiwf/precond1 [ Functions ]
[ Top ] [ m_chebfiwf ] [ Functions ]
NAME
precond1
FUNCTION
This routine applies a preconditionning to a block of memory
SIDE EFFECTS
W <type(xgBlock_t)>= memory block
SOURCE
574 subroutine precond1(W) 575 576 implicit none 577 578 !Arguments ------------------------------------ 579 type(xgBlock_t), intent(inout) :: W 580 581 !Local variables------------------------------- 582 !scalars 583 integer :: ispinor 584 585 ! ********************************************************************* 586 587 !Precondition resid_vec 588 do ispinor = 1,l_nspinor 589 call xgBlock_colwiseMul(W,l_pcon,l_icplx*l_npw*(ispinor-1)) 590 end do 591 592 end subroutine precond1