1! 2! Copyright (C) 2002-2010 Quantum ESPRESSO groups 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8 9subroutine gtable( ipol, ctable) 10 11 ! this subroutine prepares the correspondence array to 12 ! compute the operator exp(iG_ipol.r) 13 14 ! ctable : output correspondence table 15 ! in (ig,1) correspondence for g+1 16 ! in (ig,2) correspondence for (-g)+1 17 ! we use the rule: if non point ngw+1 18 ! if found positive = normal 19 ! negative = conjugate 20 ! ipol : input polarization direction 21 ! a orthorombic primitive cell is supposed 22 use kinds, only: dp 23 use gvecw, only: ngw 24 use gvect, only: mill 25 use mp, only: mp_sum 26 use io_global, only: ionode, stdout 27 use mp_global, only: intra_bgrp_comm 28 29 implicit none 30 integer :: ipol, ctable(ngw,2) 31 !local variables 32 integer :: i,j,k, ig, jg 33 logical :: found 34 real(dp) :: test 35 36 test=0.d0 37 do ig=1,ngw!loop on g vectors 38 ! first +g 39 i = mill(1,ig) 40 j = mill(2,ig) 41 k = mill(3,ig) 42 if(ipol.eq.1) i=i+1 43 if(ipol.eq.2) j=j+1 44 if(ipol.eq.3) k=k+1 45 46 found = .false. 47 48 do jg=1,ngw 49 if(mill(1,jg).eq.i .and. mill(2,jg).eq.j .and. mill(3,jg).eq.k) then 50 found=.true. 51 ctable(ig,1)=jg 52 endif 53 enddo 54 55 if(.not. found) then 56 do jg=1,ngw 57 if(-mill(1,jg).eq.i .and. -mill(2,jg).eq.j .and. -mill(3,jg).eq.k) then 58 found=.true. 59 ctable(ig,1)=-jg 60 endif 61 enddo 62 if(.not. found) then 63 ctable(ig,1)= ngw+1 64 test=test+1.d0 65 endif 66 endif 67 68 ! now -g 69 i = -mill(1,ig) 70 j = -mill(2,ig) 71 k = -mill(3,ig) 72 if(ipol.eq.1) i=i+1 73 if(ipol.eq.2) j=j+1 74 if(ipol.eq.3) k=k+1 75 76 found = .false. 77 78 do jg=1,ngw 79 if (-mill(1,jg).eq.i .and. -mill(2,jg).eq.j .and. -mill(3,jg).eq.k)then 80 found=.true. 81 ctable(ig,2)=-jg 82 endif 83 enddo 84 85 if(.not.found) then 86 do jg=1,ngw 87 if(mill(1,jg).eq.i .and. mill(2,jg).eq.j .and. mill(3,jg).eq.k)then 88 found=.true. 89 ctable(ig,2)=jg 90 endif 91 enddo 92 if(.not.found) then 93 ctable(ig,2)=ngw+1 94 test=test+1.d0 95 endif 96 endif 97 enddo 98 99 call mp_sum(test, intra_bgrp_comm) 100 if(ionode) write(stdout,*) '#not found, gtable: ', test 101 102 return 103end subroutine gtable 104 105subroutine gtablein( ipol, ctabin) 106 107 ! this subroutine prepare the inverse correspondence array to 108 ! compute the operator exp(iG_ipol.r) 109 110 ! ctabin(ngw,2) : output correspondence table 111 ! if negative to take complex conjugate, 1 g'+1, 2 g' -1 112 ! if not found = ngw+1 113 ! ipol : input polarization direction 114 ! a orthorombic primitive cell is supposed 115 116 use kinds, only: dp 117 use gvecw, only: ngw 118 use gvect, only: mill 119 use mp, only: mp_sum 120 use io_global, only: ionode, stdout 121 use mp_global, only: intra_bgrp_comm 122 123 implicit none 124 125 integer :: ipol, ctabin(ngw,2) 126 127 !local variables 128 integer :: i,j,k, ig, jg 129 logical :: found 130 real(dp) :: test 131 132 test=0.d0 133 134 do ig=1,ngw!loop on g vectors 135 i = mill(1,ig) 136 j = mill(2,ig) 137 k = mill(3,ig) 138 if(ipol.eq.1) i=i+1 139 if(ipol.eq.2) j=j+1 140 if(ipol.eq.3) k=k+1 141 found = .false. 142 143 do jg=1,ngw 144 if(i.eq.mill(1,jg).and. j.eq.mill(2,jg) .and. k.eq.mill(3,jg))then 145 found = .true. 146 ctabin(ig,1)=jg 147 else if(i.eq.-mill(1,jg).and. j.eq.-mill(2,jg) .and. k.eq.-mill(3,jg))then 148 found=.true. 149 ctabin(ig,1)=-jg 150 endif 151 enddo 152 if(.not.found) then 153 ctabin(ig,1)=ngw+1 154 test=test+1 155 endif 156 enddo 157 158 do ig=1,ngw!loop on g vectors 159 i = mill(1,ig) 160 j = mill(2,ig) 161 k = mill(3,ig) 162 if(ipol.eq.1) i=i-1 163 if(ipol.eq.2) j=j-1 164 if(ipol.eq.3) k=k-1 165 found = .false. 166 167 do jg=1,ngw 168 if(i.eq.mill(1,jg).and. j.eq.mill(2,jg) .and. k.eq.mill(3,jg))then 169 found = .true. 170 ctabin(ig,2)=jg 171 else if(i.eq.-mill(1,jg).and. j.eq.-mill(2,jg) .and. k.eq.-mill(3,jg))then 172 found=.true. 173 ctabin(ig,2)=-jg 174 endif 175 enddo 176 if(.not.found) then 177 ctabin(ig,2)=ngw+1 178 test=test+1 179 endif 180 enddo 181 182 call mp_sum(test, intra_bgrp_comm) 183 if(ionode) write(stdout,*) '#not found, gtabin: ', test 184 185 return 186 187end subroutine gtablein 188 189 190 191subroutine find_whose_is_g 192!this subroutine set the correspondence G-->Proc 193 194 USE gvecw, ONLY : ngw, ngw_g 195 USE gvect, ONLY : ig_l2g, mill_g, mill 196 USE mp, ONLY : mp_sum 197 USE io_global, ONLY : stdout 198 USE mp_global, ONLY : me_bgrp, nproc_bgrp, intra_bgrp_comm 199 USE efield_module, ONLY : whose_is_g 200 201 implicit none 202 203 INTEGER :: ig 204 205 whose_is_g(:)=0 206 207 208 do ig=1,ngw 209 if(ig_l2g(ig) > ngw_g) then 210 write(stdout,*) 'find_whose_is_g: too large' 211 stop 212 endif 213 whose_is_g(ig_l2g(ig))=me_bgrp+1 214 enddo 215 call mp_sum(whose_is_g,intra_bgrp_comm) 216 whose_is_g(:)=whose_is_g(:)-1 217 ! mill_g is used in gtable_missing and re-initialized here 218 ! workaround by PG to avoid a large array like mill_g allocated all the time 219 220 allocate ( mill_g(3,ngw_g) ) 221 do ig=1,ngw 222 mill_g(:,ig_l2g(ig)) = mill(:,ig) 223 end do 224 call mp_sum(mill_g,intra_bgrp_comm) 225 226 227return 228end subroutine find_whose_is_g 229 230 231subroutine gtable_missing 232 233 USE efield_module, ONLY : ctable_missing_1, ctable_missing_2, & 234 whose_is_g,n_g_missing_p, & 235 ctable_missing_rev_1,ctable_missing_rev_2 236 USE gvecw, ONLY : ngw, ngw_g 237 USE gvect, ONLY : ig_l2g, mill_g, mill, gstart 238 USE mp, ONLY : mp_sum, mp_max, mp_alltoall 239 USE io_global, ONLY : stdout 240 USE mp_global, ONLY : me_bgrp, nproc_bgrp, intra_bgrp_comm 241 USE parallel_include 242 243 implicit none 244 245 INTEGER :: ipol, i,j,k,ig,igg, nfound_max, ip 246 LOGICAL :: found 247 INTEGER :: nfound_proc(nproc_bgrp,2) 248 INTEGER, ALLOCATABLE :: igg_found(:,:,:), ig_send(:,:,:), igg_found_snd(:,:,:) 249 INTEGER, ALLOCATABLE :: igg_found_rcv(:,:,:) 250 INTEGER :: ierr,sndint,rcvint 251 252 allocate( igg_found(ngw_g,2,nproc_bgrp), ig_send(ngw_g,2,nproc_bgrp) ) 253 do ipol=1,2 254 255 nfound_max=0 256 nfound_proc(:,:)=0 257 ig_send(:,:,:)=0 258 259 do ig=1,ngw!loop on g vectors 260 ! first +g 261 i = mill(1,ig) 262 j = mill(2,ig) 263 k = mill(3,ig) 264 if(ipol.eq.1) i=i+1 265 if(ipol.eq.2) j=j+1 266 if(ipol.eq.3) k=k+1 267 do igg=1,ngw_g 268 if( i==mill_g(1,igg) .and. j==mill_g(2,igg) .and. k==mill_g(3,igg)) then 269 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 270 nfound_max=nfound_max+1 271 nfound_proc(whose_is_g(igg)+1,1)=nfound_proc(whose_is_g(igg)+1,1)+1 272 ig_send(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=ig 273 igg_found(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=igg 274 endif 275 276 else if( i==-mill_g(1,igg) .and. j==-mill_g(2,igg) .and. k==-mill_g(3,igg)) then 277 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 278 nfound_max=nfound_max+1 279 nfound_proc(whose_is_g(igg)+1,1)=nfound_proc(whose_is_g(igg)+1,1)+1 280 ig_send(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=ig 281 igg_found(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=-igg 282 endif 283 endif 284 285 enddo 286 enddo 287 288 do ig=gstart,ngw!loop on g vectors 289 ! first +g 290 i = -mill(1,ig) 291 j = -mill(2,ig) 292 k = -mill(3,ig) 293 if(ipol.eq.1) i=i+1 294 if(ipol.eq.2) j=j+1 295 if(ipol.eq.3) k=k+1 296 do igg=1,ngw_g 297 if( i==mill_g(1,igg) .and. j==mill_g(2,igg) .and. k==mill_g(3,igg)) then 298 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 299 nfound_max=nfound_max+1 300 nfound_proc(whose_is_g(igg)+1,2)=nfound_proc(whose_is_g(igg)+1,2)+1 301 ig_send(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=ig 302 igg_found(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=igg 303 endif 304 305 else if( i==-mill_g(1,igg) .and. j==-mill_g(2,igg) .and. k==-mill_g(3,igg)) then 306 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 307 nfound_max=nfound_max+1 308 nfound_proc(whose_is_g(igg)+1,2)=nfound_proc(whose_is_g(igg)+1,2)+1 309 ig_send(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=ig 310 igg_found(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=-igg 311 endif 312 endif 313 314 enddo 315 enddo 316 317! determine the largest nfound for processor and set it as dimension 318! for ctable_missing and ctable_missing_rev 319! copy ig_send to ctable_missing 320 321 call mp_sum(nfound_max, intra_bgrp_comm) 322 write(stdout,*) 'Additional found:', nfound_max 323 324 n_g_missing_p(ipol)=maxval(nfound_proc(:,:)) 325 326 call mp_max(n_g_missing_p(ipol), intra_bgrp_comm) 327 328 329 if(ipol==1) then 330 allocate(ctable_missing_1(n_g_missing_p(ipol),2,nproc_bgrp)) 331 ctable_missing_1(:,:,:)=0 332 do ip=1,nproc_bgrp 333 ctable_missing_1(1:nfound_proc(ip,1),1,ip)=ig_send(1:nfound_proc(ip,1),1,ip) 334 ctable_missing_1(1:nfound_proc(ip,2),2,ip)=ig_send(1:nfound_proc(ip,2),2,ip) 335 enddo 336 else 337 allocate(ctable_missing_2(n_g_missing_p(ipol),2,nproc_bgrp)) 338 ctable_missing_2(:,:,:)=0 339 do ip=1,nproc_bgrp 340 ctable_missing_2(1:nfound_proc(ip,1),1,ip)=ig_send(1:nfound_proc(ip,1),1,ip) 341 ctable_missing_2(1:nfound_proc(ip,2),2,ip)=ig_send(1:nfound_proc(ip,2),2,ip) 342 enddo 343 endif 344 345 346!mpi all to all for igg_found 347 348 allocate(igg_found_snd(n_g_missing_p(ipol),2,nproc_bgrp)) 349 allocate(igg_found_rcv(n_g_missing_p(ipol),2,nproc_bgrp)) 350 igg_found_snd(:,:,:)=0 351 do ip=1,nproc_bgrp 352 igg_found_snd(1:nfound_proc(ip,1),1,ip)=igg_found(1:nfound_proc(ip,1),1,ip) 353 igg_found_snd(1:nfound_proc(ip,2),2,ip)=igg_found(1:nfound_proc(ip,2),2,ip) 354 enddo 355 356 357 call mp_alltoall( igg_found_snd, igg_found_rcv, intra_bgrp_comm ) 358 359 if(ipol==1) then 360 allocate(ctable_missing_rev_1(n_g_missing_p(ipol),2,nproc_bgrp)) 361 ctable_missing_rev_1(:,:,:)=0 362 else 363 allocate(ctable_missing_rev_2(n_g_missing_p(ipol),2,nproc_bgrp)) 364 ctable_missing_rev_2(:,:,:)=0 365 endif 366 367 368 369 nfound_max=0 370 371 do ip=1,nproc_bgrp 372 do igg=1, n_g_missing_p(ipol) 373 if(igg_found_rcv(igg,1,ip) /= 0 ) then 374 found=.false. 375 do ig=1,ngw 376 if(igg_found_rcv(igg,1,ip)>0) then 377 if(ig_l2g(ig)==igg_found_rcv(igg,1,ip)) then 378 nfound_max=nfound_max+1 379 if(ipol==1) then 380 ctable_missing_rev_1(igg,1,ip)=ig 381 else 382 ctable_missing_rev_2(igg,1,ip)=ig 383 endif 384 found=.true. 385 endif 386 else 387 if(ig_l2g(ig)==-igg_found_rcv(igg,1,ip)) then 388 nfound_max=nfound_max+1 389 if(ipol==1) then 390 ctable_missing_rev_1(igg,1,ip)=-ig 391 else 392 ctable_missing_rev_2(igg,1,ip)=-ig 393 endif 394 found=.true. 395 endif 396 endif 397 enddo 398 if(.not.found) write(stdout,*) 'NOT FOUND:', igg_found_rcv(igg,1,ip) 399 endif 400 enddo 401 do igg=1, n_g_missing_p(ipol) 402 if(igg_found_rcv(igg,2,ip) /= 0 ) then 403 found=.false. 404 do ig=1,ngw 405 if(igg_found_rcv(igg,2,ip)>0) then 406 if(ig_l2g(ig)==igg_found_rcv(igg,2,ip)) then 407 nfound_max=nfound_max+1 408 if(ipol==1) then 409 ctable_missing_rev_1(igg,2,ip)=ig 410 else 411 ctable_missing_rev_2(igg,2,ip)=ig 412 endif 413 found=.true. 414 endif 415 else 416 if(ig_l2g(ig)==-igg_found_rcv(igg,2,ip)) then 417 nfound_max=nfound_max+1 418 if(ipol==1) then 419 ctable_missing_rev_1(igg,2,ip)=-ig 420 else 421 ctable_missing_rev_2(igg,2,ip)=-ig 422 endif 423 found=.true. 424 endif 425 endif 426 enddo 427 if(.not.found) write(stdout,*) 'NOT FOUND:', igg_found_rcv(igg,2,ip) 428 endif 429 enddo 430 431 enddo 432 call mp_sum(nfound_max, intra_bgrp_comm) 433 !write(stdout,*) 'Found check', nfound_max 434 deallocate(igg_found_snd,igg_found_rcv) 435 enddo 436 437 438 439 deallocate(igg_found, ig_send) 440return 441 442end subroutine gtable_missing 443 444 445 446 447subroutine gtable_missing_inv 448 449 USE efield_module, ONLY : ctabin_missing_1,ctabin_missing_2, whose_is_g,n_g_missing_m,& 450 & ctabin_missing_rev_1,ctabin_missing_rev_2 451 USE gvecw, ONLY : ngw, ngw_g 452 USE gvect, ONLY : ig_l2g, mill_g, mill, gstart 453 USE mp, ONLY : mp_sum, mp_max, mp_alltoall 454 USE io_global, ONLY : stdout 455 USE mp_global, ONLY : me_bgrp, nproc_bgrp, intra_bgrp_comm 456 USE parallel_include 457 458 459 implicit none 460 461 INTEGER :: ipol, i,j,k,ig,igg, nfound_max, ip 462 LOGICAL :: found 463 INTEGER :: nfound_proc(nproc_bgrp,2) 464 INTEGER, ALLOCATABLE :: igg_found(:,:,:), ig_send(:,:,:), igg_found_snd(:,:,:) 465 INTEGER, ALLOCATABLE :: igg_found_rcv(:,:,:) 466 INTEGER :: ierr,sndint,rcvint 467 468 469 470 allocate( igg_found(ngw_g,2,nproc_bgrp), ig_send(ngw_g,2,nproc_bgrp)) 471 do ipol=1,2 472 473 474 475 nfound_max=0 476 nfound_proc(:,:)=0 477 ig_send(:,:,:)=0 478 479 do ig=1,ngw!loop on g vectors 480 ! first +g 481 i = mill(1,ig) 482 j = mill(2,ig) 483 k = mill(3,ig) 484 if(ipol.eq.1) i=i+1 485 if(ipol.eq.2) j=j+1 486 if(ipol.eq.3) k=k+1 487 do igg=1,ngw_g 488 if( i==mill_g(1,igg) .and. j==mill_g(2,igg) .and. k==mill_g(3,igg)) then 489 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 490 nfound_max=nfound_max+1 491 nfound_proc(whose_is_g(igg)+1,1)=nfound_proc(whose_is_g(igg)+1,1)+1 492 ig_send(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=ig 493 igg_found(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=igg 494 endif 495 496 else if( i==-mill_g(1,igg) .and. j==-mill_g(2,igg) .and. k==-mill_g(3,igg)) then 497 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 498 nfound_max=nfound_max+1 499 nfound_proc(whose_is_g(igg)+1,1)=nfound_proc(whose_is_g(igg)+1,1)+1 500 ig_send(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=ig 501 igg_found(nfound_proc(whose_is_g(igg)+1,1),1,whose_is_g(igg)+1)=-igg 502 endif 503 endif 504 505 enddo 506 enddo 507 508 do ig=1,ngw!loop on g vectors 509 ! first +g 510 i = mill(1,ig) 511 j = mill(2,ig) 512 k = mill(3,ig) 513 if(ipol.eq.1) i=i-1 514 if(ipol.eq.2) j=j-1 515 if(ipol.eq.3) k=k-1 516 do igg=1,ngw_g 517 if( i==mill_g(1,igg) .and. j==mill_g(2,igg) .and. k==mill_g(3,igg)) then 518 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 519 nfound_max=nfound_max+1 520 nfound_proc(whose_is_g(igg)+1,2)=nfound_proc(whose_is_g(igg)+1,2)+1 521 ig_send(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=ig 522 igg_found(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=igg 523 endif 524 525 else if( i==-mill_g(1,igg) .and. j==-mill_g(2,igg) .and. k==-mill_g(3,igg)) then 526 if(whose_is_g(igg) /= -1 .and. whose_is_g(igg) /= me_bgrp) then 527 nfound_max=nfound_max+1 528 nfound_proc(whose_is_g(igg)+1,2)=nfound_proc(whose_is_g(igg)+1,2)+1 529 ig_send(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=ig 530 igg_found(nfound_proc(whose_is_g(igg)+1,2),2,whose_is_g(igg)+1)=-igg 531 endif 532 endif 533 534 enddo 535 enddo 536 537 538!determine the largest nfound for processor and set it as dimensione for ctabin_missing and ctabin_missing_rev 539!copy ig_send to ctabin_missing 540 541 call mp_sum(nfound_max, intra_bgrp_comm) 542 write(stdout,*) 'Additional found:', nfound_max 543 544 545 n_g_missing_m(ipol)=maxval(nfound_proc(:,:)) 546 call mp_max(n_g_missing_m(ipol), intra_bgrp_comm) 547 548 549 if(ipol==1) then 550 allocate(ctabin_missing_1(n_g_missing_m(ipol),2,nproc_bgrp)) 551 ctabin_missing_1(:,:,:)=0 552 do ip=1,nproc_bgrp 553 ctabin_missing_1(1:nfound_proc(ip,1),1,ip)=ig_send(1:nfound_proc(ip,1),1,ip) 554 ctabin_missing_1(1:nfound_proc(ip,2),2,ip)=ig_send(1:nfound_proc(ip,2),2,ip) 555 enddo 556 else 557 allocate(ctabin_missing_2(n_g_missing_m(ipol),2,nproc_bgrp)) 558 ctabin_missing_2(:,:,:)=0 559 do ip=1,nproc_bgrp 560 ctabin_missing_2(1:nfound_proc(ip,1),1,ip)=ig_send(1:nfound_proc(ip,1),1,ip) 561 ctabin_missing_2(1:nfound_proc(ip,2),2,ip)=ig_send(1:nfound_proc(ip,2),2,ip) 562 enddo 563 endif 564 565 566!mpi all to all for igg_found 567 568 allocate(igg_found_snd(n_g_missing_m(ipol),2,nproc_bgrp)) 569 allocate(igg_found_rcv(n_g_missing_m(ipol),2,nproc_bgrp)) 570 igg_found_snd(:,:,:)=0 571 do ip=1,nproc_bgrp 572 igg_found_snd(1:nfound_proc(ip,1),1,ip)=igg_found(1:nfound_proc(ip,1),1,ip) 573 igg_found_snd(1:nfound_proc(ip,2),2,ip)=igg_found(1:nfound_proc(ip,2),2,ip) 574 enddo 575 576 577 CALL mp_alltoall( igg_found_snd, igg_found_rcv, intra_bgrp_comm ) 578 579 if(ipol==1) then 580 allocate(ctabin_missing_rev_1(n_g_missing_m(ipol),2,nproc_bgrp)) 581 ctabin_missing_rev_1(:,:,:)=0 582 else 583 allocate(ctabin_missing_rev_2(n_g_missing_m(ipol),2,nproc_bgrp)) 584 ctabin_missing_rev_2(:,:,:)=0 585 endif 586 587 588 589 nfound_max=0 590 591 do ip=1,nproc_bgrp 592 do igg=1, n_g_missing_m(ipol) 593 if(igg_found_rcv(igg,1,ip) /= 0 ) then 594 found=.false. 595 do ig=1,ngw 596 if(igg_found_rcv(igg,1,ip)>0) then 597 if(ig_l2g(ig)==igg_found_rcv(igg,1,ip)) then 598 nfound_max=nfound_max+1 599 if(ipol==1) then 600 ctabin_missing_rev_1(igg,1,ip)=ig 601 else 602 ctabin_missing_rev_2(igg,1,ip)=ig 603 endif 604 found=.true. 605 endif 606 else 607 if(ig_l2g(ig)==-igg_found_rcv(igg,1,ip)) then 608 nfound_max=nfound_max+1 609 if(ipol==1) then 610 ctabin_missing_rev_1(igg,1,ip)=-ig 611 else 612 ctabin_missing_rev_2(igg,1,ip)=-ig 613 endif 614 found=.true. 615 endif 616 endif 617 enddo 618 if(.not.found) write(stdout,*) 'NOT FOUND:', igg_found_rcv(igg,1,ip) 619 endif 620 enddo 621 do igg=1, n_g_missing_m(ipol) 622 if(igg_found_rcv(igg,2,ip) /= 0 ) then 623 found=.false. 624 do ig=1,ngw 625 if(igg_found_rcv(igg,2,ip)>0) then 626 if(ig_l2g(ig)==igg_found_rcv(igg,2,ip)) then 627 nfound_max=nfound_max+1 628 if(ipol==1) then 629 ctabin_missing_rev_1(igg,2,ip)=ig 630 else 631 ctabin_missing_rev_2(igg,2,ip)=ig 632 endif 633 found=.true. 634 endif 635 else 636 if(ig_l2g(ig)==-igg_found_rcv(igg,2,ip)) then 637 nfound_max=nfound_max+1 638 if(ipol==1) then 639 ctabin_missing_rev_1(igg,2,ip)=-ig 640 else 641 ctabin_missing_rev_2(igg,2,ip)=-ig 642 endif 643 found=.true. 644 endif 645 endif 646 enddo 647 if(.not.found) write(stdout,*) 'NOT FOUND:', igg_found_rcv(igg,2,ip) 648 endif 649 enddo 650 651 enddo 652 call mp_sum(nfound_max, intra_bgrp_comm) 653 !write(stdout,*) 'Found check', nfound_max 654 deallocate(igg_found_snd,igg_found_rcv) 655 enddo 656 657 deallocate(igg_found, ig_send) 658 ! workaround by PG to avoid a large array like mill_g allocated all the time 659 deallocate ( mill_g ) 660 661return 662 663end subroutine gtable_missing_inv 664