1* 2* $Id$ 3* 4 5* *********************************** 6* * * 7* * Cram_Init * 8* * * 9* *********************************** 10 11 subroutine Cram_Init() 12 implicit none 13#include "errquit.fh" 14 15#include "bafdecls.fh" 16#include "cram_common.fh" 17 18 19* **** local variables **** 20 integer MASTER,taskid 21 parameter (MASTER=0) 22 23 logical value 24 integer nfft3d 25 integer i,j,k 26 integer k1,k2,k3 27 integer q,p,indx 28 integer nb 29 integer nx,ny,nz 30 integer nxh,nyh,nzh 31 integer pack(2),nidb 32 real*8 kv(3),ecut 33 34* **** external functions **** 35 logical cloak_masker 36 integer brillioun_nbrillq 37 real*8 lattice_ecut,lattice_wcut,brillioun_k 38 external cloak_masker 39 external brillioun_nbrillq 40 external lattice_ecut,lattice_wcut,brillioun_k 41 42 maxsize = brillioun_nbrillq()+1 43 44 value = BA_alloc_get(mt_int,maxsize, 45 > 'nidb_list', 46 > nidb_list(2), 47 > nidb_list(1)) 48 value = value.and. 49 > BA_alloc_get(mt_int,maxsize, 50 > 'nidb2_list', 51 > nidb2_list(2), 52 > nidb2_list(1)) 53 value = value.and. 54 > BA_alloc_get(mt_int,maxsize, 55 > 'nwaveall_list', 56 > nwaveall_list(2), 57 > nwaveall_list(1)) 58 value = value.and. 59 > BA_alloc_get(mt_int,2*maxsize, 60 > 'pack_list', 61 > pack_list(2), 62 > pack_list(1)) 63 if (.not. value) 64 > call errquit('Cram_init: out of heap memory',0,MA_ERR) 65 66 67 call cloak_init() 68 call Parallel3d_taskid_i(taskid) 69 call C3dB_nfft3d(1,nfft3d) 70 call C3dB_nx(1,nx) 71 call C3dB_ny(1,ny) 72 call C3dB_nz(1,nz) 73 nxh = nx/2 74 nyh = ny/2 75 nzh = nz/2 76 77 do nb=0,maxsize-1 78 value = BA_alloc_get(mt_int,nfft3d, 79 > 'pack',pack(2),pack(1)) 80 if (.not. value) 81 > call errquit('Cram_init: out of heap memory',0, MA_ERR) 82 83 int_mb(pack_list(1)+ 2*nb) = pack(1) 84 int_mb(pack_list(1)+ 2*nb+1) = pack(2) 85 86 if (nb.eq.0) then 87 ecut = lattice_ecut() 88 kv(1) = 0.0d0 89 kv(2) = 0.0d0 90 kv(3) = 0.0d0 91 else 92 ecut = lattice_wcut() 93 kv(1) = brillioun_k(1,nb) 94 kv(2) = brillioun_k(2,nb) 95 kv(3) = brillioun_k(3,nb) 96 !kv(1) = 0.0d0 97 !kv(2) = 0.0d0 98 !kv(3) = 0.0d0 99 end if 100 call cloak_set(kv,ecut) 101 102 nidb = 0 103 104* **** k=(k1,k2,k3) **** 105 do k=(-nzh+1),(nzh-1) 106 do j=(-nyh+1),(nyh-1) 107 do i=(-nxh+1),(nxh-1) 108 k1=i 109 k2=j 110 k3=k 111 if (k1.lt.0) k1 = k1 + nx 112 if (k2.lt.0) k2 = k2 + ny 113 if (k3.lt.0) k3 = k3 + nz 114 115 !call C3dB_ktoqp(1,k3+1,q,p) 116 call C3dB_ijktoindexp(1,k1+1,k2+1,k3+1,indx,p) 117 118 if (p.eq.taskid) then 119 !indx = (q-1)*(nx)*ny + k2*(nx) + k1 + 1 120 if (.not.cloak_masker(indx)) then 121 nidb = nidb + 1 122 int_mb(pack(1)+nidb-1) = indx 123 end if 124 end if 125 end do 126 end do 127 end do 128 129 int_mb(nidb_list(1)+nb) = nidb 130 end do 131 132 call c_Balance_Init(maxsize, 133 > int_mb(nidb_list(1)), 134 > int_mb(nidb2_list(1))) 135 136 do nb=0,maxsize-1 137 int_mb(nwaveall_list(1)+nb) = int_mb(nidb_list(1)+nb) 138 call C3dB_ISumAll(int_mb(nwaveall_list(1)+nb)) 139 end do 140 141* **** set npack_max - note that npack_max **** 142* **** can be associated with different k-points **** 143* **** on different processors **** 144 npack_max = 0 145 do nb=1,maxsize-1 146 if (int_mb(nidb2_list(1)+nb).gt.npack_max) then 147 npack_max = int_mb(nidb2_list(1)+nb) 148 end if 149 end do 150 151 return 152 end 153 154 155* *********************************** 156* * * 157* * Cram_end * 158* * * 159* *********************************** 160 161 subroutine Cram_end() 162 implicit none 163#include "errquit.fh" 164 165#include "bafdecls.fh" 166#include "cram_common.fh" 167 168* **** local variables **** 169 logical value 170 integer dum2,nb 171 172 value = .true. 173 do nb=0,maxsize-1 174 dum2 = int_mb(pack_list(1)+2*nb+1) 175 value = value.and.BA_free_heap(dum2) 176 end do 177 value = value.and.BA_free_heap(nidb_list(2)) 178 value = value.and.BA_free_heap(nidb2_list(2)) 179 value = value.and.BA_free_heap(nwaveall_list(2)) 180 value = value.and.BA_free_heap(pack_list(2)) 181 182 if (.not. value) 183 > call errquit('Cram_end:error freeing heap memory',0, MA_ERR) 184 185 call c_Balance_End() 186 call cloak_end() 187 return 188 end 189 190 191 192* *********************************** 193* * * 194* * Cram_c_pack * 195* * * 196* *********************************** 197 198 subroutine Cram_c_pack(nb,A) 199 implicit none 200 integer nb 201 complex*16 A(*) 202 203#include "bafdecls.fh" 204#include "errquit.fh" 205#include "cram_common.fh" 206 207 208* **** local variables **** 209 logical value 210 integer i,nfft3d 211 integer tmp1(2) 212 integer nidb,pack 213 214 call nwpw_timing_start(9) 215 nidb = int_mb(nidb_list(1)+nb) 216 pack = int_mb(pack_list(1)+2*nb) 217 218 call C3dB_nfft3d(1,nfft3d) 219 value = BA_push_get(mt_dcpl,nfft3d,'tmp1',tmp1(2),tmp1(1)) 220 if (.not.value) call errquit('out of stack memory',0, MA_ERR) 221 222 call dcopy(2*nfft3d,A,1,dcpl_mb(tmp1(1)),1) 223 call dcopy(2*nfft3d,0.0d0,0,A,1) 224 225 do i=1,(nidb) 226 A(i) = dcpl_mb(tmp1(1) + int_mb(pack+i-1) - 1) 227 end do 228 229 value = BA_pop_stack(tmp1(2)) 230 if (.not.value) call errquit('Cram_c_pack:error popping stack',0, 231 & MA_ERR) 232 233 call c_Balance_c_balance(nb,A) 234 235 call nwpw_timing_end(9) 236 237 return 238 end 239 240 241 242* *********************************** 243* * * 244* * Cram_c_pack_start * 245* * * 246* *********************************** 247 248 subroutine Cram_c_pack_start(nb,A,tmp1,request,reqcnt,msgtype) 249 implicit none 250 integer nb 251 complex*16 A(*) 252 complex*16 tmp1(*) 253 integer request(*),reqcnt 254 integer msgtype 255 256#include "bafdecls.fh" 257#include "errquit.fh" 258#include "cram_common.fh" 259 260* **** local variables **** 261 integer i,nfft3d 262 integer nidb,pack 263 264* **** external functions **** 265 logical control_balance 266 external control_balance 267 268 call nwpw_timing_start(9) 269 270 nidb = int_mb(nidb_list(1)+nb) 271 pack = int_mb(pack_list(1)+2*nb) 272 273 call C3dB_nfft3d(1,nfft3d) 274 275 call dcopy(2*nfft3d,A,1,tmp1,1) 276 call dcopy(2*nfft3d,0.0d0,0,A,1) 277 278 do i=1,(nidb) 279 A(i) = tmp1(int_mb(pack+i-1)) 280 end do 281 282 if (control_balance()) 283 > call c_Balance_c_balance_start(nb,A,request,reqcnt,msgtype) 284 285 call nwpw_timing_end(9) 286 return 287 end 288 289 290* *********************************** 291* * * 292* * Cram_c_pack_end * 293* * * 294* *********************************** 295 296 subroutine Cram_c_pack_end(nb,tmp1,request,reqcnt) 297 implicit none 298 integer nb 299 complex*16 tmp1(*) 300 integer request(*),reqcnt 301 302#include "bafdecls.fh" 303#include "errquit.fh" 304#include "cram_common.fh" 305 306* **** external functions **** 307 logical control_balance 308 external control_balance 309 310 call nwpw_timing_start(9) 311 312 if (control_balance()) 313 > call c_Balance_c_balance_end(nb,tmp1,request,reqcnt) 314 315 call nwpw_timing_end(9) 316 return 317 end 318 319 320 321* *********************************** 322* * * 323* * Cram_r_pack * 324* * * 325* *********************************** 326 327 subroutine Cram_r_pack(nb,A) 328 implicit none 329 integer nb 330 real*8 A(*) 331 332#include "bafdecls.fh" 333#include "cram_common.fh" 334#include "errquit.fh" 335 336 337* **** local variables **** 338 logical value 339 integer i,nfft3d 340 integer tmp1(2) 341 integer nidb,pack 342 343 call nwpw_timing_start(9) 344 nidb = int_mb(nidb_list(1)+nb) 345 pack = int_mb(pack_list(1)+2*nb) 346 347 call C3dB_nfft3d(1,nfft3d) 348 value = BA_push_get(mt_dbl,nfft3d,'tmp1',tmp1(2),tmp1(1)) 349 if (.not.value) 350 > call errquit('Cram_r_pack:out of stack memory',0, MA_ERR) 351 352 call dcopy(nfft3d,A,1,dbl_mb(tmp1(1)),1) 353 call dcopy(nfft3d,0.0d0,0,A,1) 354 do i=1,(nidb) 355 A(i) = dbl_mb(tmp1(1) + int_mb(pack+i-1) - 1) 356 end do 357 358 value = BA_pop_stack(tmp1(2)) 359 if (.not.value) 360 > call errquit('Cram_r_pack:error popping stack',0, MA_ERR) 361 362 call c_Balance_t_balance(nb,A) 363 364 call nwpw_timing_end(9) 365 return 366 end 367 368 369* *********************************** 370* * * 371* * Cram_i_pack * 372* * * 373* *********************************** 374 375 subroutine Cram_i_pack(nb,A) 376 implicit none 377#include "errquit.fh" 378 integer nb 379 integer A(*) 380 381#include "bafdecls.fh" 382#include "cram_common.fh" 383 384 385* **** local variables **** 386 logical value 387 integer i,nfft3d 388 integer tmp1(2) 389 integer nidb,pack 390 391 call nwpw_timing_start(9) 392 nidb = int_mb(nidb_list(1)+nb) 393 pack = int_mb(pack_list(1)+2*nb) 394 395 call C3dB_nfft3d(1,nfft3d) 396 value = BA_push_get(mt_int,nfft3d,'tmp1',tmp1(2),tmp1(1)) 397 if (.not.value) 398 > call errquit('Cram_i_pack:out of stack memory',0, MA_ERR) 399 400 call icopy(nfft3d,A,1,int_mb(tmp1(1)),1) 401 call icopy(nfft3d,0,0,A,1) 402 do i=1,(nidb) 403 A(i) = int_mb(tmp1(1) + int_mb(pack+i-1) - 1) 404 end do 405 406 value = BA_pop_stack(tmp1(2)) 407 if (.not.value) 408 > call errquit('Cram_i_pack:error popping stack',1, MA_ERR) 409 410 call c_Balance_i_balance(nb,A) 411 412 call nwpw_timing_end(9) 413 return 414 end 415 416 417* *********************************** 418* * * 419* * Cram_c_unpack * 420* * * 421* *********************************** 422 423 subroutine Cram_c_unpack(nb,A) 424 implicit none 425#include "errquit.fh" 426 integer nb 427 complex*16 A(*) 428 429#include "bafdecls.fh" 430#include "cram_common.fh" 431 432 433* **** local variables **** 434 logical value 435 integer i,nfft3d 436 integer tmp1(2) 437 integer nidb,pack 438 439 call nwpw_timing_start(9) 440 nidb = int_mb(nidb_list(1)+nb) 441 pack = int_mb(pack_list(1)+2*nb) 442 443 call c_Balance_c_unbalance(nb,A) 444 445 call C3dB_nfft3d(1,nfft3d) 446 447 value = BA_push_get(mt_dcpl,(nidb), 448 > 'tmp1',tmp1(2),tmp1(1)) 449 if (.not.value) 450 > call errquit('Cram_c_unpack:out of stack memory',0, MA_ERR) 451 452 call dcopy(2*(nidb),A,1,dcpl_mb(tmp1(1)),1) 453 call dcopy(2*nfft3d,0.0d0,0,A,1) 454 do i=1,(nidb) 455 A(int_mb(pack+i-1)) = dcpl_mb(tmp1(1)+i-1) 456 end do 457 value = BA_pop_stack(tmp1(2)) 458 if (.not.value) 459 > call errquit('Cram_c_unpack:error popping stack',0, MA_ERR) 460 461 call nwpw_timing_end(9) 462 return 463 end 464 465 466* *********************************** 467* * * 468* * Cram_c_unpack_start * 469* * * 470* *********************************** 471 472 subroutine Cram_c_unpack_start(nb,A,tmp1,request,reqcnt,msgtype) 473 implicit none 474 integer nb 475 complex*16 A(*) 476 complex*16 tmp1(*) 477 integer request(*),reqcnt 478 integer msgtype 479 480#include "bafdecls.fh" 481#include "errquit.fh" 482#include "cram_common.fh" 483 484* **** external functions **** 485 logical control_balance 486 external control_balance 487 488 call nwpw_timing_start(9) 489 call Cram_c_Copy(nb,A,tmp1) 490 if (control_balance()) 491 > call c_Balance_c_unbalance_start(nb,tmp1,request,reqcnt,msgtype) 492 493 call nwpw_timing_end(9) 494 495 return 496 end 497 498 499* *********************************** 500* * * 501* * Cram_c_unpack_end * 502* * * 503* *********************************** 504 505 subroutine Cram_c_unpack_end(nb,tmp1,tmp2,request,reqcnt) 506 implicit none 507 integer nb 508 complex*16 tmp1(*) 509 complex*16 tmp2(*) 510 integer request(*),reqcnt 511 512#include "bafdecls.fh" 513#include "errquit.fh" 514#include "cram_common.fh" 515 516* **** local variables **** 517 integer i,nfft3d 518 integer nidb,pack 519 520* **** external functions **** 521 logical control_balance 522 external control_balance 523 524 call nwpw_timing_start(9) 525 nidb = int_mb(nidb_list(1)+nb) 526 pack = int_mb(pack_list(1)+2*nb) 527 528 call C3dB_nfft3d(1,nfft3d) 529 530 call nwpw_timing_start(9) 531 if (control_balance()) 532 > call c_Balance_c_unbalance_end(nb,tmp1,request,reqcnt) 533 534 call dcopy(2*(nidb),tmp1,1,tmp2,1) 535 call dcopy(2*nfft3d,0.0d0,0,tmp1,1) 536 do i=1,(nidb) 537 tmp1(int_mb(pack+i-1)) = tmp2(i) 538 end do 539 540 call nwpw_timing_end(9) 541 542 return 543 end 544 545 546 547* *********************************** 548* * * 549* * Cram_npack * 550* * * 551* *********************************** 552 553 subroutine Cram_npack(nb,npack) 554 implicit none 555 integer nb 556 integer npack 557 558#include "bafdecls.fh" 559#include "cram_common.fh" 560 561 562 npack = int_mb(nidb2_list(1)+nb) 563 return 564 end 565 566 567* *********************************** 568* * * 569* * Cram_max_npack * 570* * * 571* *********************************** 572 573 subroutine Cram_max_npack(npack) 574 implicit none 575 integer npack 576 577#include "bafdecls.fh" 578#include "cram_common.fh" 579 580 npack = npack_max 581 return 582 end 583 584* *********************************** 585* * * 586* * Cram_nwave * 587* * * 588* *********************************** 589 590 integer function Cram_nwave(nb) 591 implicit none 592 integer nb 593 594#include "bafdecls.fh" 595#include "cram_common.fh" 596 597 Cram_nwave = int_mb(nidb2_list(1)+nb) 598 return 599 end 600 601* *********************************** 602* * * 603* * Cram_nwave_brdcst * 604* * * 605* *********************************** 606 integer function Cram_nwave_brdcst(nb) 607 implicit none 608 integer nb 609 610#include "bafdecls.fh" 611#include "cram_common.fh" 612 613 integer nbq,taskid_k,p 614 integer iw 615 616 if (nb.eq.0) then 617 iw = int_mb(nidb2_list(1)) 618 else 619 iw = 0 620 call Parallel3d_taskid_k(taskid_k) 621 call K1dB_ktoqp(nb,nbq,p) 622 if (p.eq.taskid_k) iw = int_mb(nidb2_list(1)+nbq) 623 call K1dB_ISumAll(iw) 624 end if 625 626 Cram_nwave_brdcst = iw 627 return 628 end 629 630 631* *********************************** 632* * * 633* * Cram_nwave_all * 634* * * 635* *********************************** 636 integer function Cram_nwave_all(nb) 637 implicit none 638 integer nb 639 640#include "bafdecls.fh" 641#include "cram_common.fh" 642 643 Cram_nwave_all = int_mb(nwaveall_list(1)+nb) 644 return 645 end 646 647* *********************************** 648* * * 649* * Cram_nwave_all_brdcst * 650* * * 651* *********************************** 652 integer function Cram_nwave_all_brdcst(nb) 653 implicit none 654 integer nb 655 656#include "bafdecls.fh" 657#include "cram_common.fh" 658 659 integer nbq,taskid_k,p 660 integer iw 661 662 if (nb.eq.0) then 663 iw = int_mb(nwaveall_list(1)) 664 else 665 iw = 0 666 call Parallel3d_taskid_k(taskid_k) 667 call K1dB_ktoqp(nb,nbq,p) 668 if (p.eq.taskid_k) iw = int_mb(nwaveall_list(1)+nbq) 669 call K1dB_ISumAll(iw) 670 end if 671 672 Cram_nwave_all_brdcst = iw 673 return 674 end 675 676 677* *********************************** 678* * * 679* * Cram_zero * 680* * * 681* *********************************** 682 683 subroutine Cram_zero(nb,zero,pzero) 684 implicit none 685 integer nb 686 integer zero,pzero 687 688 integer qzero 689 690* ********************************************************* 691* **** warning this routine assumes a specific packing **** 692* ********************************************************* 693* index = (qzero-1)*(nx/2+1)*ny + (j-1)*(nx/2+1) + i 694 zero = 1 695 !call C3dB_ktoqp(1,1,qzero,pzero) 696 call C3dB_ijktoindexp(1,1,1,1,zero,pzero) 697 698 return 699 end 700 701 702 703* *********************************** 704* * * 705* * Cram_cc_nzdot * 706* * * 707* *********************************** 708 709 subroutine Cram_cc_nzdot(nb,ne,A,B,sum) 710 implicit none 711 integer nb 712 integer ne 713 complex*16 A(*) 714 complex*16 B(*) 715 complex*16 sum(*) 716 717#include "bafdecls.fh" 718#include "cram_common.fh" 719 720* **** local variables **** 721 integer n,np 722 integer shift 723 724* **** external functions **** 725 complex*16 tzdotc 726 external tzdotc 727 728 729 call nwpw_timing_start(2) 730 call Parallel3d_np_i(np) 731 732 do n=1,ne 733 shift = (n-1)*npack_max 734 735 sum(n) = tzdotc(int_mb(nidb2_list(1)+nb), 736 > A(1+shift),1, 737 > B(1), 1) 738 end do 739 740 if (np.gt.1) then 741 call C3dB_Vector_SumAll(2*ne,sum) 742 end if 743 744 call nwpw_timing_end(2) 745 746 return 747 end 748 749* *********************************** 750* * * 751* * Cram_cc_nzdot * 752* * * 753* *********************************** 754 755 subroutine Cram_cc_nzdot2com(nb,ne,ne1,A,B,sum) 756 implicit none 757 integer nb 758 integer ne,ne1 759 complex*16 A(*) 760 complex*16 B(*) 761 complex*16 sum(*) 762 763#include "bafdecls.fh" 764#include "cram_common.fh" 765 766* **** local variables **** 767 integer n,np 768 integer shift,shift2 769 770* **** external functions **** 771 complex*16 tzdotc 772 external tzdotc 773 774 775 call nwpw_timing_start(2) 776 call Parallel3d_np_i(np) 777 778 do n=1,ne 779 shift = (n-1)*npack_max 780 shift2 = ne1*npack_max 781 sum(n) = tzdotc(int_mb(nidb2_list(1)+nb), 782 > A(1+shift),1, 783 > B(1), 1) 784 sum(n) = sum(n) + tzdotc(int_mb(nidb2_list(1)+nb), 785 > A(1+shift+shift2),1, 786 > B(1+shift2),1) 787 end do 788 789 if (np.gt.1) then 790 call C3dB_Vector_SumAll(2*ne,sum) 791 end if 792 793 call nwpw_timing_end(2) 794 795 return 796 end 797 798 799* *********************************** 800* * * 801* * Cram_cc_inzdot * 802* * * 803* *********************************** 804 805 subroutine Cram_cc_inzdot(nb,ne,A,B,sum) 806 implicit none 807 integer nb 808 integer ne 809 complex*16 A(*) 810 complex*16 B(*) 811 complex*16 sum(*) 812 813#include "bafdecls.fh" 814#include "cram_common.fh" 815 816* **** local variables **** 817 integer n,np 818 integer shift 819 820* **** external functions **** 821 complex*16 tzdotc 822 external tzdotc 823 824 825 call nwpw_timing_start(2) 826 827 do n=1,ne 828 shift = (n-1)*npack_max 829 830 sum(n) = tzdotc(int_mb(nidb2_list(1)+nb), 831 > A(1+shift),1, 832 > B(1), 1) 833 end do 834 835 836 call nwpw_timing_end(2) 837 838 return 839 end 840 841* *********************************** 842* * * 843* * Cram_cc_inzdotAdd * 844* * * 845* *********************************** 846 847 subroutine Cram_cc_inzdotAdd(nb,ne,A,B,sum) 848 implicit none 849 integer nb 850 integer ne 851 complex*16 A(*) 852 complex*16 B(*) 853 complex*16 sum(*) 854 855#include "bafdecls.fh" 856#include "cram_common.fh" 857 858* **** local variables **** 859 integer n,np 860 integer shift 861 862* **** external functions **** 863 complex*16 tzdotc 864 external tzdotc 865 866 867 call nwpw_timing_start(2) 868 869 do n=1,ne 870 shift = (n-1)*npack_max 871 872 sum(n) = sum(n)+tzdotc(int_mb(nidb2_list(1)+nb), 873 > A(1+shift),1, 874 > B(1), 1) 875 end do 876 877 878 call nwpw_timing_end(2) 879 880 return 881 end 882 883 884* *********************************** 885* * * 886* * Cram_cc_izdot * 887* * * 888* *********************************** 889 890 subroutine Cram_cc_izdot(nb,A,B,sum) 891 implicit none 892 integer nb 893 complex*16 A(*) 894 complex*16 B(*) 895 complex*16 sum 896 897#include "bafdecls.fh" 898#include "cram_common.fh" 899 900 901* **** external functions **** 902 complex*16 tzdotc 903 external tzdotc 904 905 call nwpw_timing_start(2) 906 907 sum = tzdotc(int_mb(nidb2_list(1)+nb),A,1,B,1) 908 909 call nwpw_timing_end(2) 910 911 return 912 end 913 914 915* *********************************** 916* * * 917* * Cram_cc_izdotAdd * 918* * * 919* *********************************** 920 921 subroutine Cram_cc_izdotAdd(nb,A,B,sum) 922 implicit none 923 integer nb 924 complex*16 A(*) 925 complex*16 B(*) 926 complex*16 sum 927 928#include "bafdecls.fh" 929#include "cram_common.fh" 930 931* **** local variables **** 932 integer n,np 933 integer shift 934 935* **** external functions **** 936 complex*16 tzdotc 937 external tzdotc 938 939 940 call nwpw_timing_start(2) 941 942 sum = sum+tzdotc(int_mb(nidb2_list(1)+nb),A,1,B,1) 943 944 call nwpw_timing_end(2) 945 946 return 947 end 948 949 950 951* ******************************* 952* * * 953* * Cram_cc_zdot * 954* * * 955* ******************************* 956 957 subroutine Cram_cc_zdot(nb,A,B,tsum) 958 implicit none 959 integer nb 960 complex*16 A(*) 961 complex*16 B(*) 962 complex*16 tsum 963 964#include "bafdecls.fh" 965#include "cram_common.fh" 966 967* **** local variables **** 968 integer np 969 970* **** external functions **** 971 complex*16 tzdotc 972 external tzdotc 973 974 integer i 975 call nwpw_timing_start(2) 976 977 call Parallel3d_np_i(np) 978 979 tsum = tzdotc(int_mb(nidb2_list(1)+nb),A,1,B,1) 980 981 if (np.gt.1) call C3dB_Vector_SumAll(2,tsum) 982 983 call nwpw_timing_end(2) 984 985 return 986 end 987ccccccccccccccccccccccccccccccccccccccccccccccccccccc 988 complex*16 function tzdotc(n,x,incx,y,incy) 989 implicit none 990 integer n,incx,incy 991 complex*16 x(*),y(*) 992 integer i,ix,iy 993 complex*16 tsum 994 tsum = dcmplx(0.0d0,0.0d0) 995 if ((incx.eq.1).and.(incy.eq.1)) then 996 do i=1,n 997 tsum = tsum + dconjg(x(i))*y(i) 998 end do 999 else 1000 ix=1 1001 iy=1 1002 if (incx.lt.0) ix=1+(-n+1)*incx 1003 if (incy.lt.0) iy=1+(-n+1)*incy 1004 do i=1,n 1005 tsum = tsum + dconjg(x(ix))*y(iy) 1006 ix=ix+incx 1007 iy=iy+incy 1008 end do 1009 end if 1010 tzdotc = tsum 1011 return 1012 end 1013 1014 1015* *********************************** 1016* * * 1017* * Cram_cc_idot * 1018* * * 1019* *********************************** 1020 1021 subroutine Cram_cc_idot(nb,A,B,sum) 1022 implicit none 1023 integer nb 1024 complex*16 A(*) 1025 complex*16 B(*) 1026 real*8 sum 1027 1028#include "bafdecls.fh" 1029#include "cram_common.fh" 1030 1031 1032* **** local variables **** 1033c integer np 1034 1035* **** external functions **** 1036 real*8 ddot 1037 external ddot 1038 1039 1040 call nwpw_timing_start(2) 1041 1042 1043c call Parallel3d_np_i(np) 1044 1045 sum = ddot(2*int_mb(nidb2_list(1)+nb),A,1,B,1) 1046 1047c if (np.gt.1) call C3dB_SumAll(sum) 1048 1049 call nwpw_timing_end(2) 1050 1051 return 1052 end 1053 1054 1055 1056 1057* *********************************** 1058* * * 1059* * Cram_cc_dot * 1060* * * 1061* *********************************** 1062 1063 subroutine Cram_cc_dot(nb,A,B,sum) 1064 implicit none 1065 integer nb 1066 complex*16 A(*) 1067 complex*16 B(*) 1068 real*8 sum 1069 1070#include "bafdecls.fh" 1071#include "cram_common.fh" 1072 1073* **** local variables **** 1074 integer np 1075 1076* **** external functions **** 1077 real*8 ddot 1078 external ddot 1079 1080 1081 call nwpw_timing_start(2) 1082 call Parallel3d_np_i(np) 1083 1084 sum = ddot(2*int_mb(nidb2_list(1)+nb),A,1,B,1) 1085 1086 if (np.gt.1) call C3dB_SumAll(sum) 1087 1088 call nwpw_timing_end(2) 1089 return 1090 end 1091 1092 1093 1094 1095* *********************************** 1096* * * 1097* * Cram_ccm_idgemm * 1098* * * 1099* *********************************** 1100 1101 subroutine Cram_ccm_idgemm(nb,n,m,A,B,alpha,beta,hml) 1102 implicit none 1103 integer nb,n,m 1104 complex*16 A(*) 1105 complex*16 B(*) 1106 real*8 alpha,beta 1107 real*8 hml(*) 1108 1109#include "bafdecls.fh" 1110#include "cram_common.fh" 1111 1112* **** local variables **** 1113 integer np,npack1 1114 1115 1116 call nwpw_timing_start(2) 1117 call Parallel3d_np_i(np) 1118 npack1 = int_mb(nidb2_list(1)+nb) 1119 1120 call DGEMM('T','N',n,m,2*npack1, 1121 > alpha, 1122 > A, 2*npack_max, 1123 > B, 2*npack_max, 1124 > beta, 1125 > hml,n) 1126 1127 1128 call nwpw_timing_end(2) 1129 return 1130 end 1131 1132 1133 1134* *********************************** 1135* * * 1136* * Cram_ccm_izgemm * 1137* * * 1138* *********************************** 1139 1140 subroutine Cram_ccm_izgemm(nb,n,m,A,B,alpha,beta,hml) 1141 implicit none 1142 integer nb,n,m 1143 complex*16 A(*) 1144 complex*16 B(*) 1145 complex*16 alpha,beta 1146 complex*16 hml(*) 1147 1148#include "bafdecls.fh" 1149#include "cram_common.fh" 1150 1151* **** local variables **** 1152 integer npack1 1153 1154 1155 call nwpw_timing_start(2) 1156 npack1 = int_mb(nidb2_list(1)+nb) 1157 1158 call ZGEMM('C','N',n,m,npack1, 1159 > alpha, 1160 > A, npack_max, 1161 > B, npack_max, 1162 > beta, 1163 > hml,n) 1164 1165 1166 call nwpw_timing_end(2) 1167 return 1168 end 1169 1170 1171* *********************************** 1172* * * 1173* * Cram_ccm_sym_izgemm * 1174* * * 1175* *********************************** 1176 1177 subroutine Cram_ccm_sym_izgemm(nb,n,A,B,alpha,beta,hml) 1178 implicit none 1179 integer nb,n 1180 complex*16 A(*) 1181 complex*16 B(*) 1182 complex*16 alpha,beta 1183 complex*16 hml(n,n) 1184 1185#include "bafdecls.fh" 1186#include "cram_common.fh" 1187 1188* **** local variables **** 1189 integer j,k 1190 integer npack1 1191 1192 call nwpw_timing_start(2) 1193 npack1 = int_mb(nidb2_list(1)+nb) 1194 1195 do k=1,n 1196 call ZGEMM('C','N',k,1,npack1, 1197 > alpha, 1198 > A, npack_max, 1199 > B(1+(k-1)*npack_max), npack_max, 1200 > beta, 1201 > hml(1,k),k) 1202 end do 1203 1204 do k=1,n 1205 do j=k+1,n 1206 hml(j,k) = dconjg(hml(k,j)) 1207 end do 1208 end do 1209 do k=1,n 1210 hml(k,k) = dcmplx(dble(hml(k,k)),0.0d0) 1211 end do 1212 1213c do k=1,n 1214c do j=1,n 1215c write(*,*) "sym_izgemm=",j,k,hml(j,k) 1216c end do 1217c end do 1218 1219 1220 call nwpw_timing_end(2) 1221 return 1222 end 1223 1224 1225 1226* *********************************** 1227* * * 1228* * Cram_cmc_zgemm * 1229* * * 1230* *********************************** 1231 1232 subroutine Cram_cmc_zgemm(nb,n,m,A,B,alpha,beta,hml) 1233 implicit none 1234 integer nb,n,m 1235 complex*16 A(*) 1236 complex*16 B(*) 1237 complex*16 alpha,beta 1238 complex*16 hml(*) 1239 1240#include "bafdecls.fh" 1241#include "cram_common.fh" 1242 1243* **** local variables **** 1244 integer np,npack1 1245 1246 1247 call nwpw_timing_start(2) 1248 call Parallel3d_np_i(np) 1249 npack1 = int_mb(nidb2_list(1)+nb) 1250 1251 call ZGEMM('N','N',npack1,n,m, 1252 > alpha, 1253 > A, npack_max, 1254 > hml,n, 1255 > beta, 1256 > B, npack_max) 1257 1258 1259 call nwpw_timing_end(2) 1260 return 1261 end 1262 1263 1264 1265 1266* *********************************** 1267* * * 1268* * Cram_rr_idot * 1269* * * 1270* *********************************** 1271 1272 subroutine Cram_rr_idot(nb,A,B,sum) 1273 implicit none 1274 integer nb 1275 real*8 A(*) 1276 real*8 B(*) 1277 real*8 sum 1278 1279#include "bafdecls.fh" 1280#include "cram_common.fh" 1281 1282* **** external functions **** 1283 real*8 ddot 1284 external ddot 1285 1286 call nwpw_timing_start(2) 1287 1288 sum = ddot(int_mb(nidb2_list(1)+nb),A(1),1,B(1),1) 1289 1290 call nwpw_timing_end(2) 1291 return 1292 end 1293 1294 1295 1296 1297 1298* *********************************** 1299* * * 1300* * Cram_rr_dot * 1301* * * 1302* *********************************** 1303 1304 subroutine Cram_rr_dot(nb,A,B,sum) 1305 implicit none 1306 integer nb 1307 real*8 A(*) 1308 real*8 B(*) 1309 real*8 sum 1310 1311#include "bafdecls.fh" 1312#include "cram_common.fh" 1313 1314* **** local variables **** 1315 integer np 1316 1317* **** external functions **** 1318 real*8 ddot 1319 external ddot 1320 1321 call nwpw_timing_start(2) 1322 1323 call Parallel3d_np_i(np) 1324 1325 sum = ddot(int_mb(nidb2_list(1)+nb),A(1),1,B(1),1) 1326 1327 if (np.gt.1) call C3dB_SumAll(sum) 1328 1329 call nwpw_timing_end(2) 1330 return 1331 end 1332 1333 1334* *********************************** 1335* * * 1336* * Cram_r_dsum * 1337* * * 1338* *********************************** 1339 1340 subroutine Cram_r_dsum(nb,A,sum) 1341 implicit none 1342 integer nb 1343 real*8 A(*) 1344 real*8 sum 1345 1346#include "bafdecls.fh" 1347#include "cram_common.fh" 1348 1349* **** local variables **** 1350 integer np 1351 1352* **** external functions **** 1353 real*8 dsum 1354 external dsum 1355 1356 call nwpw_timing_start(2) 1357 call Parallel3d_np_i(np) 1358 1359 sum = dsum(int_mb(nidb2_list(1)+nb),A(1),1) 1360 1361 if (np.gt.1) call C3dB_SumAll(sum) 1362 1363 call nwpw_timing_end(2) 1364 return 1365 end 1366 1367 1368 1369* *********************************** 1370* * * 1371* * Cram_c_Copy * 1372* * * 1373* *********************************** 1374 1375 subroutine Cram_c_Copy(nb,A,B) 1376 implicit none 1377 integer nb 1378 complex*16 A(*) 1379 complex*16 B(*) 1380 1381#include "bafdecls.fh" 1382#include "cram_common.fh" 1383 1384 call dcopy(2*(int_mb(nidb2_list(1)+nb)),A,1,B,1) 1385 1386 return 1387 end 1388 1389* *********************************** 1390* * * 1391* * Cram_r_Copy * 1392* * * 1393* *********************************** 1394 1395 subroutine Cram_r_Copy(nb,A,B) 1396 implicit none 1397 integer nb 1398 real*8 A(*) 1399 real*8 B(*) 1400 1401#include "bafdecls.fh" 1402#include "cram_common.fh" 1403 1404 call dcopy((int_mb(nidb2_list(1)+nb)),A,1,B,1) 1405 1406 return 1407 end 1408 1409* *********************************** 1410* * * 1411* * Cram_rc_Copy * 1412* * * 1413* *********************************** 1414 1415 subroutine Cram_rc_Copy(nb,A,B) 1416 implicit none 1417 integer nb 1418 real*8 A(*) 1419 complex*16 B(*) 1420 1421#include "bafdecls.fh" 1422#include "cram_common.fh" 1423 1424 integer i 1425 1426 do i=1,(int_mb(nidb2_list(1)+nb)) 1427 B(i) = dcmplx(A(i),0.0d0) 1428 end do 1429 return 1430 end 1431 1432* *********************************** 1433* * * 1434* * Cram_cr_Copy * 1435* * * 1436* *********************************** 1437 1438 subroutine Cram_cr_Copy(nb,A,B) 1439 implicit none 1440 integer nb 1441 complex*16 A(*) 1442 real*8 B(*) 1443 1444#include "bafdecls.fh" 1445#include "cram_common.fh" 1446 1447 integer i 1448 1449 do i=1,(int_mb(nidb2_list(1)+nb)) 1450 B(i) = dble(A(i)) 1451 end do 1452 return 1453 end 1454 1455 1456* *********************************** 1457* * * 1458* * Cram_c_Zero * 1459* * * 1460* *********************************** 1461 1462 subroutine Cram_c_Zero(nb,A) 1463 implicit none 1464 integer nb 1465 complex*16 A(*) 1466 1467#include "bafdecls.fh" 1468#include "cram_common.fh" 1469 1470 call dcopy(2*(int_mb(nidb2_list(1)+nb)),0.0d0,0,A,1) 1471 1472 return 1473 end 1474 1475 1476 1477* *********************************** 1478* * * 1479* * Cram_cc_Sum * 1480* * * 1481* *********************************** 1482 1483 subroutine Cram_cc_Sum(nb,A,B,C) 1484 implicit none 1485 integer nb 1486 complex*16 A(*) 1487 complex*16 B(*) 1488 complex*16 C(*) 1489 1490#include "bafdecls.fh" 1491#include "cram_common.fh" 1492 1493* **** local variables **** 1494 integer i 1495 1496 do i=1,(int_mb(nidb2_list(1)+nb)) 1497 C(i) = A(i) + B(i) 1498 end do 1499 1500 return 1501 end 1502 1503 1504* *********************************** 1505* * * 1506* * Cram_cc_Sum2 * 1507* * * 1508* *********************************** 1509 1510 subroutine Cram_cc_Sum2(nb,A,B) 1511 implicit none 1512 integer nb 1513 complex*16 A(*) 1514 complex*16 B(*) 1515 1516#include "bafdecls.fh" 1517#include "cram_common.fh" 1518 1519* **** local variables **** 1520 integer i 1521 1522 do i=1,(int_mb(nidb2_list(1)+nb)) 1523 B(i) = B(i) + A(i) 1524 end do 1525 1526 return 1527 end 1528 1529 1530* *********************************** 1531* * * 1532* * Cram_rr_Sum * 1533* * * 1534* *********************************** 1535 1536 subroutine Cram_rr_Sum(nb,A,B,C) 1537 implicit none 1538 integer nb 1539 real*8 A(*) 1540 real*8 B(*) 1541 real*8 C(*) 1542 1543#include "bafdecls.fh" 1544#include "cram_common.fh" 1545 1546* **** local variables **** 1547 integer i 1548 1549 do i=1,(int_mb(nidb2_list(1)+nb)) 1550 C(i) = A(i) + B(i) 1551 end do 1552 1553 return 1554 end 1555 1556 1557 1558* *********************************** 1559* * * 1560* * Cram_rr_Sum2 * 1561* * * 1562* *********************************** 1563 1564 subroutine Cram_rr_Sum2(nb,A,B) 1565 implicit none 1566 integer nb 1567 real*8 A(*) 1568 real*8 B(*) 1569 1570#include "bafdecls.fh" 1571#include "cram_common.fh" 1572 1573* **** local variables **** 1574 integer i 1575 1576 do i=1,(int_mb(nidb2_list(1)+nb)) 1577 B(i) = B(i) + A(i) 1578 end do 1579 1580 return 1581 end 1582 1583 1584* *********************************** 1585* * * 1586* * Cram_cc_Sub * 1587* * * 1588* *********************************** 1589 1590 subroutine Cram_cc_Sub(nb,A,B,C) 1591 implicit none 1592 integer nb 1593 complex*16 A(*) 1594 complex*16 B(*) 1595 complex*16 C(*) 1596 1597#include "bafdecls.fh" 1598#include "cram_common.fh" 1599 1600* **** local variables **** 1601 integer i 1602 1603 do i=1,(int_mb(nidb2_list(1)+nb)) 1604 C(i) = A(i) - B(i) 1605 end do 1606 1607 return 1608 end 1609 1610* *********************************** 1611* * * 1612* * Cram_rr_Sub * 1613* * * 1614* *********************************** 1615 1616 subroutine Cram_rr_Sub(nb,A,B,C) 1617 implicit none 1618 integer nb 1619 real*8 A(*) 1620 real*8 B(*) 1621 real*8 C(*) 1622 1623#include "bafdecls.fh" 1624#include "cram_common.fh" 1625 1626* **** local variables **** 1627 integer i 1628 1629 do i=1,(int_mb(nidb2_list(1)+nb)) 1630 C(i) = A(i) - B(i) 1631 end do 1632 1633 return 1634 end 1635 1636* *********************************** 1637* * * 1638* * Cram_rr_Sqrt * 1639* * * 1640* *********************************** 1641 1642 subroutine Cram_rr_Sqrt(nb,A,C) 1643 implicit none 1644 integer nb 1645 real*8 A(*) 1646 real*8 C(*) 1647 1648#include "bafdecls.fh" 1649#include "cram_common.fh" 1650 1651* **** local variables **** 1652 integer i 1653 1654 do i=1,(int_mb(nidb2_list(1)+nb)) 1655 C(i) = dsqrt(A(i)) 1656 end do 1657 1658 return 1659 end 1660 1661 1662* *********************************** 1663* * * 1664* * Cram_rr_Sqrt1 * 1665* * * 1666* *********************************** 1667 1668 subroutine Cram_rr_Sqrt1(nb,A) 1669 implicit none 1670 integer nb 1671 real*8 A(*) 1672 1673#include "bafdecls.fh" 1674#include "cram_common.fh" 1675 1676* **** local variables **** 1677 integer i 1678 1679 do i=1,(int_mb(nidb2_list(1)+nb)) 1680 A(i) = dsqrt(A(i)) 1681 end do 1682 1683 return 1684 end 1685 1686 1687* *********************************** 1688* * * 1689* * Cram_cr_Sqr * 1690* * * 1691* *********************************** 1692 1693 subroutine Cram_cr_Sqr(nb,A,C) 1694 implicit none 1695 integer nb 1696 complex*16 A(*) 1697 real*8 C(*) 1698 1699#include "bafdecls.fh" 1700#include "cram_common.fh" 1701 1702 1703* **** local variables **** 1704 integer i 1705 1706 do i=1,(int_mb(nidb2_list(1)+nb)) 1707 C(i) = dble(A(i))**2 + dimag(A(i))**2 1708 end do 1709 1710 return 1711 end 1712 1713 1714 1715* *********************************** 1716* * * 1717* * Cram_c_SMul * 1718* * * 1719* *********************************** 1720 1721 subroutine Cram_c_SMul(nb,alpha,A,C) 1722 implicit none 1723 integer nb 1724 real*8 alpha 1725 complex*16 A(*) 1726 complex*16 C(*) 1727 1728#include "bafdecls.fh" 1729#include "cram_common.fh" 1730 1731* **** local variables **** 1732 integer i 1733 1734 do i=1,(int_mb(nidb2_list(1)+nb)) 1735 C(i) = alpha*A(i) 1736 end do 1737 1738 return 1739 end 1740 1741 1742 1743* *********************************** 1744* * * 1745* * Cram_c_SMul1 * 1746* * * 1747* *********************************** 1748 1749 subroutine Cram_c_SMul1(nb,alpha,A) 1750 implicit none 1751 integer nb 1752 real*8 alpha 1753 complex*16 A(*) 1754 1755#include "bafdecls.fh" 1756#include "cram_common.fh" 1757 1758* **** local variables **** 1759 integer i 1760 1761 do i=1,(int_mb(nidb2_list(1)+nb)) 1762 A(i) = alpha*A(i) 1763 end do 1764 1765 return 1766 end 1767 1768 1769 1770* *********************************** 1771* * * 1772* * Cram_c_ZMul * 1773* * * 1774* *********************************** 1775 1776 subroutine Cram_c_ZMul(nb,alpha,A,C) 1777 implicit none 1778 integer nb 1779 complex*16 alpha 1780 complex*16 A(*) 1781 complex*16 C(*) 1782 1783#include "bafdecls.fh" 1784#include "cram_common.fh" 1785 1786* **** local variables **** 1787 integer i 1788 1789 do i=1,(int_mb(nidb2_list(1)+nb)) 1790 C(i) = alpha*A(i) 1791 end do 1792 1793 return 1794 end 1795 1796 1797 1798* *********************************** 1799* * * 1800* * Cram_r_SMul * 1801* * * 1802* *********************************** 1803 1804 subroutine Cram_r_SMul(nb,alpha,A,C) 1805 implicit none 1806 integer nb 1807 real*8 alpha 1808 real*8 A(*) 1809 real*8 C(*) 1810 1811#include "bafdecls.fh" 1812#include "cram_common.fh" 1813 1814* **** local variables **** 1815 integer i,nn 1816 1817 nn = int_mb(nidb2_list(1)+nb) 1818 1819 do i=1,nn 1820 C(i) = alpha*A(i) 1821 end do 1822 1823 return 1824 end 1825 1826 1827* *********************************** 1828* * * 1829* * Cram_r_SMul1 * 1830* * * 1831* *********************************** 1832 1833 subroutine Cram_r_SMul1(nb,alpha,A) 1834 implicit none 1835 integer nb 1836 real*8 alpha 1837 real*8 A(*) 1838 1839#include "bafdecls.fh" 1840#include "cram_common.fh" 1841 1842* **** local variables **** 1843 integer i 1844 1845 do i=1,(int_mb(nidb2_list(1)+nb)) 1846 A(i) = alpha*A(i) 1847 end do 1848 1849 return 1850 end 1851 1852 1853 1854* *********************************** 1855* * * 1856* * Cram_cc_daxpy * 1857* * * 1858* *********************************** 1859 1860 subroutine Cram_cc_daxpy(nb,alpha,A,B) 1861 implicit none 1862 integer nb 1863 real*8 alpha 1864 complex*16 A(*) 1865 complex*16 B(*) 1866 1867#include "bafdecls.fh" 1868#include "cram_common.fh" 1869 1870 call daxpy(2*(int_mb(nidb2_list(1)+nb)),alpha,A,1,B,1) 1871 1872 return 1873 end 1874 1875 1876 1877* *********************************** 1878* * * 1879* * Cram_rr_daxpy * 1880* * * 1881* *********************************** 1882 1883 subroutine Cram_rr_daxpy(nb,alpha,A,C) 1884 implicit none 1885 integer nb 1886 real*8 alpha 1887 real*8 A(*) 1888 real*8 C(*) 1889 1890#include "bafdecls.fh" 1891#include "cram_common.fh" 1892 1893* **** local variables **** 1894 integer i 1895 1896 do i=1,(int_mb(nidb2_list(1)+nb)) 1897 C(i) = C(i) + alpha*A(i) 1898 end do 1899 1900 return 1901 end 1902 1903 1904 1905* *********************************** 1906* * * 1907* * Cram_cc_zaxpy * 1908* * * 1909* *********************************** 1910 1911 subroutine Cram_cc_zaxpy(nb,alpha,A,C) 1912 implicit none 1913 integer nb 1914 complex*16 alpha 1915 complex*16 A(*) 1916 complex*16 C(*) 1917 1918#include "bafdecls.fh" 1919#include "cram_common.fh" 1920 1921* **** local variables **** 1922 integer i 1923 1924 do i=1,(int_mb(nidb2_list(1)+nb)) 1925 C(i) = C(i) + alpha*A(i) 1926 end do 1927 1928 return 1929 end 1930 1931 1932* *********************************** 1933* * * 1934* * Cram_rc_Mul * 1935* * * 1936* *********************************** 1937 1938 subroutine Cram_rc_Mul(nb,A,B,C) 1939 implicit none 1940 integer nb 1941 real*8 A(*) 1942 complex*16 B(*) 1943 complex*16 C(*) 1944 1945#include "bafdecls.fh" 1946#include "cram_common.fh" 1947 1948* **** local variables **** 1949 integer i 1950 1951 do i=1,( int_mb(nidb2_list(1)+nb)) 1952 C(i) = A(i)*B(i) 1953 end do 1954 1955 return 1956 end 1957 1958 1959* *********************************** 1960* * * 1961* * Cram_rc_Mul2 * 1962* * * 1963* *********************************** 1964 1965 subroutine Cram_rc_Mul2(nb,A,B) 1966 implicit none 1967 integer nb 1968 real*8 A(*) 1969 complex*16 B(*) 1970 1971#include "bafdecls.fh" 1972#include "cram_common.fh" 1973 1974* **** local variables **** 1975 integer i 1976 1977 do i=1,( int_mb(nidb2_list(1)+nb)) 1978 B(i) = B(i)*A(i) 1979 end do 1980 1981 return 1982 end 1983 1984* *********************************** 1985* * * 1986* * Cram_irc_Mul * 1987* * * 1988* *********************************** 1989 1990 subroutine Cram_irc_Mul(nb,A,B,C) 1991 implicit none 1992 integer nb 1993 real*8 A(*) 1994 complex*16 B(*) 1995 complex*16 C(*) 1996 1997#include "bafdecls.fh" 1998#include "cram_common.fh" 1999 2000* **** local variables **** 2001 integer i 2002 2003 do i=1,(int_mb(nidb2_list(1)+nb)) 2004 C(i) = A(i) * dcmplx(-dimag(B(i)),dble(B(i))) 2005 end do 2006 2007 return 2008 end 2009 2010 2011* *********************************** 2012* * * 2013* * Cram_icc_Mul * 2014* * * 2015* *********************************** 2016 2017 subroutine Cram_icc_Mul(nb,A,B,C) 2018 implicit none 2019 integer nb 2020 complex*16 A(*) 2021 complex*16 B(*) 2022 complex*16 C(*) 2023 2024#include "bafdecls.fh" 2025#include "cram_common.fh" 2026 2027* **** local variables **** 2028 integer i 2029 2030 do i=1,(int_mb(nidb2_list(1)+nb)) 2031 C(i) = A(i) * dcmplx(-dimag(B(i)),dble(B(i))) 2032 end do 2033 2034 return 2035 end 2036 2037 2038 2039* *********************************** 2040* * * 2041* * Cram_cc_Mul * 2042* * * 2043* *********************************** 2044 2045 subroutine Cram_cc_Mul(nb,A,B,C) 2046 implicit none 2047 integer nb 2048 complex*16 A(*) 2049 complex*16 B(*) 2050 complex*16 C(*) 2051 2052#include "bafdecls.fh" 2053#include "cram_common.fh" 2054 2055* **** local variables **** 2056 integer i 2057 2058 do i=1,(int_mb(nidb2_list(1)+nb)) 2059 C(i) = A(i)*B(i) 2060 end do 2061 2062 return 2063 end 2064 2065 2066* *********************************** 2067* * * 2068* * Cram_cc_Mul2 * 2069* * * 2070* *********************************** 2071 2072 subroutine Cram_cc_Mul2(nb,A,B) 2073 implicit none 2074 integer nb 2075 complex*16 A(*) 2076 complex*16 B(*) 2077 2078#include "bafdecls.fh" 2079#include "cram_common.fh" 2080 2081* **** local variables **** 2082 integer i 2083 2084 do i=1,(int_mb(nidb2_list(1)+nb)) 2085 B(i) = B(i)*A(i) 2086 end do 2087 2088 return 2089 end 2090 2091 2092* *********************************** 2093* * * 2094* * Cram_ccr_conjgMul * 2095* * * 2096* *********************************** 2097 subroutine Cram_ccr_conjgMul(nb,A,B,C) 2098 implicit none 2099 integer nb 2100 complex*16 A(*) 2101 complex*16 B(*) 2102 real*8 C(*) 2103 2104#include "bafdecls.fh" 2105#include "cram_common.fh" 2106 2107* **** local variables **** 2108 integer i 2109 2110 do i=1,( int_mb(nidb2_list(1)+nb)) 2111 C(i) = dble(dconjg(A(i))*B(i)) 2112 end do 2113 2114 return 2115 end 2116 2117 2118* *********************************** 2119* * * 2120* * Cram_cc_conjgMul * 2121* * * 2122* *********************************** 2123 2124 subroutine Cram_cc_conjgMul(nb,A,B,C) 2125 implicit none 2126 integer nb 2127 complex*16 A(*) 2128 complex*16 B(*) 2129 complex*16 C(*) 2130 2131#include "bafdecls.fh" 2132#include "cram_common.fh" 2133 2134* **** local variables **** 2135 integer i 2136 2137 do i=1,(int_mb(nidb2_list(1)+nb)) 2138 C(i) = dconjg(A(i))*B(i) 2139 end do 2140 2141 return 2142 end 2143 2144 2145 2146 2147* *********************************** 2148* * * 2149* * Cram_rc_iMul * 2150* * * 2151* *********************************** 2152 2153 subroutine Cram_rc_iMul(nb,A,B,C) 2154 implicit none 2155 integer nb 2156 real*8 A(*) 2157 complex*16 B(*) 2158 complex*16 C(*) 2159 2160#include "bafdecls.fh" 2161#include "cram_common.fh" 2162 2163* **** local variables **** 2164 integer i 2165 2166 do i=1,(int_mb(nidb2_list(1)+nb)) 2167 C(i) = dcmplx(0.0d0,A(i)) * B(i) 2168 end do 2169 2170 return 2171 end 2172 2173 2174 2175* *********************************** 2176* * * 2177* * Cram_rr_Mul * 2178* * * 2179* *********************************** 2180 2181 subroutine Cram_rr_Mul(nb,A,B,C) 2182 implicit none 2183 integer nb 2184 real*8 A(*) 2185 real*8 B(*) 2186 real*8 C(*) 2187 2188#include "bafdecls.fh" 2189#include "cram_common.fh" 2190 2191* **** local variables **** 2192 integer i 2193 2194 do i=1,(int_mb(nidb2_list(1)+nb)) 2195 C(i) = A(i)*B(i) 2196 end do 2197 2198 return 2199 end 2200 2201* *********************************** 2202* * * 2203* * Cram_rr_Mul2 * 2204* * * 2205* *********************************** 2206 2207 subroutine Cram_rr_Mul2(nb,A,B) 2208 implicit none 2209 integer nb 2210 real*8 A(*) 2211 real*8 B(*) 2212 2213#include "bafdecls.fh" 2214#include "cram_common.fh" 2215 2216* **** local variables **** 2217 integer i 2218 2219 do i=1,(int_mb(nidb2_list(1)+nb)) 2220 B(i) = B(i)*A(i) 2221 end do 2222 2223 return 2224 end 2225 2226 2227 2228* *********************************** 2229* * * 2230* * Cram_zccr_Mulitply2 * 2231* * * 2232* *********************************** 2233* 2234* This routine used by cpsp force routines computes 2235* C(i) = Imag (alpha*A(i)*dconjg(B(i))) 2236* 2237 subroutine Cram_zccr_Multiply2(nb,alpha,A,B,C) 2238 implicit none 2239 integer nb 2240 complex*16 alpha 2241 complex*16 A(*) 2242 complex*16 B(*) 2243 real*8 C(*) 2244 2245#include "bafdecls.fh" 2246#include "cram_common.fh" 2247 2248* **** local variables **** 2249 integer i 2250 2251 do i=1,(int_mb(nidb2_list(1)+nb)) 2252 C(i) = dimag(alpha*A(i)*dconjg(B(i))) 2253 end do 2254 2255 return 2256 end 2257 2258* *********************************** 2259* * * 2260* * Cram_zccr_Mulitply2Add * 2261* * * 2262* *********************************** 2263* 2264* This routine used by cpsp force routines computes 2265* C(i) = Imag (alpha*A(i)*dconjg(B(i))) 2266* 2267 subroutine Cram_zccr_Multiply2Add(nb,alpha,A,B,C) 2268 implicit none 2269 integer nb 2270 complex*16 alpha 2271 complex*16 A(*) 2272 complex*16 B(*) 2273 real*8 C(*) 2274 2275#include "bafdecls.fh" 2276#include "cram_common.fh" 2277 2278* **** local variables **** 2279 integer i 2280 2281 do i=1,(int_mb(nidb2_list(1)+nb)) 2282 C(i) = C(i)+dimag(alpha*A(i)*dconjg(B(i))) 2283 end do 2284 2285 return 2286 end 2287 2288 2289 2290