1 subroutine argos_space_trvl(xw,vw,xwcr,gw,iwl,iwlp, 2 + xs,vs,gs,isl,islp, 3 + boxsiz,ibownr,ipl,ndx,itmp,rtmp,lenx) 4c 5 implicit none 6c 7#include "argos_space_common.fh" 8#include "global.fh" 9#include "msgids.fh" 10c 11 real*8 xw(mwm,3,mwa),vw(mwm,3,mwa),xwcr(mwm,3) 12 real*8 xs(msa,3),vs(msa,3) 13 real*8 gw(mwm,3,mwa),gs(msa,3) 14 integer iwl(mwm,miw2),iwlp(mwm,npackw) 15 integer isl(msa,mis2),islp(msa,npack) 16 integer lenx 17 real*8 boxsiz(maxbox,3) 18 integer ipl(mbox,mip2),ibownr(maxbox,3) 19 integer ndx(lenx),itmp(lenx) 20 real*8 rtmp(lenx) 21 logical lrec(27) 22c 23 integer i,indexw,indexs,j,k,ibx,iby,ibz,ipx,ipy,ipz 24 integer isbox,isnod,nrbox 25 integer ilp,ihp,jlp,jhp 26 integer il,ih,jl,jh,ilw,ihw,jlw,jhw,ils,ihs,jls,jhs 27 integer iliw,ihiw,jliw,jhiw 28 integer ilis,ihis,jlis,jhis 29 integer iwm,iwstay,jwstay,lwstay,nwgo,nwgosm 30 integer nwgtsm 31 integer nwnew,nwstay,iwmloc,jwmloc,lwmloc,irw 32 integer isa,jsa,isstay,jsstay,icsgm,ifsgm,ilsgm 33 integer nsnew,nsstay,lsstay,isaloc,jsaloc,lsaloc,irs 34 integer nsgo,jnode,iwfr,iwto,isfr,isto 35 real*8 factor,xscx,xscy,xscz,boxi(3) 36 integer itemps,nfold 37 logical lend 38 character*255 string 39c 40 boxi(1)=one/box(1) 41 boxi(2)=one/box(2) 42 boxi(3)=one/box(3) 43 nfold=0 44 lpbc9=.false. 45c 46 nwstay=0 47c 48c order the solvent molecules 49c 50 if(nwmloc.gt.0) then 51 do 1 i=1,nwmloc 52 ndx(i)=i 53 1 continue 54 endif 55 if(nwmloc.gt.1) then 56 lwmloc=nwmloc/2+1 57 irw=nwmloc 58 2 continue 59 if(lwmloc.gt.1) then 60 lwmloc=lwmloc-1 61 itemps=ndx(lwmloc) 62 else 63 itemps=ndx(irw) 64 ndx(irw)=ndx(1) 65 irw=irw-1 66 if(irw.eq.1) then 67 ndx(1)=itemps 68 goto 3 69 endif 70 endif 71 iwmloc=lwmloc 72 jwmloc=lwmloc+lwmloc 73 4 continue 74 if(jwmloc.le.irw) then 75 if(jwmloc.lt.irw) then 76 if((iwl(ndx(jwmloc),lwnod).eq.iwl(ndx(jwmloc+1),lwnod).and. 77 + iwl(ndx(jwmloc),lwbox).le.iwl(ndx(jwmloc+1),lwbox)).or. 78 + ((iwl(ndx(jwmloc),lwnod).eq.me.or. 79 + (iwl(ndx(jwmloc),lwnod).ne.me.and. 80 + iwl(ndx(jwmloc),lwnod).le.iwl(ndx(jwmloc+1),lwnod))).and. 81 + iwl(ndx(jwmloc+1),lwnod).ne.me)) jwmloc=jwmloc+1 82 endif 83 if((iwl(itemps,lwnod).eq.iwl(ndx(jwmloc),lwnod).and. 84 + iwl(itemps,lwbox).le.iwl(ndx(jwmloc),lwbox)).or. 85 + ((iwl(itemps,lwnod).eq.me.or. (iwl(itemps,lwnod).ne.me.and. 86 + iwl(itemps,lwnod).le.iwl(ndx(jwmloc),lwnod))).and. 87 + iwl(ndx(jwmloc),lwnod).ne.me)) then 88 ndx(iwmloc)=ndx(jwmloc) 89 iwmloc=jwmloc 90 jwmloc=jwmloc+jwmloc 91 else 92 jwmloc=irw+1 93 endif 94 goto 4 95 endif 96 ndx(iwmloc)=itemps 97 goto 2 98 3 continue 99c 100 do 5 k=1,3 101 do 8 i=1,nwmloc 102 rtmp(i)=xwcr(i,k) 103 8 continue 104 do 9 i=1,nwmloc 105 xwcr(i,k)=rtmp(ndx(i)) 106 9 continue 107 do 10 j=1,nwa 108 do 11 i=1,nwmloc 109 rtmp(i)=xw(i,k,j) 110 11 continue 111 do 12 i=1,nwmloc 112 xw(i,k,j)=rtmp(ndx(i)) 113 12 continue 114 do 13 i=1,nwmloc 115 rtmp(i)=vw(i,k,j) 116 13 continue 117 do 14 i=1,nwmloc 118 vw(i,k,j)=rtmp(ndx(i)) 119 14 continue 120 if(iguide.gt.0) then 121 do 113 i=1,nwmloc 122 rtmp(i)=gw(i,k,j) 123 113 continue 124 do 114 i=1,nwmloc 125 gw(i,k,j)=rtmp(ndx(i)) 126 114 continue 127 endif 128 10 continue 129 5 continue 130 do 18 k=1,miw2 131 do 19 i=1,nwmloc 132 itmp(i)=iwl(i,k) 133 19 continue 134 do 20 i=1,nwmloc 135 iwl(i,k)=itmp(ndx(i)) 136 20 continue 137 18 continue 138 endif 139c 140 if(nwmloc.gt.0) then 141 do 21 iwm=1,nwmloc 142 if(iwl(iwm,lwnod).eq.me) then 143 nwstay=iwm 144 else 145c 146c check if moving atoms go to neighboring processor 147c 148 do 222 k=1,27 149 if(iwl(iwm,lwnod).eq.neighb(k,1)) goto 223 150 222 continue 151 write(string,'(a,i4,a,i4,a,i4,9f6.2)') 152 + 'argos_space_travel: solvent molecule ', 153 + iwl(iwm,lwgmn),' moving to non-neighbor ',iwl(iwm,lwnod), 154 + ' from ',me,((xw(iwm,i,j),i=1,3),j=1,3) 155 call md_abort(string,me) 156 223 continue 157c 158 endif 159c 160c testcode 161c 162 if(iand(idebug,8).eq.8) then 163 if(iwl(iwm,lwnod).ne.me) write(lfndbg,'(a,3i5)') 164 + 'Travel w fnd ',me,iwl(iwm,lwnod),iwl(iwm,lwgmn) 165 endif 166c 167c end test code 168c 169 21 continue 170 endif 171c 172c order the solute atoms 173c 174c isl(isa,lsbox) : box 175c isl(isa,lsnod) : node 176c isl(isa,lssgm) : segment 177c 178 nsstay=0 179 if(nsaloc.gt.0) then 180 do 22 i=1,nsaloc 181 ndx(i)=i 182 22 continue 183 endif 184c 185 if(nsaloc.gt.1) then 186 lsaloc=nsaloc/2+1 187 irs=nsaloc 188 23 continue 189 if(lsaloc.gt.1) then 190 lsaloc=lsaloc-1 191 itemps=ndx(lsaloc) 192 else 193 itemps=ndx(irs) 194 ndx(irs)=ndx(1) 195 irs=irs-1 196 if(irs.eq.1) then 197 ndx(1)=itemps 198 goto 24 199 endif 200 endif 201 isaloc=lsaloc 202 jsaloc=lsaloc+lsaloc 203 25 continue 204 if(jsaloc.le.irs) then 205 if(jsaloc.lt.irs) then 206 if((isl(ndx(jsaloc),lsnod).eq.isl(ndx(jsaloc+1),lsnod).and. 207 + (isl(ndx(jsaloc),lsbox).lt.isl(ndx(jsaloc+1),lsbox).or. 208 + (isl(ndx(jsaloc),lsbox).eq.isl(ndx(jsaloc+1),lsbox).and. 209 + isl(ndx(jsaloc),lssgm).le.isl(ndx(jsaloc+1),lssgm)))).or. 210 + ((isl(ndx(jsaloc),lsnod).eq.me.or. 211 + (isl(ndx(jsaloc),lsnod).ne.me.and. 212 + isl(ndx(jsaloc),lsnod).le.isl(ndx(jsaloc+1),lsnod))).and. 213 + isl(ndx(jsaloc+1),lsnod).ne.me)) jsaloc=jsaloc+1 214 endif 215 if((isl(itemps,lsnod).eq.isl(ndx(jsaloc),lsnod).and. 216 + (isl(itemps,lsbox).lt.isl(ndx(jsaloc),lsbox).or. 217 + (isl(itemps,lsbox).eq.isl(ndx(jsaloc),lsbox).and. 218 + isl(itemps,lssgm).le.isl(ndx(jsaloc),lssgm)))).or. 219 + ((isl(itemps,lsnod).eq.me.or. (isl(itemps,lsnod).ne.me.and. 220 + isl(itemps,lsnod).le.isl(ndx(jsaloc),lsnod))).and. 221 + isl(ndx(jsaloc),lsnod).ne.me)) then 222 ndx(isaloc)=ndx(jsaloc) 223 isaloc=jsaloc 224 jsaloc=jsaloc+jsaloc 225 else 226 jsaloc=irs+1 227 endif 228 goto 25 229 endif 230 ndx(isaloc)=itemps 231 goto 23 232 24 continue 233c 234 do 26 k=1,3 235 do 27 i=1,nsaloc 236 rtmp(i)=xs(i,k) 237 27 continue 238 do 28 i=1,nsaloc 239 xs(i,k)=rtmp(ndx(i)) 240 28 continue 241 do 29 i=1,nsaloc 242 rtmp(i)=vs(i,k) 243 29 continue 244 do 30 i=1,nsaloc 245 vs(i,k)=rtmp(ndx(i)) 246 30 continue 247 if(iguide.gt.0) then 248 do 2129 i=1,nsaloc 249 rtmp(i)=gs(i,k) 250 2129 continue 251 do 2130 i=1,nsaloc 252 gs(i,k)=rtmp(ndx(i)) 253 2130 continue 254 endif 255 26 continue 256 do 40 k=1,mis2 257 do 41 i=1,nsaloc 258 itmp(i)=isl(i,k) 259 41 continue 260 do 42 i=1,nsaloc 261 isl(i,k)=itmp(ndx(i)) 262 42 continue 263 40 continue 264 endif 265c 266 if(nsa.gt.0) then 267 do 43 isa=1,nsaloc 268 if(isl(isa,lsnod).eq.me) then 269 nsstay=isa 270 else 271c 272c check if moving atoms go to neighboring processor 273c 274 do 444 k=1,27 275 if(isl(isa,lsnod).eq.neighb(k,1)) goto 445 276 444 continue 277 write(string,'(a,i4,a,i4,a,i4,3f6.2)') 278 + 'argos_space_travel: solute segment ', 279 + isl(isa,lssgm),' moving to non-neighbor ',isl(isa,lsnod), 280 + ' from ',me,(xs(isa,i),i=1,3) 281 call md_abort(string,me) 282 445 continue 283c 284 endif 285 43 continue 286 endif 287c 288c make packages ready for shipment 289c 290c loop over all neighboring nodes 291c 292 call ga_distribution(ga_iw,me,iliw,ihiw,jliw,jhiw) 293 call ga_distribution(ga_w,me,ilw,ihw,jlw,jhw) 294 call ga_distribution(ga_is,me,ilis,ihis,jlis,jhis) 295 call ga_distribution(ga_s,me,ils,ihs,jls,jhs) 296c 297 indexw=0 298 indexs=0 299 nwgosm=0 300c 301 do 70 i=1,27 302 jnode=neighb(i,1) 303 if(jnode.ge.0.and.jnode.ne.me) then 304c 305c for the solvent 306c 307 iwfr=0 308 iwto=0 309 do 71 iwm=nwstay+1,nwmloc 310 if(iwl(iwm,lwnod).eq.jnode) then 311 if(iwfr.eq.0) iwfr=iwm 312 iwto=iwm 313c 314c testcode 315c 316 if(iand(idebug,8).eq.8) then 317 if(iwl(iwm,lwnod).ne.me) write(lfndbg,'(a,3i5)') 318 + 'Travel w snd ',me,iwl(iwm,lwnod),iwl(iwm,lwgmn) 319 endif 320c 321c end test code 322c 323 endif 324 71 continue 325c 326c if molecules need to travel copy coordinates etc into global array 327c 328 nwgo=iwto-iwfr+1 329 if(iwfr.eq.0) nwgo=0 330 ipl(1,1)=0 331 ipl(1,2)=0 332c 333 if(nwgo.gt.0) then 334 nwgosm=nwgosm+nwgo 335 il=iliw+indexw 336 ih=il+nwgo-1 337 if(npackw.eq.0) then 338 call ga_put(ga_iw,il,ih,jliw,jhiw,iwl(iwfr,1),mwm) 339 else 340 call argos_space_packw(ih-il+1,iwl(iwfr,1),iwlp(iwfr,1)) 341 call ga_put(ga_iw,il,ih,jliw,jliw+npackw-1,iwlp(iwfr,1),mwm) 342 endif 343 il=ilw+indexw 344 ih=il+nwgo-1 345 call ga_put(ga_w,il,ih,jlw,jlw+3*mwa-1,xw(iwfr,1,1),mwm) 346 call ga_put(ga_w,il,ih,jlw+3*mwa,jlw+6*mwa-1,vw(iwfr,1,1),mwm) 347 call ga_put(ga_w,il,ih,jlw+6*mwa,jlw+6*mwa+2,xwcr(iwfr,1),mwm) 348 if(iguide.gt.0) then 349 call ga_put(ga_w,il,ih,jlw+6*mwa+3,jlw+9*mwa+2,gw(iwfr,1,1),mwm) 350 endif 351 ipl(1,1)=indexw+1 352 ipl(1,2)=indexw+nwgo 353 indexw=indexw+nwgo 354 endif 355c 356c for the solute 357c 358 isfr=0 359 isto=0 360 do 72 isa=nsstay+1,nsaloc 361 if(isl(isa,lsnod).eq.jnode) then 362 if(isfr.eq.0) isfr=isa 363 isto=isa 364 endif 365 72 continue 366 nsgo=isto-isfr+1 367 if(isfr.eq.0) nsgo=0 368 ipl(1,3)=0 369 ipl(1,4)=0 370 if(nsgo.gt.0) then 371 il=ilis+indexs 372 ih=il+nsgo-1 373 if(npack.eq.0) then 374 call ga_put(ga_is,il,ih,jlis,jhis,isl(isfr,1),msa) 375 else 376 call argos_space_pack(ih-il+1,isl(isfr,1),islp(isfr,1)) 377 call ga_put(ga_is,il,ih,jlis,jlis+npack-1,islp(isfr,1),msa) 378 endif 379 call ga_put(ga_s,il,ih,jls,jls+2,xs(isfr,1),msa) 380 call ga_put(ga_s,il,ih,jls+3,jls+5,vs(isfr,1),msa) 381 if(iguide.gt.0) then 382 call ga_put(ga_s,il,ih,jls+6,jls+8,gs(isfr,1),msa) 383 endif 384 ipl(1,3)=indexs+1 385 ipl(1,4)=indexs+nsgo 386 indexs=indexs+nsgo 387 endif 388c 389c inform other node of number of molecules to get 390c 391 if(ipl(1,1).gt.0.or.ipl(1,3).gt.0) then 392 call ga_distribution(ga_ip,jnode,ilp,ihp,jlp,jhp) 393 ilp=ilp+2+i 394 call ga_put(ga_ip,ilp,ilp,jlp,jhp,ipl,mbox) 395 endif 396 endif 397 lrec(i)=.false. 398 70 continue 399c 400 call ga_sync() 401c 402c receive molecules from other nodes 403c 404 nwgtsm=0 405c 406 call ga_distribution(ga_ip,me,ilp,ihp,jlp,jhp) 407 call ga_get(ga_ip,ilp,ilp+30,jlp,jhp,ipl,mbox) 408c 409 do 74 i=1,27 410 jnode=neighb(i,2) 411 if(jnode.ge.0.and.jnode.ne.me.and..not.lrec(i)) then 412c 413 iwfr=ipl(3+i,1) 414 iwto=ipl(3+i,2) 415 isfr=ipl(3+i,3) 416 isto=ipl(3+i,4) 417c 418 nwnew=iwto-iwfr+1 419 nsnew=isto-isfr+1 420c 421 if(iwfr.eq.0) nwnew=0 422 if(isfr.eq.0) nsnew=0 423c 424 if(nwstay+nwnew.gt.mwm) then 425 write(string,'(a,i7,a,i7)') 426 + 'Travel: mwm needs increase with ',nwnew,' to ',nwstay+nwnew 427 call md_abort(string,me) 428 endif 429 if(nsstay+nsnew.gt.msa) then 430 write(string,'(a,i7,a,i7)') 431 + 'Travel: msa needs increase with ',nsnew,' to ',nsstay+nsnew 432 call md_abort(string,me) 433 endif 434c 435 lrec(i)=.true. 436c 437 if(iwfr.gt.0) then 438 nwgtsm=nwgtsm+nwnew 439 iwto=ipl(3+i,2) 440 call ga_distribution(ga_iw,jnode,iliw,ihiw,jliw,jhiw) 441 call ga_distribution(ga_w,jnode,ilw,ihw,jlw,jhw) 442c 443c get data for additional molecules 444c 445 il=iliw+iwfr-1 446 ih=iliw+iwto-1 447 if(npackw.eq.0) then 448 call ga_get(ga_iw,il,ih,jliw,jhiw,iwl(nwstay+1,1),mwm) 449 else 450 call ga_get(ga_iw,il,ih,jliw,jliw+npackw-1,iwlp(nwstay+1,1),mwm) 451 call argos_space_unpackw(ih-il+1,iwl(nwstay+1,1),iwlp(nwstay+1,1)) 452 endif 453 call ga_get(ga_w,il,ih,jlw,jlw+3*mwa-1,xw(nwstay+1,1,1),mwm) 454 call ga_get(ga_w,il,ih,jlw+3*mwa,jlw+6*mwa-1,vw(nwstay+1,1,1),mwm) 455 call ga_get(ga_w,il,ih,jlw+6*mwa,jlw+6*mwa+2,xwcr(nwstay+1,1),mwm) 456 if(iguide.gt.0) then 457 call ga_get(ga_w,il,ih,jlw+6*mwa+3,jlw+9*mwa+2, 458 + gw(nwstay+1,1,1),mwm) 459 endif 460c 461c testcode 462c 463 if(iand(idebug,8).eq.8) then 464 write(lfndbg,'(a,3i5)') 465 + ('Travel w rcv ',me,jnode,iwl(nwstay+k,lwgmn),k=1,nwnew) 466 endif 467c 468c end test code 469c 470c 471c update number of local solvent molecules 472c 473 nwstay=nwstay+nwnew 474c 475 endif 476c 477c for the solute 478c 479 if(isfr.gt.0) then 480 call ga_distribution(ga_is,jnode,ilis,ihis,jlis,jhis) 481 call ga_distribution(ga_s,jnode,ils,ihs,jls,jhs) 482 il=ilis+isfr-1 483 ih=ilis+isto-1 484 jl=jlis 485 jh=jhis 486 if(npack.eq.0) then 487 call ga_get(ga_is,il,ih,jlis,jhis,isl(nsstay+1,1),msa) 488 else 489 call ga_get(ga_is,il,ih,jlis,jlis+npack-1,islp(nsstay+1,1),msa) 490 call argos_space_unpack(ih-il+1,isl(nsstay+1,1),islp(nsstay+1,1)) 491 endif 492 call ga_get(ga_s,il,ih,jls,jls+2,xs(nsstay+1,1),msa) 493 call ga_get(ga_s,il,ih,jls+3,jls+5,vs(nsstay+1,1),msa) 494 if(iguide.gt.0) then 495 call ga_get(ga_s,il,ih,jls+6,jls+8,gs(nsstay+1,1),msa) 496 endif 497c 498 nsstay=nsstay+nsnew 499 endif 500c 501 endif 502c 503c reset the pointers to zero 504c 505 ipl(3+i,1)=0 506 ipl(3+i,2)=0 507 ipl(3+i,3)=0 508 ipl(3+i,4)=0 509c 510 74 continue 511c 512c reset ipl in global array 513c 514 call ga_put(ga_ip,ilp,ilp+30,jlp,jhp,ipl,mbox) 515c 516c order the solvent molecules according to subbox and 517c store indices into ip 518c 519c ip(1,1) : number of boxes on this node 520c ip(1,2) : number of solvent molecules on this node 521c ip(2,2) : number of solute atoms on this node 522c 523c ip(3+i,1) : index for solvents to be moved to the i-th neighbor 524c 525c ip(30+i,1) : number of i-th box on this node 526c ip(30+i,2) : index to first solvent in i-th box 527c ip(30+i,3) : index to lasst solvent in i-th box 528c 529 if(nwstay.gt.0.and.(nwgosm.gt.0.or.nwgtsm.gt.0)) then 530 do 81 i=1,nwstay 531 ndx(i)=i 532 81 continue 533 if(nwstay.gt.1) then 534 lwstay=nwstay/2+1 535 irw=nwstay 536 82 continue 537 if(lwstay.gt.1) then 538 lwstay=lwstay-1 539 itemps=ndx(lwstay) 540 else 541 itemps=ndx(irw) 542 ndx(irw)=ndx(1) 543 irw=irw-1 544 if(irw.eq.1) then 545 ndx(1)=itemps 546 goto 83 547 endif 548 endif 549 iwstay=lwstay 550 jwstay=lwstay+lwstay 551 84 continue 552 if(jwstay.le.irw) then 553 if(jwstay.lt.irw) then 554 if(iwl(ndx(jwstay),lwbox).le.iwl(ndx(jwstay+1),lwbox)) 555 + jwstay=jwstay+1 556 endif 557 if(iwl(itemps,lwbox).le.iwl(ndx(jwstay),lwbox)) then 558 ndx(iwstay)=ndx(jwstay) 559 iwstay=jwstay 560 jwstay=jwstay+jwstay 561 else 562 jwstay=irw+1 563 endif 564 goto 84 565 endif 566 ndx(iwstay)=itemps 567 goto 82 568 83 continue 569c 570 do 85 k=1,3 571 do 88 i=1,nwstay 572 rtmp(i)=xwcr(i,k) 573 88 continue 574 do 89 i=1,nwstay 575 xwcr(i,k)=rtmp(ndx(i)) 576 89 continue 577 do 90 j=1,nwa 578 do 91 i=1,nwstay 579 rtmp(i)=xw(i,k,j) 580 91 continue 581 do 92 i=1,nwstay 582 xw(i,k,j)=rtmp(ndx(i)) 583 92 continue 584 do 93 i=1,nwstay 585 rtmp(i)=vw(i,k,j) 586 93 continue 587 do 94 i=1,nwstay 588 vw(i,k,j)=rtmp(ndx(i)) 589 94 continue 590 if(iguide.gt.0) then 591 do 193 i=1,nwstay 592 rtmp(i)=gw(i,k,j) 593 193 continue 594 do 194 i=1,nwstay 595 gw(i,k,j)=rtmp(ndx(i)) 596 194 continue 597 endif 598 90 continue 599 85 continue 600 do 98 k=1,miw2 601 do 99 i=1,nwstay 602 itmp(i)=iwl(i,k) 603 99 continue 604 do 100 i=1,nwstay 605 iwl(i,k)=itmp(ndx(i)) 606 100 continue 607 98 continue 608c 609 endif 610 endif 611c 612c order the solute according to segment 613c 614 if(nsstay.gt.0) then 615 do 122 i=1,nsstay 616 ndx(i)=i 617 122 continue 618 if(nsstay.gt.1) then 619 lsstay=nsstay/2+1 620 irs=nsstay 621 123 continue 622 if(lsstay.gt.1) then 623 lsstay=lsstay-1 624 itemps=ndx(lsstay) 625 else 626 itemps=ndx(irs) 627 ndx(irs)=ndx(1) 628 irs=irs-1 629 if(irs.eq.1) then 630 ndx(1)=itemps 631 goto 124 632 endif 633 endif 634 isstay=lsstay 635 jsstay=lsstay+lsstay 636 125 continue 637 if(jsstay.le.irs) then 638 if(jsstay.lt.irs) then 639 if(isl(ndx(jsstay),lssgm).le.isl(ndx(jsstay+1),lssgm)) 640 + jsstay=jsstay+1 641 endif 642 if(isl(itemps,lssgm).le.isl(ndx(jsstay),lssgm)) then 643 ndx(isstay)=ndx(jsstay) 644 isstay=jsstay 645 jsstay=jsstay+jsstay 646 else 647 jsstay=irs+1 648 endif 649 goto 125 650 endif 651 ndx(isstay)=itemps 652 goto 123 653 124 continue 654 endif 655c 656c for each segment : 1. determine box number 657c 2. assign box number to each atom 658c 3. when box not owned by node: 659c a. assign box number 660c b. assign correct node number 661c 662 goto 666 663 icsgm=isl(ndx(1),lssgm) 664 ifsgm=1 665 ilsgm=1 666 do 126 isa=2,nsstay+1 667c 668c if isa is first atom of a new segment or very last atom 669c 670 671 if(isa.le.nsstay) then 672 lend=isl(ndx(isa),lssgm).ne.icsgm 673 else 674 lend=.true. 675 endif 676 if(lend) then 677 if(isa.gt.nsstay) ilsgm=nsstay 678 if(ifsgm.gt.0.and.ilsgm.ge.ifsgm) then 679 xscx=zero 680 xscy=zero 681 xscz=zero 682 do 127 jsa=ifsgm,ilsgm 683 xscx=xscx+xs(ndx(jsa),1) 684 xscy=xscy+xs(ndx(jsa),2) 685 xscz=xscz+xs(ndx(jsa),3) 686 127 continue 687 factor=one/dble(ilsgm-ifsgm+1) 688 xscx=factor*xscx 689 xscy=factor*xscy 690 xscz=factor*xscz 691 if(npbtyp.ne.0) then 692 if(abs(xscx).gt.boxh(1)) then 693 xscx=xscx-nint(xscx*boxi(1))*box(1) 694 nfold=1 695 endif 696 if(abs(xscy).gt.boxh(2)) then 697 xscy=xscy-nint(xscy*boxi(2))*box(2) 698 nfold=1 699 endif 700 if(abs(xscz).gt.boxh(3)) then 701 xscz=xscz-nint(xscz*boxi(3))*box(3) 702 nfold=1 703 endif 704 endif 705c 706c determine the box number 707c 708 ibx=0 709 iby=0 710 ibz=0 711 do 128 i=1,nbx-1 712 if(xscx+boxh(1).gt.boxsiz(i,1)) ibx=i 713 128 continue 714 do 129 i=1,nby-1 715 if(xscy+boxh(2).gt.boxsiz(i,2)) iby=i 716 129 continue 717 do 1130 i=1,nbz-1 718 if(xscz+boxh(3).gt.boxsiz(i,3)) ibz=i 719 1130 continue 720 if(npbtyp.gt.0) then 721 if(ibx.ge.nbx) ibx=ibx-nbx 722 if(iby.ge.nby) iby=iby-nby 723 if(ibx.lt.0) ibx=ibx+nbx 724 if(iby.lt.0) iby=iby+nby 725 if(npbtyp.eq.1) then 726 if(ibz.ge.nbz) ibz=ibz-nbz 727 if(ibz.lt.0) ibz=ibz+nbz 728 else 729 if(ibz.ge.nbz) ibz=nbz-1 730 if(ibz.lt.0) ibz=0 731 endif 732 else 733 if(ibx.ge.nbx) ibx=nbx-1 734 if(iby.ge.nby) iby=nby-1 735 if(ibz.ge.nbz) ibz=nbz-1 736 if(ibx.lt.0) ibx=0 737 if(iby.lt.0) iby=0 738 if(ibz.lt.0) ibz=0 739 endif 740 ipx=ibownr(ibx+1,1) 741 ipy=ibownr(iby+1,2) 742 ipz=ibownr(ibz+1,3) 743 isbox=(ibz*nby+iby)*nbx+ibx 744 isnod=(ipz*npy+ipy)*npx+ipx 745c 746c assign box and node numbers 747c 748 do 1131 jsa=ifsgm,ilsgm 749 isl(ndx(jsa),lsbox)=isbox 750 isl(ndx(jsa),lsnod)=isnod 751 1131 continue 752c 753 endif 754 if(isa.le.nsstay) icsgm=isl(ndx(isa),lssgm) 755 ifsgm=isa 756 else 757 ilsgm=isa 758 endif 759 126 continue 760 666 continue 761c 762c order solute according to box, segment, charge group, atom number 763c 764 if(nsstay.gt.1) then 765 lsstay=nsstay/2+1 766 irs=nsstay 767 132 continue 768 if(lsstay.gt.1) then 769 lsstay=lsstay-1 770 itemps=ndx(lsstay) 771 else 772 itemps=ndx(irs) 773 ndx(irs)=ndx(1) 774 irs=irs-1 775 if(irs.eq.1) then 776 ndx(1)=itemps 777 goto 133 778 endif 779 endif 780 isstay=lsstay 781 jsstay=lsstay+lsstay 782 134 continue 783 if(jsstay.le.irs) then 784 if(jsstay.lt.irs) then 785 if(isl(ndx(jsstay),lsbox).lt.isl(ndx(jsstay+1),lsbox).or. 786 + (isl(ndx(jsstay),lsbox).eq.isl(ndx(jsstay+1),lsbox).and. 787 + (isl(ndx(jsstay),lssgm).lt.isl(ndx(jsstay+1),lssgm).or. 788 + (isl(ndx(jsstay),lssgm).eq.isl(ndx(jsstay+1),lssgm).and. 789 + (isl(ndx(jsstay),lsgrp).lt.isl(ndx(jsstay+1),lsgrp).or. 790 + (isl(ndx(jsstay),lsgrp).eq.isl(ndx(jsstay+1),lsgrp).and. 791 + isl(ndx(jsstay),lsgan).le.isl(ndx(jsstay+1),lsgan))))))) 792 + jsstay=jsstay+1 793 endif 794 if(isl(itemps,lsbox).lt.isl(ndx(jsstay),lsbox).or. 795 + (isl(itemps,lsbox).eq.isl(ndx(jsstay),lsbox).and. 796 + (isl(itemps,lssgm).lt.isl(ndx(jsstay),lssgm).or. 797 + (isl(itemps,lssgm).eq.isl(ndx(jsstay),lssgm).and. 798 + (isl(itemps,lsgrp).lt.isl(ndx(jsstay),lsgrp).or. 799 + (isl(itemps,lsgrp).eq.isl(ndx(jsstay),lsgrp).and. 800 + isl(itemps,lsgan).le.isl(ndx(jsstay),lsgan))))))) then 801 ndx(isstay)=ndx(jsstay) 802 isstay=jsstay 803 jsstay=jsstay+jsstay 804 else 805 jsstay=irs+1 806 endif 807 goto 134 808 endif 809 ndx(isstay)=itemps 810 goto 132 811 133 continue 812 endif 813c 814 do 135 k=1,3 815 do 136 i=1,nsstay 816 rtmp(i)=xs(i,k) 817 136 continue 818 do 137 i=1,nsstay 819 xs(i,k)=rtmp(ndx(i)) 820 137 continue 821 do 138 i=1,nsstay 822 rtmp(i)=vs(i,k) 823 138 continue 824 do 139 i=1,nsstay 825 vs(i,k)=rtmp(ndx(i)) 826 139 continue 827 if(iguide.gt.0) then 828 do 1138 i=1,nsstay 829 rtmp(i)=gs(i,k) 830 1138 continue 831 do 1139 i=1,nsstay 832 gs(i,k)=rtmp(ndx(i)) 833 1139 continue 834 endif 835 135 continue 836 do 149 k=1,mis2 837 do 150 i=1,nsstay 838 itmp(i)=isl(i,k) 839 150 continue 840 do 151 i=1,nsstay 841 isl(i,k)=itmp(ndx(i)) 842 151 continue 843 149 continue 844c 845 endif 846c 847 do 200 i=1,ipl(1,1) 848 ipl(30+i,2)=0 849 ipl(30+i,3)=0 850 ipl(30+i,4)=0 851 ipl(30+i,5)=0 852 200 continue 853c 854 do 201 i=1,ipl(1,1) 855 nrbox=ipl(30+i,1) 856 if(nwstay.gt.0) then 857 do 202 iwm=1,nwstay 858 if(iwl(iwm,lwbox).eq.nrbox) then 859 if(ipl(30+i,2).eq.0) ipl(30+i,2)=iwm 860 ipl(30+i,3)=iwm 861 endif 862 202 continue 863 endif 864 if(nsstay.gt.0) then 865 do 203 isa=1,nsstay 866 if(isl(isa,lsbox).eq.nrbox) then 867 if(ipl(30+i,4).eq.0) ipl(30+i,4)=isa 868 ipl(30+i,5)=isa 869 endif 870 203 continue 871 endif 872 201 continue 873c 874 nwmloc=nwstay 875 ipl(1,2)=nwmloc 876 nsaloc=nsstay 877 ipl(2,2)=nsaloc 878c 879 call ga_igop(msp_23,nfold,1,'+') 880 lpbc9=nfold.gt.0 881c 882 return 883 end 884c $Id$ 885