API tools faq
paste
Advertisement
Guest User

Untitled

a guest
Oct 2nd, 2018
116
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Fortran 9.25 KB | None | 0 0
raw download clone embed print report
  1.       SUBROUTINE DORBAX (F,NP1,NP2,BBXS,PWP,NBB,NPW,XLFT,XRGT,NAQRS,
  2.      $ NND,NNE,MAXLOT,igam)
  3. c
  4.       implicit double precision (a-h,o-z)
  5.       parameter (NHT=7,nh4=4*nht-3)
  6.       parameter (nhq=nht*nht*nht*nht)
  7.       parameter (zero=0.d0)
  8.       DIMENSION F(igam),BBXS(3*nbb),PWP(NPW,3),XLFT(NND,MAXLOT,3),
  9.      # XRGT(NNE,MAXLOT,3),BX(NAQRS*(NHKTA+NHKTB+NHKTC+NHKTD))
  10.       COMMON /XA/ XAND(3,8,nh4)
  11.       COMMON /CONI/ NUCT(4),NRCT(4),JSTT(4),
  12.      1 IFD1X,IFD2X,NUCAQ,NUCRS,NXXXX,NCQRS,KQ,KR,KS,IVA
  13.       LOGICAL IFPL,NP1(3),NP2(3)
  14.       COMMON /TAB/ NNN
  15.       COMMON / /FACT(nh4),RFACT(nh4),FACTM(nh4),RFACTM(nh4)
  16.      & ,MAA,IFD1,IFD2,KCD,KBCD,NHCD,NHBCD,NNC,IFPL(3),
  17.      & NN1,NN2,NHKTA,NHKTB,NHKTC,NHKTD,KHKTA,KHKTB,KHKTC,KHKTD
  18.      $ ,NNB
  19.       COMMON /DOR/ INDXY(nhq),INDZ(nhq)
  20. C
  21. C     THE DIMENSION OF INDXY AND INDZ IS (L+1)**4
  22. C
  23.       NUC12=NUCT(1)*NUCT(2)
  24.       NUMM=0
  25.       NUML=0
  26.       NUM=0
  27.       NUMJ=0
  28.       DO 111 NA=1,NHKTA
  29.       DO 1 NB=1,NHKTB
  30.       I12=NA+NB-1
  31.       NUMI=0
  32.       DO 341 NC=1,NHKTC
  33.       DO 340 ND=1,NHKTD
  34.       I34=NC+ND-1
  35.       NUML=NUML+1
  36.       JT=I12+I34-1
  37.       ijt = (jt-1)*naqrs+1
  38.       IM=NNN-JT+1
  39.       INDXY(NUML)=NUM*NCQRS+1
  40.       INDZ(NUML)=((NUML-1)*(NNN+1)-NUM-1)*NCQRS+1
  41.       IB4= NCQRS * (NUM+1)
  42.       IB6 = NUCAQ * (NUMJ+I12-1)
  43. c
  44.       IF (JT.GT.2) THEN
  45. C
  46. C.... GENERAL CODE.  USE WHEN SUM OF Q.N.  (JT) >2
  47. C
  48.         DO 230 IXY=1,3
  49.         IB44=IB4 + (IXY-1)*NBB
  50.         do i=1,NAQRS*(JT-1)
  51.           bx(naqrs+i) = zero
  52.         enddo
  53.         DO 52 IP=2,JT
  54.         IF(NP1(IXY) .AND. NP2(IXY) .AND. (iand(1,IP+JT).EQ.1))
  55.      1    GO TO 52
  56.         MIN9=MAX(1,IP-I12+1)
  57.         MAX9=MIN(IP,I34)
  58.         DO 51 I9=MIN9,MAX9
  59.         IF((NP2(IXY) .AND. (iand(1,I9+I34).EQ.1)) .OR.
  60.      1   (NP1(IXY) .AND. (iand(1,I9+I12-IP).EQ.0))) GO TO 51
  61. C...  NO SIMPLIFICATIONS DUE TO PLANARITY
  62.         CALL MPPROD (BX,XLFT(NUCAQ*(NUMJ+IP-I9)+1,KQ,IXY),NUCAQ,
  63.      $  XRGT(NUC12*(NUMI+I9-1)+1,IVA,IXY),NUCRS,
  64.      $  XAND(IXY,KR,I34-I9+1)*FACT(IP))
  65.         ii = (ip-1)*naqrs
  66.         DO 53 I=1,NAQRS
  67.         BX(I+ii) = BX(I+ii) + BX(I)
  68.  53     continue
  69.    51   CONTINUE
  70.    52   CONTINUE
  71.         CALL CMPRI(BX(naqrs+1),BX,JT-1)
  72. C
  73. C...  FIRST DO THE CONTRIBUTIONS FOR I = 1&2
  74. C
  75.        IF(IXY.LT.3)  THEN
  76. c
  77. C...     DO  X  AND  Y  COMPONENTS.
  78.          IF((.NOT.NP1(IXY) .OR. (iand(1,I12).EQ.1)) .AND.
  79.      1   (.NOT.NP2(IXY) .OR. (iand(1,I34).EQ.1))) THEN
  80. C...       NO SIMPLIFICATIONS DUE TO PLANARITY
  81.            CALL MPPROD (BX(ijt),XLFT(NUMJ*NUCAQ+1,KQ,IXY),
  82.      *     NUCAQ, XRGT(NUMI*NUC12+1,IVA,IXY),NUCRS,XAND(IXY,KR,I34))
  83.            CALL CMPR1(BX(iJT),BBXS(NCQRS*NUM+1+(IXY-1)*NBB))
  84.            DO 10 JPQRS=1,NCQRS
  85.            BBXS(JPQRS+IB44) = BX(JPQRS+NCQRS) + BX(JPQRS)*
  86.      *     PWP(JPQRS,IXY)
  87.  10        continue
  88. C
  89.          ELSE
  90. C     PLANAR AB OR CD FRAGMENT
  91.           do i=1,ncqrs
  92.             bbxs(NCQRS*NUM+(IXY-1)*NBB+i) = zero
  93.           enddo
  94.           IF(IFPL(IXY))THEN
  95. c            CALL MOVE (BX(NCQRS+1),BBXS(IB44+1),NCQRS)
  96.             do jpqrs=1,ncqrs
  97.               bbxs(ib44+jpqrs) = bx(ncqrs+jpqrs)
  98.             enddo
  99.           ELSE
  100.             DO 13 JPQRS=1,NCQRS
  101.             BBXS(JPQRS+IB44) = BX(JPQRS+NCQRS) + BX(JPQRS)*
  102.      *      PWP(JPQRS,IXY)
  103.  13         continue
  104.           ENDIF
  105. C
  106.         ENDIF
  107. C
  108.       ELSE
  109. c
  110. C...  DO  Z  COMPONENTS
  111.         IF((.NOT.NP1(3) .OR. (iand(1,I12).EQ.1)) .AND.
  112.      1     (.NOT.NP2(3) .OR. (iand(1,I34).EQ.1))) THEN
  113.           CALL MPPROD (BX(iJT) ,XLFT(NUMJ*NUCAQ+1,KQ,3),NUCAQ ,
  114.      #    XRGT(NUMI*NUC12+1,IVA,3) , NUCRS, XAND(3,KR,I34))
  115.           DO 91 IJ=1,IM
  116.           IFB = NCQRS*IJ
  117.           IB2 = NCQRS*(NUMM+IJ-1) + 2*NBB
  118.           CALL CMPR1 (BX(iJT),BX(ijt+naqrs))
  119.           DO 90 JPQRS=1,NCQRS
  120.           BBXS(JPQRS+IB2) = BX(JPQRS+ijt+naqrs-1)*F(JPQRS+IFB-NCQRS)+
  121.      $    F(JPQRS+IFB)*(BX(JPQRS+NCQRS)+BX(JPQRS)*PWP(JPQRS,3))
  122.    90     CONTINUE
  123.  91       continue
  124. C
  125. C
  126.         ELSE IF (IFPL(3)) THEN
  127. C     PLANAR ABCD GEOMETRY
  128.           DO 961 IJ=1,IM
  129.           IFB = NCQRS*IJ
  130.           IB2 = NCQRS * (NUMM+IJ-1)  +  2*NBB
  131.           DO 96 JPQRS=1,NCQRS
  132.           BBXS(JPQRS+IB2) = F(JPQRS+IFB) * BX(JPQRS+NCQRS)
  133.  96       continue
  134.  961      continue
  135. C
  136.         ELSE
  137. C     PLANAR AB OR CD FRAGMENT
  138.           DO 971 IJ=1,IM
  139.           IFB = NCQRS*IJ
  140.           IB2 = NCQRS * (NUMM+IJ-1)  +  2*NBB
  141.           DO 97 JPQRS=1,NCQRS
  142.           BBXS(JPQRS+IB2) =
  143.      1    F(JPQRS+IFB)*(BX(JPQRS+NCQRS)+BX(JPQRS)*PWP(JPQRS,3))
  144.  97       continue
  145.  971      continue
  146.         ENDIF
  147. C
  148.       ENDIF
  149. C
  150. C...  NOW DO CONTRIBUTIONS FOR I>2
  151. C...
  152.       IF(IFPL(IXY)) THEN
  153.         IMAX = (JT+1)/2
  154.       ELSE
  155.         IMAX = JT
  156.       ENDIF
  157. C
  158.       DO 220 I=3,IMAX
  159.       IF(IFPL(IXY)) THEN
  160. C...    PLANAR ABCD GEOMETRY
  161.         IIX = (I+I-3)*NCQRS
  162.         DO 125 JPQRS=1,NCQRS
  163.         BX(JPQRS+ijt-1) = RFACT(I)*BX(JPQRS+IIX)
  164.  125    continue
  165.       ELSE
  166.         IRUM = MIN(I,JT-I+1)
  167.         IIX = NCQRS * (I-2)
  168.         DO 110 JPQRS=1,NCQRS
  169.         BX(JPQRS+iJT-1) = RFACT(I)*BX(JPQRS+IIX)*PWP(JPQRS+IIX,IXY)
  170.  110    continue
  171.         if(irum.lt.2) go to 160
  172.         DO 131 IRU=3,IRUM
  173.         IIX = (I+IRU-3) * NCQRS
  174.         IF(IRU.EQ.I) THEN
  175.           DO 126 JPQRS=1,NCQRS
  176.           BX(JPQRS+iJT-1) = BX(JPQRS+iJT-1) + RFACT(I)*BX(JPQRS+IIX)
  177.  126      continue
  178.         ELSE
  179.           IIY= (I-IRU-1) * NCQRS
  180.           DO 130 JPQRS=1,NCQRS
  181.           BX(JPQRS+iJT-1) = BX(JPQRS+iJT-1)+RFACT(IRU)*RFACT(I-IRU+1)*
  182.      1    PWP(JPQRS+IIY,IXY)*BX(JPQRS+IIX)
  183.  130      continue
  184.         ENDIF
  185.   131   CONTINUE
  186.         IIX = NCQRS *(I-1)
  187.         IIY = NCQRS * (I-3)
  188.         DO 150 JPQRS=1,NCQRS
  189.         BX(JPQRS+iJT-1) = BX(JPQRS+iJT-1) + RFACT(I-1)*
  190.      1  BX(JPQRS+IIX)*PWP(JPQRS+IIY,IXY)
  191.  150    continue
  192.       ENDIF
  193. C
  194.   160 IF(IXY .LT. 3)  THEN
  195. C...  X AND Y COMPONENTS
  196.         IB5 = NCQRS * (NUM+I-1) + (IXY-1)*NBB
  197.         DO 170 JPQRS=1,NCQRS
  198.         BBXS(JPQRS + IB5) = BX(JPQRS+iJT-1)
  199.  170    continue
  200.         ELSE
  201. C...  Z COMPONENTS
  202.           DO 211 IJ=1,IM
  203.           IFB = NCQRS*(IJ+I-2)
  204.           IB2 = NCQRS * (NUMM+IJ-1)  +  2*NBB
  205.           DO 210 JPQRS=1,NCQRS
  206.           BBXS(JPQRS+IB2) = BBXS(JPQRS+IB2)+F(JPQRS+IFB)*BX(JPQRS+iJT-1)
  207.  210      continue
  208.  211      continue
  209.         ENDIF
  210. C
  211.   220   CONTINUE
  212.   230   CONTINUE
  213. C
  214.       ELSE IF (JT.EQ.2) THEN
  215. C...  X AND Y COMPONENTS
  216.         IB8 = NUC12 * (NUMI+I34-1)
  217.         DO 270 IXY=1,2
  218.         IB44=IB4 + (IXY-1)*NBB
  219.         IF(IFPL(IXY)) GO TO 270
  220.         IF((.NOT.NP1(IXY) .OR. (I12.EQ.1)) .AND.
  221.      1   (.NOT.NP2(IXY) .OR. (I34.EQ.1))) THEN
  222.         CALL MVPROD(BX,XLFT(IB6+1,KQ,IXY),NUCAQ,XRGT(IB8+1,IVA,IXY),
  223.      #  NUCRS)
  224.         CALL MPPROD(BX(naqrs+1),XLFT(NUMJ*NUCAQ+1,KQ,IXY),NUCAQ,
  225.      #  XRGT(NUMI*NUC12+1,IVA,IXY),NUCRS,XAND(IXY,KR,I34))
  226.         CALL CMPR2(BX,BX(naqrs+1),BX(2*naqrs+1),
  227.      *             BBXS(NCQRS*NUM+1+(IXY-1)*NBB))
  228.         DO 260 JPQRS=1,NCQRS
  229.         BBXS(JPQRS+IB44) = BX(JPQRS+2*naqrs)*PWP(JPQRS,IXY)
  230.  260    continue
  231. C
  232.       ELSE
  233. C     PLANAR AB OR CD FRAGMENT
  234.         CALL MVPROD(BX(naqrs+1),XLFT(IB6+1,KQ,IXY),NUCAQ,
  235.      #  XRGT(IB8+1,IVA,IXY),NUCRS)
  236.         CALL CMPR1(BX(naqrs+1),BX)
  237.         do i=1,ncqrs
  238.           BBXS(NCQRS*NUM+(IXY-1)*NBB+i) = zero
  239.         enddo
  240.         DO 261 JPQRS=1,NCQRS
  241.         BBXS(JPQRS + IB44) = BX(JPQRS)*PWP(JPQRS,IXY)
  242.  261    continue
  243.       ENDIF
  244. C
  245.   270 CONTINUE
  246. C...  Z COMPONENTS
  247.       IF( .NOT.IFPL(3)) THEN
  248.       IF((.NOT.NP1(3) .OR. (I12.EQ.1)) .AND.
  249.      1   (.NOT.NP2(3) .OR. (I34.EQ.1))) THEN
  250.       CALL MVPROD(BX(3*naqrs+1),XLFT(IB6+1,KQ,3),NUCAQ,
  251.      * XRGT(IB8+1,IVA,3),NUCRS)
  252.       CALL MPPROD(BX(2*naqrs+1),XLFT(NUMJ*NUCAQ+1,KQ,3),NUCAQ,
  253.      $ XRGT(NUMI*NUC12+1,IVA,3),NUCRS,XAND(3,KR,I34))
  254.       CALL CMPR2(BX(2*naqrs+1),BX(3*naqrs+1),BX,BX(naqrs+1))
  255.       DO 285 JPQRS=1,NCQRS
  256.       BX(JPQRS+naqrs) = BX(JPQRS+naqrs) * PWP(JPQRS,3)
  257.  285  continue
  258.       DO 291 IJ=1,IM
  259.       IFB = NCQRS*IJ
  260.       IB2 = NCQRS *(NUMM+IJ-1)  +  2*NBB
  261.       DO 290 JPQRS=1,NCQRS
  262.       BBXS(JPQRS + IB2) = BX(JPQRS+naqrs) * F(JPQRS+IFB) +
  263.      1 BX(JPQRS) * F(JPQRS+IFB-NCQRS)
  264.  290  continue
  265.  291  continue
  266. C
  267.       ELSE
  268. C     PLANAR AB OR CD FRAGMENT
  269.       CALL MVPROD (BX(3*naqrs+1),XLFT(IB6+1,KQ,3),NUCAQ,
  270.      * XRGT(IB8+1,IVA,3),NUCRS)
  271.       CALL CMPR1(BX(3*naqrs+1),BX(naqrs+1))
  272.       DO 286 JPQRS=1,NCQRS
  273.       BX(JPQRS+naqrs) = BX(JPQRS+naqrs) * PWP(JPQRS,3)
  274.  286  continue
  275. c      GO TO 298
  276. c      DO 297 IJ=1,IM
  277. c      IFB = NCQRS*IJ
  278. c      IB2 = NCQRS * (NUMM+IJ-1)  +  2*NBB
  279. c      DO 291 JPQRS=1,NCQRS
  280. c  291 BBXS(JPQRS + IB2) = BX(JPQRS,2) * F(JPQRS+IFB)
  281. c  297 CONTINUE
  282.   298 CALL PROD(BBXS(NCQRS*NUMM+2*NBB+1),BX(naqrs+1),F(NCQRS+1),NCQRS,
  283.      * IM)
  284. C
  285. C     HAD TO DO THIS 'CAUSE OF A COMPILER BUG
  286. C     PROD  DOES WHAT THE LOOP ->297 WAS SUPPOSEDTO
  287. C     J.A. JULY -84
  288. C
  289.       ENDIF
  290. C
  291.         ENDIF
  292.       ELSE
  293. C...  SPECIAL CASE FOR JT=1
  294.         DO 320 JPQRS=1,NCQRS
  295.         BBXS(JPQRS) = 1.d0
  296.         BBXS(JPQRS+NBB) = 1.d0
  297.  320    continue
  298. c        CALL MOVE(F,BBXS(2*NBB+1),IM*NCQRS)
  299.         do i=1,im*ncqrs
  300.           bbxs(2*nbb+i) = f(i)
  301.         enddo
  302.       ENDIF
  303. C
  304.       NUMM=NUMM +IM
  305.       NUM=NUM+JT
  306.       NUMI=NUMI+I34
  307.  340  continue
  308.  341  continue
  309.       NUMJ=NUMJ+I12
  310.  1    continue
  311.  111  continue
  312. c
  313.       RETURN
  314.       END
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!