C N. R. BADNELL UoS v2.13 - QUB v1.2 27/08/09 C C SUBROUTINE BAKSUB(N,X,B,U,V,W) C C SOLVES FOR X THE VECTOR EQUATION UPP X=B C WHERE UPP IS AN UPPER TRIANGULAR MATRIX WITH ONLY THREE NON-ZERO C DIAGONALS, X AND B ARE COLUMN VECTORS. THE ARRAYS U, V AND W ARE C AS DESCRIBED IN SUBROUTINE VECTOR. C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION X(N),B(N),U(N),V(N),W(N) C X(N) = B(N)/U(N) X(N-1) = (B(N-1)-V(N-1)*X(N))/U(N-1) IF (N.EQ.2) GOTO 20 N1 = N - 1 DO 10 I = 2,N1 X(N-I) = (B(N-I)-V(N-I)*X(N-I+1)-W(N-I)*X(N-I+2))/U(N-I) 10 CONTINUE C 20 CONTINUE C END C C C BLOCK DATA BDLIB1 C*********************************************************************** C C BLOCK DATA C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C COMMON /BPSIZE/MXLR1,MXLR2,MXNC2,MXNR1,MXOCC COMMON /CONSTS/ZERO,TENTH,HALF,ONE,TWO,THREE,FOUR,SEVEN,ELEVEN,EPS COMMON /TERMS/NROWS,I(18),J(18),N(189) C C SET COMPILE PARAMETERS USED FOR STG2 AND RECUPD C DATA MXLR1,MXLR2,MXNC2,MXNR1,MXOCC/MZLR1,MZLR2,MZNC2,MZNR1,MZOCC/ C C SET GLOBAL REAL CONSTANTS C DATA ZERO,TENTH,HALF,ONE,TWO,THREE,FOUR,SEVEN,ELEVEN,EPS/ A 0.0D00,0.1D00,0.5D00,1.0D00, B 2.0D00,3.0D00,4.0D00,7.0D00, C 1.1D01,1.0D-05/ C C --- READS IN QUANTUM NUMBERS OF TERMS WHICH CAN BE FORMED FROM C CONFIGURATIONS L**Q . ONLY THE FIRST HALF OF THAT PART OF THE C TABLE, CORRESPONDING TO A GIVEN L, IS INCLUDED, BECAUSE OF THE C SYMMETRY OF THE TABLE. E.G. D**7 FORMS THE SAME TERMS AS D**3 C C THE ARRAYS I,J,N CORRESPOND TO THE ARRAYS ITAB,JTAB,NTAB C DATA NROWS/18/ DATA I(1),I(2),I(3),I(4),I(5),I(6)/1,1,1,3,3,1/ DATA I(7),I(8),I(9),I(10),I(11),I(12)/5,8,16,16,1,1/ DATA I(13),I(14),I(15),I(16),I(17),I(18)/1,1,1,1,1,1/ DATA J(1),J(2),J(3),J(4),J(5),J(6)/0,3,6,9,18,27/ DATA J(7),J(8),J(9),J(10),J(11),J(12)/30,45,69,117,165,168/ DATA J(13),J(14),J(15),J(16),J(17),J(18)/171,174,177,180,183,186/ DATA N(1),N(2),N(3),N(4),N(5),N(6)/1,1,2,0,1,1/ DATA N(7),N(8),N(9),N(10),N(11),N(12)/1,3,2,0,1,1/ DATA N(13),N(14),N(15),N(16),N(17),N(18)/2,5,1,2,3,3/ DATA N(19),N(20),N(21),N(22),N(23),N(24)/1,3,2,3,5,2/ DATA N(25),N(26),N(27),N(28),N(29),N(30)/3,1,4,1,5,2/ DATA N(31),N(32),N(33),N(34),N(35),N(36)/0,1,1,2,5,1/ DATA N(37),N(38),N(39),N(40),N(41),N(42)/2,9,1,2,3,3/ DATA N(43),N(44),N(45),N(46),N(47),N(48)/2,7,3,1,5,2/ DATA N(49),N(50),N(51),N(52),N(53),N(54)/3,3,2,3,5,2/ DATA N(55),N(56),N(57),N(58),N(59),N(60)/3,7,2,3,9,2/ DATA N(61),N(62),N(63),N(64),N(65),N(66)/3,11,2,3,3,4/ DATA N(67),N(68),N(69),N(70),N(71),N(72)/3,7,4,0,1,1/ DATA N(73),N(74),N(75),N(76),N(77),N(78)/2,5,1,2,9,1/ DATA N(79),N(80),N(81),N(82),N(83),N(84)/2,3,3,2,7,3/ DATA N(85),N(86),N(87),N(88),N(89),N(90)/4,1,1,4,5,1/ DATA N(91),N(92),N(93),N(94),N(95),N(96)/4,7,1,4,9,1/ DATA N(97),N(98),N(99),N(100),N(101),N(102)/4,13,1,4,3,3/ DATA N(103),N(104),N(105),N(106),N(107),N(108)/4,5,3,4,7,3/ DATA N(109),N(110),N(111),N(112),N(113),N(114)/4,9,3,4,11,3/ DATA N(115),N(116),N(117),N(118),N(119),N(120)/4,5,5,1,5,2/ DATA N(121),N(122),N(123),N(124),N(125),N(126)/3,3,2,3,5,2/ DATA N(127),N(128),N(129),N(130),N(131),N(132)/3,7,2,3,9,2/ DATA N(133),N(134),N(135),N(136),N(137),N(138)/3,11,2,3,3,4/ DATA N(139),N(140),N(141),N(142),N(143),N(144)/3,7,4,5,1,2/ DATA N(145),N(146),N(147),N(148),N(149),N(150)/5,5,2,5,7,2/ DATA N(151),N(152),N(153),N(154),N(155),N(156)/5,9,2,5,13,2/ DATA N(157),N(158),N(159),N(160),N(161),N(162)/5,5,4,5,9,4/ DATA N(163),N(164),N(165),N(166),N(167),N(168)/5,1,6,1,7,2/ DATA N(169),N(170),N(171)/1,9,2/ DATA N(172),N(173),N(174),N(175),N(176),N(177)/1,11,2,1,13,2/ DATA N(178),N(179),N(180),N(181),N(182),N(183)/1,15,2,1,17,2/ DATA N(184),N(185),N(186),N(187),N(188),N(189)/1,19,2,1,21,2/ C END C C C BLOCK DATA BDLIB2 C*********************************************************************** C C INSERT THE PROGRAM PSHELL C F P BY D C S ALLISON CATALOGUE C NUMBER ACQB (CPC 1 (1969) 15) AND THE PROGRAM A NEW D SHELL C F P C BY A T CHIVERS CATALOGUE NUMBER ACRN (CPC 6 (1973) 88) C C BLOCK DATA,CFPP,CFPD C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C COMMON /FRPAR2/I(719) C C BLOCK DATA FOR CFPD SUBROUTINE C DATA I(1),I(2),I(3),I(4),I(5),I(6),I(7),I(8),I(9),I(10),I(11), A I(12),I(13),I(14),I(15),I(16),I(17),I(18),I(19),I(20),I(21), B I(22),I(23),I(24),I(25),I(26),I(27),I(28),I(29),I(30),I(31), C I(32),I(33),I(34),I(35),I(36),I(37),I(38),I(39),I(40),I(41), D I(42),I(43),I(44),I(45),I(46),I(47),I(48),I(49),I(50),I(51), E I(52),I(53),I(54),I(55),I(56),I(57),I(58),I(59),I(60),I(61), F I(62),I(63),I(64),I(65),I(66),I(67),I(68),I(69),I(70),I(71), G I(72),I(73),I(74),I(75),I(76),I(77),I(78),I(79),I(80)/1,5,8, H 16,16,1,2,3,4,5,0,2,3,4,5,0,2,3,4,3,0,2,3,2,5,0,0,3,4,3,0,0, I 1,4,5,0,0,3,2,3,0,0,3,4,3,0,0,0,4,5,0,0,0,2,3,0,0,0,4,5,0,0, J 0,4,1,0,0,0,2,3,0,0,0,4,5,0,0,0,0,3/ DATA I(81),I(82),I(83),I(84),I(85),I(86),I(87),I(88),I(89),I(90), A I(91),I(92),I(93),I(94),I(95),I(96),I(97),I(98),I(99),I(100), B I(101),I(102),I(103),I(104),I(105),I(106),I(107),I(108), C I(109),I(110),I(111),I(112),I(113),I(114),I(115),I(116), D I(117),I(118),I(119),I(120),I(121),I(122),I(123),I(124), E I(125),I(126),I(127),I(128),I(129),I(130),I(131),I(132), F I(133),I(134),I(135),I(136),I(137),I(138),I(139),I(140), G I(141),I(142),I(143),I(144),I(145)/0,0,0,4,5,2,3,3,2,0,0,1,1, H 5,4,0,4,5,4,3,0,2,4,3,2,0,0,3,3,1,0,0,2,2,6,0,0,2,1,5,0,0,1, I 1,4,0,0,0,6,4,0,0,0,4,3,0,0,0,4,3,0,0,0,3,2/ DATA I(146),I(147),I(148),I(149),I(150),I(151),I(152),I(153), A I(154),I(155),I(156),I(157),I(158),I(159),I(160),I(161), B I(162),I(163),I(164),I(165),I(166),I(167),I(168),I(169), C I(170),I(171),I(172),I(173),I(174),I(175),I(176),I(177), D I(178),I(179),I(180),I(181),I(182),I(183),I(184),I(185), E I(186),I(187),I(188),I(189),I(190),I(191),I(192),I(193), F I(194),I(195),I(196),I(197),I(198),I(199),I(200),I(201), G I(202),I(203),I(204),I(205),I(206),I(207),I(208),I(209), H I(210),I(211),I(212),I(213),I(214),I(215),I(216),I(217), I I(218),I(219),I(220),I(221),I(222),I(223),I(224),I(225)/0,0, J 0,2,2,0,0,0,2,2,0,0,0,0,1,0,0,0,0,0,2,3,4,5,6,0,3,4,3,4,0,1, K 2,3,4,0,1,2,3,4,0,1,2,3,4,0,0,2,3,2,0,0,2,3,2,0,0,2,3,2,0,0, L 0,1,2,0,0,0,1,2,0,0,0,1,2,0,0,0,1,2/ DATA I(226),I(227),I(228),I(229),I(230),I(231),I(232),I(233), A I(234),I(235),I(236),I(237),I(238),I(239),I(240),I(241), B I(242),I(243),I(244),I(245),I(246),I(247),I(248),I(249), C I(250),I(251),I(252),I(253),I(254),I(255),I(256),I(257), D I(258),I(259),I(260),I(261),I(262),I(263),I(264),I(265), E I(266),I(267),I(268),I(269),I(270),I(271),I(272),I(273), F I(274),I(275),I(276),I(277),I(278),I(279),I(280),I(281), G I(282),I(283),I(284),I(285),I(286),I(287),I(288),I(289), H I(290)/0,0,0,1,2,0,0,0,1,2,0,0,0,1,2,0,0,0,1,2,1,1,1,1,1,4, I -7,-1,21,7,-21,21,-8,-1,-8,0,0,28,-9,-49,7,0,0,1,11,-25,-9, J -25,0,0,0,0,-10,-10,-5,45,15,0,0,0,0,0,16,0,0/ DATA I(291),I(292),I(293),I(294),I(295),I(296),I(297),I(298), A I(299),I(300),I(301),I(302),I(303),I(304),I(305),I(306), B I(307),I(308),I(309),I(310),I(311),I(312),I(313),I(314), C I(315),I(316),I(317),I(318),I(319),I(320),I(321),I(322), D I(323),I(324),I(325),I(326),I(327),I(328),I(329),I(330), E I(331),I(332),I(333),I(334),I(335),I(336),I(337),I(338), F I(339),I(340),I(341),I(342),I(343),I(344),I(345),I(346), G I(347),I(348),I(349),I(350),I(351),I(352),I(353),I(354), H I(355),I(356),I(357),I(358),I(359),I(360),I(361),I(362), I I(363),I(364),I(365),I(366),I(367),I(368),I(369),I(370)/7,20, J -560,224,-112,-21,-56,16,0,0,0,0,0,0,0,0,3,0,0,-56,-448,49, K -64,-14,0,0,0,0,0,0,0,0,0,26,308,110,220,0,0,0,7,-154,-28, L -132,0,0,0,0,0,-9,297,90,-405,45,0,0,3,66,-507,-3,-60,15,0,0, M 0,5,315,-14,-175,-21,-56,-25,0,70,385,-105,28,63,0,0/ DATA I(371),I(372),I(373),I(374),I(375),I(376),I(377),I(378), A I(379),I(380),I(381),I(382),I(383),I(384),I(385),I(386), B I(387),I(388),I(389),I(390),I(391),I(392),I(393),I(394), C I(395),I(396),I(397),I(398),I(399),I(400),I(401),I(402), D I(403),I(404),I(405),I(406),I(407),I(408),I(409),I(410), E I(411),I(412),I(413),I(414),I(415),I(416),I(417),I(418), F I(419),I(420),I(421),I(422),I(423),I(424),I(425),I(426), G I(427),I(428),I(429),I(430),I(431),I(432),I(433),I(434), H I(435)/0,0,0,315,0,0,135,0,0,189,0,0,105,0,1,0,0,0,200,15, I 120,60,-35,10,0,-25,88,200,45,20,0,1,0,0,0,16,-200,-14,-14, J 25,0,0,0,120,-42,42,0,0,1,-105,-175,-175,-75,0,0,0,0,0,0,0,0, K 0,0,0,0/ DATA I(436),I(437),I(438),I(439),I(440),I(441),I(442),I(443), A I(444),I(445),I(446),I(447),I(448),I(449),I(450),I(451), B I(452),I(453),I(454),I(455),I(456),I(457),I(458),I(459), C I(460),I(461),I(462),I(463),I(464),I(465),I(466),I(467), D I(468),I(469),I(470),I(471),I(472),I(473),I(474),I(475), E I(476),I(477),I(478),I(479),I(480),I(481),I(482),I(483), F I(484),I(485),I(486),I(487),I(488),I(489),I(490),I(491), G I(492),I(493),I(494),I(495),I(496),I(497),I(498),I(499), H I(500),I(501),I(502),I(503),I(504),I(505),I(506),I(507), I I(508),I(509),I(510),I(511),I(512),I(513),I(514),I(515)/154, J -110,0,0,231,286,924,-308,220,-396,0,0,0,0,0,0,-66,-90,180,0, K 99,-99,891,-5577,-405,-9,0,45,45,0,0,0,0,224,0,-56,0,-220, L 1680,0,112,0,-21,21,0,-16,0,0,-70,14,-84,56,0,55,945,4235, M -175,-315,0,-21,189,-25,0,0,25,-15,-135,35,0,0,600,968,120, N 600,0,60,60,10,3,0/ DATA I(516),I(517),I(518),I(519),I(520),I(521),I(522),I(523), A I(524),I(525),I(526),I(527),I(528),I(529),I(530),I(531), B I(532),I(533),I(534),I(535),I(536),I(537),I(538),I(539), C I(540),I(541),I(542),I(543),I(544),I(545),I(546),I(547), D I(548),I(549),I(550),I(551),I(552),I(553),I(554),I(555), E I(556),I(557),I(558),I(559),I(560),I(561),I(562),I(563), F I(564),I(565),I(566),I(567),I(568),I(569),I(570),I(571), G I(572),I(573),I(574),I(575),I(576),I(577),I(578),I(579), H I(580)/0,-56,0,-64,0,0,0,0,448,0,-9,-49,0,14,0,0,0,-16,126, I 14,0,0,0,0,-200,360,0,-14,126,25,0,0,0,0,0,0,-175,182,-728, J -2184,0,0,0,0,0,0,0,0,0,0,0,0,0,220,880,0,-400,0,-9,-25,0,0, K 0,0,0/ DATA I(581),I(582),I(583),I(584),I(585),I(586),I(587),I(588), A I(589),I(590),I(591),I(592),I(593),I(594),I(595),I(596), B I(597),I(598),I(599),I(600),I(601),I(602),I(603),I(604), C I(605),I(606),I(607),I(608),I(609),I(610),I(611),I(612), D I(613),I(614),I(615),I(616),I(617),I(618),I(619),I(620), E I(621),I(622),I(623),I(624),I(625),I(626),I(627),I(628), F I(629),I(630),I(631),I(632),I(633),I(634),I(635),I(636), G I(637),I(638),I(639),I(640),I(641),I(642),I(643),I(644), H I(645),I(646),I(647),I(648),I(649),I(650),I(651),I(652), I I(653),I(654),I(655),I(656),I(657),I(658),I(659),I(660)/0,0, J 0,-45,-5,845,-1215,275,495,0,-11,99,0,0,0,0,0,0,0,0,33,-7, K -2541,105,-525,0,35,35,-15,0,0,0,0,0,0,0,0,-800,0,-160,0,-5, L 45,0,30,0,0,0,0,0,0,0,0,-100,1452,180,-100,0,-10,90,15,-2,0, M 0,0,0,0,0,0,0,0,0,0,6,0,0,0,0,0,0/ DATA I(661),I(662),I(663),I(664),I(665),I(666),I(667),I(668), A I(669),I(670),I(671),I(672),I(673),I(674),I(675),I(676), B I(677),I(678),I(679),I(680),I(681),I(682),I(683),I(684), C I(685),I(686),I(687),I(688),I(689),I(690),I(691),I(692), D I(693),I(694),I(695),I(696),I(697),I(698),I(699),I(700), E I(701),I(702),I(703),I(704),I(705),I(706),I(707),I(708), F I(709),I(710),I(711),I(712),I(713),I(714),I(715),I(716), G I(717),I(718),I(719)/0,0,0,0,0,0,0,0,0,0,-14,-56,0,0,1,1,1,1, H 1,5,15,2,42,70,60,140,30,10,60,1680,840,1680,210,360,90,10, I 504,1008,560,280,140,1,1,1,420,700,700,300,550,1100,8400, J 18480,2800,2800,50,350,700,150,5/ C END C C C SUBROUTINE BUBBLE(JPOL,FAIL) C C REDUCES A CIRCUIT OF ORDER 2,GIVING DELTA FUNCTION AND PHASE C FACTORS. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C CHARACTER*6 NAME,NAMSUB COMMON /NAM/NAMSUB DATA NAME/'BUBBLE'/ C NAMSUB = NAME K2 = 2 K23 = 3 I1 = 1 I2 = 1 IT1 = NPOINT(1) IT2 = NPOINT(2) IF (IT2.NE.ILAST) GOTO 20 IF (IT1.EQ.IFIRST) GOTO 10 IT2 = IT1 IT1 = ILAST 10 CONTINUE I1 = -1 K23 = 2 I2 = 2 C C PHASE(IT1,JDIAG,KFL2B): C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C 20 CONTINUE J7(J7C+1) = JDIAG(IT1,1) J7(J7C+2) = JDIAG(IT1,2) J7C = J7C + 3 J7(J7C) = JDIAG(IT1,3) C K = ABS((3*ARR(IT2,1)+2*ARR(IT2,2)+ARR(IT2,3))/2) + 1 IF (K.NE.4) THEN C C PHASE2(JDIAG(IT2,K)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JDIAG(IT2,K) C ENDIF C IF (NBNODE.EQ.2) GOTO 80 IL1 = IL(IT2) + I1 IT = IH(IL1) ARR(IT,K23) = ARR(IT1,K23) L = JDIAG(IT1,K23) L1 = JDIAG(IT,K23) JDIAG(IT,K23) = L IF (JPOL.NE.1) GOTO 30 MP = MP - 1 JW(2,JWC) = L J6(J6C-1) = L J6(J6C) = L IF (K.EQ.2) J8(J8C) = L GOTO 40 C 30 CONTINUE CALL DELTA(L,L1,FAIL) IF (FAIL) GOTO 80 40 CONTINUE TAB1(L,I2) = IT IF (IT1.EQ.ILAST) GOTO 70 IF (IT2.NE.ILAST) GOTO 50 TAB1(L,1) = IH(2) IL1 = 2 K2 = 1 50 CONTINUE DO 60 I = IL1,NBNODE IT = IH(I) IL(IT) = I - K2 IH(I-K2) = IT 60 CONTINUE 70 CONTINUE J9(J9C+1) = L J9C = J9C + 2 J9(J9C) = L 80 CONTINUE C END C C C SUBROUTINE CFP(LIJ,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) C*********************************************************************** C C CFP C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C C === CHOOSES APPROPRIATE FRACTIONAL PARENTAGE SUBROUTINE C C IF F-SHELL OR G-SHELL COEFFICIENT-OF-FRACTIONAL-PARENTAGE ROUTINES C ARE INCLUDED, THIS COMPUTED GOTO NEEDS MODIFYING TO ACCOUNT FOR C THIS C COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 C C --- FALSE CALL FOR S-SHELLS C COEFP = 1.0D0 C C --- P-SHELLS C IF (LIJ.EQ.1) CALL CFPP(N,ILI,ISI,ILJ,ISJ,COEFP) C C --- D-SHELLS C IF (LIJ.EQ.2) CALL CFPD(N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) C C --- F-SHELLS, G-SHELLS ETC. WITH UP TO TWO ELECTRONS C C IF(LIJ.GE.3) CALL CFPF(N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) C IF (IBUG1.GT.1) WRITE (IWRITE,3000) COEFP C 3000 FORMAT (' COEFP =',F15.9) END C C C SUBROUTINE CFPD(N,IVI,LI,ISI,IVJ,LJ,ISJ,COEFP) IMPLICIT REAL*8 (A-H,O-Z) C C THIS SUBROUTINE EVALUATES THE COEFFICIENTS OF FRACTIONAL PARENTAGE C FOR EQUIVALENT D SHELL ELECTRONS FROM TABLES GIVEN IN J.C.SLATER C QUANTUM THEORY OF ATOMIC STRUCTURE,VOLUME2,P350(1960) C IN THE SUBROUTINE LIST N,THE NO.OF ELECTRONS,V THE SENIORITY QUAN C TUM NO.,L THE ANGULAR MOMENTUM QUANTUM NO.,(2S+1) THE SPIN QUANTUM C NO. OF BOTH THE STATE IN QUESTION AND ITS PARENT STATE ARE INPUT C PARAMETERS THE RESULT IS OUTPUT AS COEFP C COMMON /FRPAR2/K(5),IV(5,16),IL(5,16),IS(5,16),ITAB1(5,1), A ITAB2(8,5),ITAB3(16,8),ITAB4(16,16),NORM1(5),NORM2(8), B NORM3(16),NORM4(16) C C TEST IF N IS IN THE FIRST HALF OF SHELL C COEFP = 1.0D0 IF (N.GE.6) GOTO 30 IF (N.EQ.1) THEN IF (IVJ.EQ.0 .AND. LJ.EQ.0 .AND. ISJ.EQ.1) GOTO 80 GOTO 70 C ENDIF C C TEST IF STATE IN QUESTION IS ALLOWED C IF IT IS, IDENTIFY THE ROW OF THE TABLE BY J1 C DO 10 J = 1,16 J1 = J IF (IV(N,J).EQ.IVI .AND. IL(N,J).EQ.LI .AND. A IS(N,J).EQ.ISI) GOTO 15 10 CONTINUE GOTO 70 C C TEST IF PARENT STATE IS ALLOWED C IF IT IS, IDENTIFY THE COLUMN OF THE TABLE BY J2 C 15 CONTINUE DO 20 J = 1,16 J2 = J IF (IV(N-1,J).EQ.IVJ .AND. IL(N-1,J).EQ.LJ .AND. A IS(N-1,J).EQ.ISJ) GOTO 60 20 CONTINUE GOTO 70 C C SIMILAR SETTING OF J1 AND J2 IF N IS IN SECOND HALF OF SHELL C 30 CONTINUE M = 10 - N IF (N.EQ.10) THEN IF (IVI.EQ.0 .AND. LI.EQ.0 .AND. ISI.EQ.1) GOTO 80 GOTO 70 C ENDIF C DO 40 J = 1,16 J2 = J IF (IV(M,J).EQ.IVI .AND. IL(M,J).EQ.LI .AND. A IS(M,J).EQ.ISI) GOTO 45 40 CONTINUE GOTO 70 C 45 CONTINUE DO 50 J = 1,16 J1 = J IF (IV(M+1,J).EQ.IVJ .AND. IL(M+1,J).EQ.LJ .AND. A IS(M+1,J).EQ.ISJ) GOTO 60 50 CONTINUE GOTO 70 C C IDENTIFY THE F.P.C AS A UNIQUE ELEMENT OF ITABN(J1,J2) C 60 CONTINUE M = MIN(N,11-N) IF (M.EQ.2) THEN IN = ITAB1(J1,J2) AN = NORM1(J1) C ELSE IF (M.EQ.3) THEN IN = ITAB2(J1,J2) AN = NORM2(J1) C ELSE IF (M.EQ.4) THEN IN = ITAB3(J1,J2) AN = NORM3(J1) C ELSE IF (M.EQ.5) THEN IN = ITAB4(J1,J2) AN = NORM4(J1) C ELSE GOTO 80 C ENDIF C COEFP = ABS(IN) IF (IN.EQ.0) GOTO 80 COEFP = SQRT(COEFP/AN) IF (IN.LT.0) COEFP = -COEFP IF (N.LE.5) GOTO 80 C C USE RECURRENCE RELATION EQUATION (19) OF RACAH FOR SECOND HALF OF C SHELL C ISIGN = (-1)** ((ISI+ISJ-7)/2+LI+LJ) FACTOR = SQRT( A DBLE((11-N)*ISJ* (2*LJ+1)) / B DBLE(N*ISI* (2*LI+1)) ) COEFP = COEFP*ISIGN*FACTOR IF (N.EQ.6 .AND. MOD((IVJ-1)/2,2).NE.0) COEFP = -COEFP GOTO 80 C C AN UNALLOWED STATE C FOR AN UNALLOWED STATE THE F.P. COEFFICIENT IS SET EQUAL TO AN C ERRONEOUS VALUE.BY REPLACING THE 11 COEFP=9.9 STATEMENT BY THE C FOLLWING 3 CARDS THE USER CAN TERMINATE THE RUN WHEN AN C UNALLOWED STATE OCCURS C 106 FORMAT(' FAIL IN CFPD AT 11 UNALLOWED STATE') C 11 WRITE(IWRITE,106) C STOP C 70 CONTINUE COEFP = 9.9D0 80 CONTINUE C END C C C SUBROUTINE CFPP(N,LI,ISI,LJ,ISJ,COEFP) IMPLICIT REAL*8 (A-H,O-Z) C C THIS SUBROUTINE EVALUATES THE COEFFICIENTS OF FRACTIONAL PARENTAGE C FOR EQUIVALENT P SHELL ELECTRONS FROM TABLES GIVEN IN J.C.SLATER C QUANTUM THEORY OF ATOMIC STRUCTURE,VOLUME2,P350(1960) C IN THE SUBROUTINE LIST N,THE NO. OF ELECTRONS,L THE ANGULAR C MOMENTUM QUANTUM NO.,(2S+1) THE SPIN QUANTUM NO. OF BOTH THE STATE C IN QUESTION AND ITS PARENT STATE ARE INPUT PARAMETERS.THE RESULT C IS OUTPUT AS COEFP C DIMENSION IL(3,3),IS(3,3),ITAB1(3,1),ITAB2(3,3),NORM1(3),NORM2(3) C C SET UP P SHELL PARAMETERS AND TABLES C DATA IL(1,1),IL(2,1),IL(2,2),IL(2,3),IL(3,1),IL(3,2),IL(3,3)/1,1, A 2,0,0,2,1/ DATA IS(1,1),IS(2,1),IS(2,2),IS(2,3),IS(3,1),IS(3,2),IS(3,3)/2,3, A 1,1,4,2,2/ DATA ITAB1(1,1),ITAB1(2,1),ITAB1(3,1)/1,1,1/ DATA ITAB2(1,1),ITAB2(1,2),ITAB2(1,3),ITAB2(2,1),ITAB2(2,2), A ITAB2(2,3),ITAB2(3,1),ITAB2(3,2),ITAB2(3,3)/1,0,0,1,-1,0,-9, B -5,4/ DATA NORM1(1),NORM1(2),NORM1(3)/1,1,1/ DATA NORM2(1),NORM2(2),NORM2(3)/1,2,18/ C C TEST IF N IS IN THE FIRST HALF OF SHELL C COEFP = 1.0D0 IF (N.GE.4) GOTO 30 IF (N.EQ.1) THEN IF (LJ.EQ.0 .AND. ISJ.EQ.1) GOTO 90 GOTO 80 C ENDIF C C TEST IF STATE IN QUESTION IS ALLOWED C IF IT IS, IDENTIFY THE ROW OF THE TABLE BY J1 C DO 10 J = 1,3 J1 = J IF (IL(N,J).EQ.LI .AND. IS(N,J).EQ.ISI) GOTO 15 10 CONTINUE GOTO 80 C C TEST IF PARENT STATE IS ALLOWED C IF IT IS, IDENTIFY THE COLUMN OF THE TABLE BY J2 C 15 CONTINUE DO 20 J = 1,3 J2 = J IF (IL(N-1,J).EQ.LJ .AND. IS(N-1,J).EQ.ISJ) GOTO 70 20 CONTINUE GOTO 80 C C SIMILAR SETTING OF J1 AND J2 IF N IS IN SECOND HALF OF SHELL C 30 CONTINUE M = 6 - N IF (N.EQ.6) THEN IF (LI.EQ.0 .AND. ISI.EQ.1) GOTO 90 GOTO 80 C ENDIF C DO 40 J = 1,3 J2 = J IF (IL(M,J).EQ.LI .AND. IS(M,J).EQ.ISI) GOTO 50 40 CONTINUE GOTO 80 C 50 CONTINUE DO 60 J = 1,3 J1 = J IF (IL(M+1,J).EQ.LJ .AND. IS(M+1,J).EQ.ISJ) GOTO 70 60 CONTINUE GOTO 80 C C IDENTIFY THE F.P.C AS A UNIQUE ELEMENT OF ITABN(J1,J2) C 70 CONTINUE M = MIN(N,7-N) IF (M.EQ.2) THEN IN = ITAB1(J1,J2) AN = NORM1(J1) C ELSE IF (M.EQ.3) THEN IN = ITAB2(J1,J2) AN = NORM2(J1) C ELSE GOTO 90 C ENDIF C COEFP = ABS(IN) IF (IN.EQ.0) GOTO 90 COEFP = SQRT(COEFP/AN) IF (IN.LT.0) COEFP = -COEFP IF (N.LE.3) GOTO 90 C C USE RECURRENCE RELATION EQUATION (19) OF RACAH FOR SECOND HALF OF C SHELL C ISIGN = (-1)** ((ISI+ISJ-5)/2+LI+LJ) FACTOR = A DBLE((7-N)*ISJ* (2*LJ+1)) / B DBLE(N*ISI* (2*LI+1)) COEFP = COEFP*ISIGN*SQRT(FACTOR) IF (N.EQ.4 .AND. LJ.NE.1) COEFP = -COEFP GOTO 90 C C FOR AN UNALLOWED STATE THE F.P.COEFFICIENT IS SET EQUAL TO AN C ERRONEOUS VALUE. THIS STATEMENT COULD BE REPLACED BY AN ERROR C STATEMENT C 80 CONTINUE COEFP = 9.9D0 90 CONTINUE C END C C C FUNCTION CG(J1,J2,J3,M1,M2,M3) C*********************************************************************** C C CG C C*********************************************************************** C C CALCULATES A CLEBSCH-GORDAN COEFFICIENT. C C J1,J2,J3 ARE THE J VALUES OF THE COEFFICIENT C M1,M2,M3 ARE THE M VALUES OF THE COEFFICIENT C C FORMULA FOR VECTOR COUPLING COEFFICIENT C (CONDON AND SHORTLEY PAGE 75 FORMULA (3.14.5)) C IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C COMMON /FACT/GAMMA(MZFAC) C CG = 0.0D0 IF (ABS(M1).GT.J1 .OR. ABS(M2).GT.J2 .OR. A ABS(M3).GT.J3) GOTO 20 IF (M1+M2.NE.M3) GOTO 20 IF (2*MAX(J1,J2,J3).GT.J1+J2+J3) GOTO 20 C XNUM = (GAMMA(-J1+J2+J3+1)/GAMMA(J2+M2+1))* A (GAMMA(J1-J2+J3+1)/GAMMA(J1+M1+1))* B (GAMMA(J1+J2-J3+1)/GAMMA(J1-M1+1))* C (GAMMA(J3-M3+1)/GAMMA(J2-M2+1))* D ((2*J3+1)*GAMMA(J3+M3+1)/GAMMA(J1+J2+J3+2)) XNUM = SQRT(XNUM) IA = -J1 + M1 IB = J2 - J1 + M3 IC = J2 + J3 + M1 ID = -J1 + J2 + J3 IE = J3 + M3 NUMIN = MAX(0,IB) NUMAX = MIN(IC,ID,IE) J = J2 + M2 C DO 10 N = NUMIN,NUMAX G = (GAMMA(N-IA+1)*GAMMA(IC-N+1))/ A (GAMMA(ID-N+1)*GAMMA(IE-N+1)*GAMMA(N-IB+1)*GAMMA(N+1)) IF (MOD(J+N,2).NE.0) G = -G CG = CG + G 10 CONTINUE C CG = CG*XNUM 20 CONTINUE C END C C C SUBROUTINE CHANGE(L,K) C C EXCHANGES THE FREE ENDS IN EITHER FIRST OR LAST TRIAD OF JDIAG. C C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4,KFL2A=2*KFL2,KFL2B=4*KFL2, A KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C C PHASE(L,JDIAG,KFL2B): C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C J7(J7C+1) = JDIAG(L,1) J7(J7C+2) = JDIAG(L,2) J7C = J7C + 3 J7(J7C) = JDIAG(L,3) C JP = JDIAG(L,K) JDIAG(L,K) = JDIAG(L,1) JDIAG(L,1) = JP JAR = ARR(L,K) ARR(L,K) = ARR(L,1) ARR(L,1) = JAR C END C C C SUBROUTINE CHOP C*********************************************************************** C C CHOP C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) COMMON /REMOVE/ICHOP(MXORB2) COMMON /DIAGNL/IDIAG,JA,JB C C --- ZEROIZE THE OUTPUT ARRAY C DO 10 I = 1,IHSH ICHOP(I) = 0 10 CONTINUE C C NO AVERAGE ENERGY TERMS FOR OFF-DIAGONAL MATRIX ELEMENTS C IF (IDIAG.EQ.0) RETURN JSTO = 0 ICOUNT = 0 DO 20 J = 1,IHSH NFULL = 4*LJ(J) + 2 I2 = NOSH1(J) C C IS THE SHELL FULL OR EMPTY C IF (I2.NE.NFULL .AND. I2.NE.0) THEN C C IF NOT, DOES IT CONTAIN ONLY ONE ELECTRON, OR ONLY ONE =HOLE= C IF (I2.EQ.1 .OR. I2.EQ. (NFULL-1)) JSTO = J C ELSE C C ICHOP SET UNITY FOR CLOSED SHELLS C ICHOP(J) = 1 ICOUNT = ICOUNT + 1 ENDIF C 20 CONTINUE C C IF ALL BUT ONE SHELL IS CLOSED, AND THIS CONTAINS ONE ELECTRON OR C =HOLE= , THEN IT CAN BE TREATED PURELY BY AVERAGE ENERGY C IF (ICOUNT.NE. (IHSH-1) .OR. JSTO.EQ.0) RETURN ICHOP(JSTO) = 1 C END C C C SUBROUTINE CHVAR(JP,NBC,KBC,JT,JINV,NSUM) C C CHANGES THE ORDER OF SUMMATION VARIABLE TO BE ABLE TO PERFORM C SEPARATELY THE SUMMATIONS IN GENSUM. C LOGICAL JT(NSUM) DIMENSION JP(NBC),JINV(NSUM) C KB = KBC + 1 IF (KB.GT.NBC) GOTO 20 DO 10 I = KB,NBC JK = JP(I) IF (.NOT.JT(JK)) GOTO 10 KBC = KBC + 1 JP(I) = JP(KBC) JP(KBC) = JINV(JK) 10 CONTINUE 20 CONTINUE C END C C C SUBROUTINE CUT1L(FAIL) C C CUT ON ONE LINE,THAT WAS LEFT AS A FREE END IN JDIAG.PUTS C CORRESPONDING DELTA IN J23. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C CHARACTER*6 NAME DATA NAME/'CUT1L '/ C IT = ITFREE(1) J0 = JDIAG(IT,1) CALL DELTA(J0,M,FAIL) IF (FAIL) GOTO 50 CALL DELTA(JDIAG(IT,3),JDIAG(IT,2),FAIL) IF (FAIL) GOTO 50 JDIAG(IT+1,3) = JDIAG(IT,3) IF (ARR(IT,2)-ARR(IT,3)) 20,10,30 10 CONTINUE ARR(IT+1,3) = 1 ARR(IT-1,2) = -1 GOTO 30 C 20 CONTINUE ARR(IT+1,3) = -1 ARR(IT-1,2) = 1 30 CONTINUE J9C = J9C + 1 J9(J9C) = JDIAG(IT,3) J = 2 CALL ZERO(J,J0,FAIL) IF (FAIL) GOTO 50 IL1 = IL(IT+1) DO 40 I = IL1,NBNODE IT = IH(I) ILP = I - 1 IL(IT) = ILP IH(ILP) = IT 40 CONTINUE NBNODE = NBNODE - 1 50 CONTINUE CALL PRINTJ(NAME,12) C END C C C SUBROUTINE CUT2L(FAIL) C C CUT ON TWO LINES THAT WERE LEFT AS FREE ENDS IN JDIAG.PUTS C CORRESPONDING DELTA IN J23. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C LOGICAL TABS INTEGER ARROW COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C CHARACTER*6 NAME DATA NAME/'CUT2L '/ C IT1 = ITFREE(1) IT2 = ITFREE(2) JT1 = JDIAG(IT1,1) JT2 = JDIAG(IT2,1) CALL DELTA(JT1,JT2,FAIL) IF (FAIL) GOTO 10 IF (ARR(IT1,1).EQ.ARR(IT2,1)) THEN C C PHASE2(JT1): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JT1 C ENDIF C ARR(IT2,1) = -ARR(IT1,1) JDIAG(IT2,1) = JT1 TAB1(JT1,2) = IT2 J9(J9C+1) = JT1 J9C = J9C + 2 J9(J9C) = JT1 C C OTHERJ(0,JT1,L1,LC1,K1): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C L1 = LINE(JT1,1) IF (L1.EQ.0 .OR. TABS(L1)) THEN K1 = 1 L1 = LINE(JT1,2) LC1 = LCOL(JT1,2) C ELSE K1 = 2 LC1 = LCOL(JT1,1) ENDIF C C OTHERJ(0,JT2,L2,LC2,K2): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C L2 = LINE(JT2,1) IF (L2.EQ.0 .OR. TABS(L2)) THEN L2 = LINE(JT2,2) LC2 = LCOL(JT2,2) C ELSE LC2 = LCOL(JT2,1) ENDIF C J23(L2,LC2) = JT1 LINE(JT1,K1) = L2 LCOL(JT1,K1) = LC2 ARROW(L2,LC2) = -ARROW(L1,LC1) 10 CONTINUE CALL PRINTJ(NAME,12) C END C C C SUBROUTINE CUTNL(FAIL) C C THIS SUBROUTINE EXAMINES THE CASE WHERE THERE ARE MORE THAN C TWO FREE ENDS,BUT THEY ARE CONTIGUOUS,SO THAT THE GRAPH CAN C BE CUT WITHOUT DESTROYING THE FLAT STRUCTURE. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C INTEGER ARROW LOGICAL TABS COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C COMMON /KEEP/JKP(2,3),JARR(2,3),IT2,IT3,IT5 C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C LOGICAL FAIL C CHARACTER*6 NAME DATA NAME/'CUTNL '/ C NTF = ITFREE(NFREE) - ITFREE(1) IF (NTF.GT.NFREE) GOTO 90 IT2 = ITFREE(1) IT3 = ITFREE(NFREE) IT1 = IT2 - 1 IT4 = IT3 + 1 IF (NTF.EQ.NFREE) GOTO 20 JT = JDIAG(IT2,3) CALL DELTA(JT,JDIAG(IT3,2),FAIL) IF (FAIL) GOTO 80 IF (ARR(IT2,3).NE.ARR(IT3,2)) GOTO 10 C C PHASE2(JT): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JT C ARR(IT2,3) = -ARR(IT2,3) ARR(IT1,2) = -ARR(IT1,2) 10 CONTINUE JDIAG(IT3,2) = JT JDIAG(IT4,3) = JT J9(J9C+1) = JT J9C = J9C + 2 J9(J9C) = JT NBTR = NBTR + NFREE IT5 = 0 GOTO 60 C 20 CONTINUE NFR = 0 DO 30 IT5 = IT2,IT3 NFR = NFR + 1 IF (ITFREE(NFR).GT.IT5) GOTO 40 30 CONTINUE 40 CONTINUE JKP(1,1) = JDIAG(IT5,1) JARR(1,1) = -ARR(IT5,1) JKP(1,2) = JDIAG(IT2,3) JARR(1,2) = -ARR(IT2,3) JKP(1,3) = JDIAG(IT3,2) JARR(1,3) = -ARR(IT3,2) DO 50 J = 1,3 JKP(2,J) = JDIAG(IT5,J) JARR(2,J) = ARR(IT5,J) 50 CONTINUE JDIAG(IT5,2) = JDIAG(IT3,2) ARR(IT5,2) = ARR(IT3,2) JDIAG(IT5,3) = JDIAG(IT2,3) ARR(IT5,3) = ARR(IT2,3) ILP = IL(IT2) IL(IT5) = ILP IH(ILP) = IT5 NBTR = NBTR + NFREE + 2 C C PHASE(IT5,JDIAG,KFL2B): C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C J7(J7C+1) = JDIAG(IT5,1) J7(J7C+2) = JDIAG(IT5,2) J7C = J7C + 3 J7(J7C) = JDIAG(IT5,3) C K = ABS((3*ARR(IT5,1)+2*ARR(IT5,2)+ARR(IT5,3))/2) + 1 IF (K.NE.4) THEN C C PHASE2(JDIAG(IT5,K)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JDIAG(IT5,K) C ENDIF C 60 CONTINUE IL1 = IL(IT4) DO 70 I = IL1,NBNODE IT = IH(I) ILP = I - NFREE IL(IT) = ILP IH(ILP) = IT 70 CONTINUE NBNODE = NBNODE - NFREE NFIN = 0 GOTO 90 C 80 CONTINUE FAIL = .TRUE. 90 CONTINUE CALL PRINTJ(NAME,8) C END C C C SUBROUTINE DELTA(JA,JB,FAIL) C C TEST FOR DELTA(JA,JB).IF THEY ARE SUMMATION VARIABLES,THE SECOND C IS CHANGED INTO THE FIRST EVERYWHERE.IF THEY ARE FIXED,THEIR C VALUE IS CHECKED,AND FAIL PUT TO .TRUE. IF THEY DIFFER. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C LOGICAL CUT COMMON /CUTDIG/CUT LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 C COMMON /DIM/J6CC,J7CC,J8CC,J9CC,JWCC,JDELC C LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C IF (IBUG3.EQ.1) PRINT 3000,JA,SUMVAR(JA),JB,SUMVAR(JB) IF (SUMVAR(JA) .AND. SUMVAR(JB)) GOTO 20 IF (FREE(JA) .OR. FREE(JB)) THEN JDEL = JDEL + 1 LDEL(JDEL,1) = JA LDEL(JDEL,2) = JB SUMVAR(JA) = .FALSE. SUMVAR(JB) = .FALSE. GOTO 160 ENDIF IF (J1(JA).NE.J1(JB)) FAIL = .TRUE. CUT = .TRUE. GOTO 160 C 20 CONTINUE IF (J6C.NE.J6CC) THEN DO 30 II = J6CC + 1,J6C IF (J6(II).EQ.JB) J6(II) = JA 30 CONTINUE ENDIF IF (J7C.NE.J7CC) THEN DO 50 II = J7CC + 1,J7C IF (J7(II).EQ.JB) J7(II) = JA 50 CONTINUE ENDIF IF (J8C.NE.J8CC) THEN DO 70 II = J8CC + 1,J8C IF (J8(II).EQ.JB) J8(II) = JA 70 CONTINUE ENDIF IF (J9C.NE.J9CC) THEN DO 90 II = J9CC + 1,J9C IF (J9(II).EQ.JB) J9(II) = JA 90 CONTINUE ENDIF IF (JWC.NE.JWCC) THEN DO 120 II = JWCC + 1,JWC IF (JW(1,II).EQ.JB) JW(1,II) = JA IF (JW(2,II).EQ.JB) JW(2,II) = JA IF (JW(3,II).EQ.JB) JW(3,II) = JA IF (JW(4,II).EQ.JB) JW(4,II) = JA IF (JW(5,II).EQ.JB) JW(5,II) = JA IF (JW(6,II).EQ.JB) JW(6,II) = JA 120 CONTINUE ENDIF IF (JDEL.NE.JDELC) THEN DO 150 II = JDELC + 1,JDEL IF (LDEL(II,1).EQ.JB) LDEL(II,1) = JA IF (LDEL(II,2).EQ.JB) LDEL(II,2) = JA 150 CONTINUE SUMVAR(JB) = .FALSE. ENDIF 160 CONTINUE C 3000 FORMAT (/' FROM DELTA JA=',I2,L2,5X,'JB=',I2,L2) END C C C SUBROUTINE DIAGRM(JUMP) C C THIS SUBROUTINE BUILDS UP A FLAT DIAGRAM FROM THE TRIADS J23 AND C PLACES THEM IN JDIAG.ARROWS ARE IN ARR (INTEGER).THE DIAGRAM IS C BUILT SO AS TO MAXIMIZE THE NUMBER OF TRIADS INVOLVED,WITHIN A C ONE-STEP-FORWARD-CHECK PROCESS.IF THE DIAGRAM DOES NOT C INCLUDE ALL THE NBTR TRIADS,IT WILL HAVE 'FREE ENDS'.JDIAG HAS C DIMENSION DOUBLE THAT OF J23,BECAUSE THE PATH MAY PROCEED EITHER C WAY. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP LOGICAL TABS INTEGER ARROW COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C COMMON /BUILD/IAL(KFL2B),IF1,IF2,NODE C CHARACTER*6 NAME C SAVE NB C DATA NAME/'DIAGRM'/ C C INITIALIZATION C IF (JUMP-2) 10,20,50 10 CONTINUE NB = 0 20 CONTINUE NB = NB + 1 IF (TABS(NB)) GOTO 20 NODE = NBTR ILAST = NBTR DO 30 J = 1,3 JDIAG(NODE,J) = J23(NB,J) ARR(NODE,J) = ARROW(NB,J) 30 CONTINUE TABS(NB) = .TRUE. DO 40 I = 1,MP IAL(I) = 0 40 CONTINUE IF1 = JDIAG(NODE,1) IF2 = JDIAG(NODE,3) IAL(IF1) = 1 IAL(IF2) = 1 50 CONTINUE NTIME = 0 I1 = 1 K1 = 1 K2 = 2 K3 = 3 60 CONTINUE JB = JDIAG(NODE,K2) C C OTHERJ(0,JB,L,LC,KP): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C L = LINE(JB,1) IF (L.EQ.0 .OR. TABS(L)) THEN L = LINE(JB,2) LC = LCOL(JB,2) C ELSE LC = LCOL(JB,1) ENDIF C C NEIBOR(LC,L1,L2): C GIVES THE POSITIONS OF THE OTHER TWO ARGUMENTS IN THE TRIAD. C IF (LC.LT.2) THEN L1 = 2 L2 = 3 C ELSE IF (LC.EQ.2) THEN L1 = 3 L2 = 1 C ELSE L1 = 1 L2 = 2 ENDIF C CALL WAY(L,L1,L2,ICH,ND) NODE = NODE + I1 TABS(L) = .TRUE. JDIAG(NODE,K3) = J23(L,LC) ARR(NODE,K3) = ARROW(L,LC) ICT = ICH*I1 IF (ICH) 70,70,80 70 CONTINUE LP = L1 L1 = L2 L2 = LP 80 CONTINUE IF (ICT.LE.0) THEN C C PHASE(L,J23,KFL2A) C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C J7(J7C+1) = J23(L,1) J7(J7C+2) = J23(L,2) J7C = J7C + 3 J7(J7C) = J23(L,3) C ENDIF C JDIAG(NODE,K1) = J23(L,L1) ARR(NODE,K1) = ARROW(L,L1) JDIAG(NODE,K2) = J23(L,L2) ARR(NODE,K2) = ARROW(L,L2) J = J23(L,L1) IAL(J) = IAL(J) + 1 J = J23(L,L2) IAL(J) = IAL(J) + 1 IF (ND.LT.1) GOTO 60 NTIME = NTIME + 1 ILAST = MAX(NODE,ILAST) IFIRST = MIN(NODE,NBTR) NBP = IAL(IF1) + IAL(IF2) IF (NBP-3) 90,90,140 90 CONTINUE IF (NTIME-1) 100,100,140 100 CONTINUE IF (NBP-2) 130,130,110 110 CONTINUE IF (IAL(IF1)-IAL(IF2)) 120,120,130 120 CONTINUE JT = JDIAG(NBTR,1) JAR = ARR(NBTR,1) JDIAG(NBTR,1) = JDIAG(NBTR,3) ARR(NBTR,1) = ARR(NBTR,3) JDIAG(NBTR,3) = JT ARR(NBTR,3) = JAR C C PHASE(NBTR,JDIAG,KFL2B): C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C J7(J7C+1) = JDIAG(NBTR,1) J7(J7C+2) = JDIAG(NBTR,2) J7C = J7C + 3 J7(J7C) = JDIAG(NBTR,3) C 130 CONTINUE NODE = NBTR I1 = -1 K2 = 3 K3 = 2 GOTO 60 C 140 CONTINUE NBNODE = ILAST - IFIRST + 1 NBTR = NBTR - NBNODE C C DEFINITION OF FREE ENDS AND OTHER QUANTITIES. C CALL INTAB CALL PRINTJ(NAME,12) C END C C C SUBROUTINE DRACAH(I,J,K,L,M,N,RAC) C*********************************************************************** C C DRACAH C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C C SUBROUTINE TO CALCULATE RACAH COEFFICIENTS C THE ARGUMENTS I,J,K,L,M,N SHOULD BE TWICE THEIR ACTUAL VALUE C WORKS FOR INTEGER AND HALF-INTEGER VALUES OF ANGULAR MOMENTA. C THE ROUTINE MAKES USE OF THE GAM ARRAY, THUS SUBROUTINE FACTT C MUST BE CALLED BEFORE THIS ROUTINE IS USED. C WRITTEN BY N S SCOTT; CHECK IF...PRINT6 ADDED WE'89JUN19TH C PARAMETER (MXFCT=500) C COMMON /FACTS/GAM(MXFCT) C DATA ZERO,ONE,HALF/0.0D0,1.0D0,0.5D0/ C RAC = ZERO J1 = I + J + M J2 = K + L + M J3 = I + K + N J4 = J + L + N IF (2*MAX(I,J,M).GT.J1 .OR. MOD(J1,2).NE.0) GOTO 30 IF (2*MAX(K,L,M).GT.J2 .OR. MOD(J2,2).NE.0) GOTO 30 IF (2*MAX(I,K,N).GT.J3 .OR. MOD(J3,2).NE.0) GOTO 30 IF (2*MAX(J,L,N).GT.J4 .OR. MOD(J4,2).NE.0) GOTO 30 J5 = (I+J+K+L)/2 + 2 J6 = (I+L+M+N)/2 + 2 J7 = (J+K+M+N)/2 + 2 NUMAX = MIN(J5,J6,J7) - 1 IF (NUMAX.GE.MXFCT) GOTO 40 RAC = ONE J1 = J1/2 J2 = J2/2 J3 = J3/2 J4 = J4/2 NUMIN = MAX(J1,J2,J3,J4) + 1 IF (NUMIN.EQ.NUMAX) GOTO 20 KF = NUMIN + 1 DO 10 KI = NUMAX,KF,-1 RAC = - (KI* (J5-KI)* (J6-KI)* (J7-KI))*RAC/ A ((KI-1-J1)* (KI-1-J2)* (KI-1-J3)* (KI-1-J4)) + ONE 10 CONTINUE 20 CONTINUE RAC = (2*MOD(J5+NUMIN,2)-1)*EXP(GAM(NUMIN+1)-GAM(NUMIN-J1)- A GAM(NUMIN-J2)-GAM(NUMIN-J3)-GAM(NUMIN-J4)-GAM(J5-NUMIN)- B GAM(J6-NUMIN)-GAM(J7-NUMIN)+ (GAM(J1+1-I)+GAM(J1+1-J)+ C GAM(J1+1-M)-GAM(J1+2)+GAM(J2+1-K)+GAM(J2+1-L)+GAM(J2+1-M)- D GAM(J2+2)+GAM(J3+1-I)+GAM(J3+1-K)+GAM(J3+1-N)-GAM(J3+2)+ E GAM(J4+1-J)+GAM(J4+1-L)+GAM(J4+1-N)-GAM(J4+2))*HALF)*RAC 30 CONTINUE RETURN C 40 CONTINUE PRINT 3000,NUMAX STOP C 3000 FORMAT (' STOP IN DRACAH, MXFCT >',I4, A ' NEEDED FOR FACTORIAL ARRAY.') END C C C SUBROUTINE EIGEN(N,EIG,EPSI,P,R,POLY,BETA) C C ACCEPTS THE ARRAYS R AND P OF MAIN AND SUPER DIAGONAL ELEMENTS C RESPECTIVELY. USING THE STURM SEQUENCE PROPERTY C A BISECTION METHOD IS APPLIED TO DETERMINE THE EIGENVALUES C (STORED IN THE ARRAY EIG ON RETURN) TO AN ACCURACY C SPECIFIED BY EPSI. N IS AS DEFINED IN THE SUBROUTINE HSLDR. C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION EIG(N),P(N),R(N),POLY(N),BETA(N) DATA ZERO,HALF,ONE/0.0D0,0.5D0,1.0D0/, A TINY/1.0D-5/ C C CALCULATE THE AVERAGE OF THE GREATEST AND SMALLEST MAIN C DIAGONAL ELEMENTS STORING THE RESULT IN AMID. C ASMALL = R(1) ALARG = R(1) DO 10 I = 2,N IF (R(I).GT.ALARG) ALARG = R(I) IF (R(I).LT.ASMALL) ASMALL = R(I) 10 CONTINUE AMID = (ALARG+ASMALL)*HALF C C REDUCE EACH MAIN DIAGONAL ELEMENT BY AMID AND CALCUATE, USING C THE STURM SEQUENCE PROPERTY, THE EIGENVALUES OF THE CORRESPONDING C REDUCED TRI-DIAGONAL MATRIX. C DO 20 I = 1,N R(I) = R(I) - AMID 20 CONTINUE C C CALCULATE THE MAXIMUM INFINITY NORM G OF THE MATRIX. THE C EIGENVALUES LIE IN THE RANGE -G TO +G. C G = ZERO DO 30 I = 1,N G1 = ABS(R(I)) IF (I.GT.1) G1 = G1 + ABS(P(I-1)) IF (I.LT.N) G1 = G1 + ABS(P(I)) IF (G.LT.G1) G = G1 30 CONTINUE C C CALCULATE THE SQUARES OF THE SUPER DIAGONAL ELEMENTS AND STORE C THESE IN THE ARRAY BETA. C N1 = N - 1 DO 40 I = 1,N1 BETA(I) = P(I)*P(I) 40 CONTINUE C C THIS LOOP DETERMINES THE EIGENVALUES ONE AT A TIME IN ORDER OF C ALGEBRAIC SIZE DOWNWARDS. C DO 100 K = 1,N AL = -G BL = G C C ONCE THROUGH THIS LOOP IS ONE BISECTION OF THE RANGE. CL1 IS THE C CURRENT ESTIMATE, CL THE IMMEDIATELY PREVIOUS ESTIMATE OF THE C EIGENVALUE. C DO 80 J = 1,100 CL1 = (AL+BL)*HALF IF (J.EQ.1) GOTO 50 C C IF THE EIGENVALUE HAS BEEN DETERMINED TO A SPECIFIED ACCURACY C EPSI, THE CALCULATION IS COMPLETE. C IF (ABS(CL1-CL).LT.EPSI) GOTO 90 C C LSUM STORES THE NUMBER OF AGREEMENTS IN SIGN IN THE STURM C SEQUENCE. C 50 CONTINUE LSUM = 0 DO 60 I = 1,N POLY(I) = R(I) - CL1 60 CONTINUE X = POLY(1) IF (POLY(1).GT.ZERO) LSUM = 1 C C THIS LOOP CALCULATES ALL THE REMAINING MEMBERS OF THE STURM C SEQUENCE. THE NUMBER OF AGREEMENTS IN SIGN IS ALSO DETERMINED. C DO 70 I = 2,N IF (X.EQ.ZERO) THEN X = POLY(I) - ABS(P(I-1))/ (ONE+TINY) C ELSE X = POLY(I) - BETA(I-1)/X ENDIF C IF (X.GT.ZERO) LSUM = LSUM + 1 70 CONTINUE CL = CL1 C C THE NEW RANGE FOR THE EIGENVALUE (DEPENDENT ON THE VALUE OF LSUM) C IS DETERMINED. C IF (LSUM.LT.K) THEN BL = CL1 C ELSE AL = CL1 ENDIF C 80 CONTINUE C C THE EIGENVALUE IS STORED IN THE ARRAY EIG. C 90 CONTINUE EIG(K) = CL1 C C RETURN TO CALCULATE THE NEXT EIGENVALUE. C 100 CONTINUE C C THE ELEMENTS OF THE ORIGINAL TRI-DIAGONAL MATRIX ARE REGAINED C AND ITS EIGENVALUES FOUND. C DO 110 I = 1,N R(I) = R(I) + AMID EIG(I) = EIG(I) + AMID 110 CONTINUE C END C C C SUBROUTINE EIGVEC(N,A,LENGTH,X,P) C C TAKES THE EIGENVECTOR OF THE TRI-DIAGONAL MATRIX C STORED IN X AND DETAILS OF THE MATRICES USED IN TRANSFORMING C THE ORIGINAL MATRIX TO TRI-DIAGONAL FORM, STORED IN A, C AND OBTAINS THE CORRESPONDING EIGENVECTOR OF THE ORIGINAL C MATRIX. N AND LENGTH ARE AS DEFINED IN THE SUBROUTINE HSLDR. C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION A(LENGTH),X(N),P(N) DATA ZERO/0.0D0/ C N2 = N - 2 N22 = (N* (N+1))/2 C C N2=N-2 TRANSFORMATIONS TO OBTAIN EACH EIGENVECTOR. C DO 30 K = 1,N2 C C K1 IS THE NUMBER OF ELEMENTS IN THE FIRST (K-1) ROWS OF THE C UPPER TRIANGLE STORED IN A. C K1 = N22 - ((K+2)* (K+3))/2 + 1 SOP = ZERO KP1 = K + 1 NK1 = N - K - 1 DO 10 I = 1,KP1 SOP = SOP + A(K1+I)*X(NK1+I) 10 CONTINUE IF (SOP.EQ.ZERO) GOTO 30 C C FROM INFORMATION STORED IN THE ARRAY A BKR IS DETERMINED AS IN C SUBROUTINE HOUSE. C BKR = - (P(NK1)*A(K1+1)) SOP = SOP/BKR DO 20 J = 1,KP1 X(NK1+J) = X(NK1+J) - A(K1+J)*SOP 20 CONTINUE 30 CONTINUE C C THE EIGENVECTOR OF THE ORIGINAL MATRIX IS NORMALISED. C CALL NORM(N,X) C END C C C SUBROUTINE FACTT C*********************************************************************** C C FACTT C C*********************************************************************** IMPLICIT REAL*8 (G,L,O,X) C C CALCULATES THE LOGS OF FACTORIALS REQUIRED BY THE RACAH C COEFFICIENT ROUTINE DRACAH C WRITTEN BY N.S. SCOTT; MXFCT CODE ADDED WE'89JUN19TH. C INCLUDE 'PARAM' PARAMETER (MXFCT=500) C COMMON /FACTS/GAM(MXFCT) C ONE = 1.D0 GAM(1) = ONE X = ONE M = MIN(MXFCT,MZFAC) DO 10 I = 2,M GAM(I) = GAM(I-1)*X X = X + ONE 10 CONTINUE DO 20 I = 1,M GAM(I) = LOG(GAM(I)) 20 CONTINUE M = M + 1 DO 30 I = M,MXFCT GAM(I) = GAM(I-1) + LOG(X) X = X + ONE 30 CONTINUE C END SUBROUTINE FANO(IRHO,ISIG,IRHOP,ISIGP) C C C*********************************************************************** C C FANO C C*********************************************************************** C C ANGULAR MOMENTUM INTEGRALS IN ATOMIC STRUCTURE BY A HIBBERT. C CATALOGUE NUMBER ACQV (CPC 2 (1971) 180) MODIFIED BY THE C CORRECTION DECK CATALOGUE NUMBER ACQV000A (CPC 6 (1973) 59) AND C THE ADAPTATION DECK CATALOGUE NUMBER ACQV0001 (CPC 7 (1974) ) C C *** DEFINITION OF DIMENSION LIST C C RMEDIR(K),K=KD1,KD2,2 - DIRECT REDUCED MATRIX ELEMENT PRODUCT C RMEEX(K),K=KE1,KE2,2 - EXCHANGE REDUCED MATRIX ELEMENT PRODUCT C KD1,KE1 ARE ALWAYS .GE. 1 C KD2,KE2 ARE .LE. 1+2*MAX(L-VALUE) C NBAR(I), I=1,IHSH - NUMBER OF SPECTATOR ELECTRONS IN EACH SHELL C THE OTHER ARRAYS ARE DEFINED IN NJGRAF. C IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXCTAB=MZLR2*MZLR2/2*MZLR2/2+MZLR2*MZLR2/2+MZLR2) PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) PARAMETER (MX2LR2=2*MZLR2) C PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20,KFLS=12, A KFLN=10,KFLV=40) C PARAMETER (MXCOUP=16,MXIJK=MXCOUP*MXCOUP) C LOGICAL FAILSD,FAILSE,FAILAD,FAILAE,FREE C DIMENSION IJKST3(MXIJK),IJKST4(MXIJK) DIMENSION RMEDIR(99),RMEEX(99),NBAR(21) DIMENSION K6SD(KFL6),K7SD(KFL7),K8SD(KFL8),K9SD(KFL9), A KWSD(6,KFLW),LDELSD(KFLW,2) DIMENSION K6SE(KFL6),K7SE(KFL7),K8SE(KFL8),K9SE(KFL9), A KWSE(6,KFLW),LDELSE(KFLW,2) DIMENSION K6AD(KFL6),K7AD(KFL7),K8AD(KFL8),K9AD(KFL9), A KWAD(6,KFLW),LDELAD(KFLW,2) DIMENSION K6AE(KFL6),K7AE(KFL7),K8AE(KFL8),K9AE(KFL9), A KWAE(6,KFLW),LDELAE(KFLW,2) DIMENSION J6PSD(KFLV),J7PSD(KFLV),J8PSD(KFLV),J9PSD(KFLV), A JWRDSD(6,KFLW),NBJSD(KFLN),NB6JSD(KFLN),K6CPSD(KFLN), B K7CPSD(KFLN),K8CPSD(KFLN),K9CPSD(KFLN),JSM6SD(KFLS), C JSM4SD(KFLS,KFLW),JSM5SD(KFLS,KFLW),IN6JSD(KFLW) DIMENSION J6PSE(KFLV),J7PSE(KFLV),J8PSE(KFLV),J9PSE(KFLV), A JWRDSE(6,KFLW),NBJSE(KFLN),NB6JSE(KFLN),K6CPSE(KFLN), B K7CPSE(KFLN),K8CPSE(KFLN),K9CPSE(KFLN),JSM6SE(KFLS), C JSM4SE(KFLS,KFLW),JSM5SE(KFLS,KFLW),IN6JSE(KFLW) DIMENSION J6PAD(KFLV),J7PAD(KFLV),J8PAD(KFLV),J9PAD(KFLV), A JWRDAD(6,KFLW),NBJAD(KFLN),NB6JAD(KFLN),K6CPAD(KFLN), B K7CPAD(KFLN),K8CPAD(KFLN),K9CPAD(KFLN),JSM6AD(KFLS), C JSM4AD(KFLS,KFLW),JSM5AD(KFLS,KFLW),IN6JAD(KFLW) DIMENSION J6PAE(KFLV),J7PAE(KFLV),J8PAE(KFLV),J9PAE(KFLV), A JWRDAE(6,KFLW),NBJAE(KFLN),NB6JAE(KFLN),K6CPAE(KFLN), B K7CPAE(KFLN),K8CPAE(KFLN),K9CPAE(KFLN),JSM6AE(KFLS), C JSM4AE(KFLS,KFLW),JSM5AE(KFLS,KFLW),IN6JAE(KFLW) C COMMON /CONSTS/ZERO,TENTH,HALF,ONE,TWO,THREE,FOUR,SEVEN,ELEVEN,EPS COMMON /CSTORE/CTABLE(MXCTAB),KPOINT(MZLR2,MZLR2),LRANG3 COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /KRON/IDEL(MXORB2,MXORB2) COMMON /TERMS/NROWS,ITAB(18),JTAB(18),NTAB(189) COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 COMMON /XATION/AMULT(MX2LR2),BMULT(MX2LR2),KD1,KD2,KE1,KE2,MULTD, A MULTE,LNOEX COMMON /NJLJ/NRHO,LRHO,NSIG,LSIG,NRHOP,LRHOP,NSIGP,LSIGP COMMON /INTERM/J1BAR1(MXORB2,3),J1BAR2(MXORB2,3),J1TLD1(3), A J1TLD2(3) COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) C C === DETERMINES POTENTIAL MATRIX ELEMENT FOR GIVEN RHO,SIG,RHOP,SIGP C C --- SET USEFUL CONSTANTS C M1 = ISIG - IRHO M2 = ISIGP - IRHOP M19 = IRHO - IRHOP M20 = ISIG - ISIGP MULTD = 0 MULTE = 0 C C JSNDIR,JANGDI=0 IMPLIES APPROPRIATE J2,J3 ARRAYS FOR CALL OF C NJGRAF HAVE NOT BEEN SET C JSNDIR = 0 JANGDI = 0 ISPIND = 0 ISPINE = 0 IANGD = 0 IANGE = 0 C C SET THE FAIL PARAMETERS .FALSE. INITIALLY SO THAT NJGRAF CAN C BE CALLED. C FAILSD = .FALSE. FAILSE = .FALSE. FAILAD = .FALSE. FAILAE = .FALSE. C C --- SET N,L VALUES OF INTERACTING SHELLS C NRHO = NJ(IRHO) LRHO = LJ(IRHO) NSIG = NJ(ISIG) LSIG = LJ(ISIG) NRHOP = NJ(IRHOP) LRHOP = LJ(IRHOP) NSIGP = NJ(ISIGP) LSIGP = LJ(ISIGP) IF (IBUG1.GE.2) WRITE (IWRITE,3040) NRHO,LRHO,NSIG,LSIG,NRHOP, A LRHOP,NSIGP,LSIGP C C --- EVALUATE MULTIPLICATIVE FACTORS C IL = IDEL(IRHO,ISIG) IR = IDEL(IRHOP,ISIGP) Q = NOSH1(IRHO)* (NOSH1(ISIG)-IL)*NOSH2(IRHOP)* (NOSH2(ISIGP)-IR) XMULT = SQRT(Q)*HALF ADIRCT = (1+ (1-IL)* (1-IR))/SQRT((LSIG+LSIG+ONE)* A (LRHOP+LRHOP+ONE)) IEXCHG = 2 - IL - IR AEXCHG = IEXCHG/SQRT((LSIG+LSIG+ONE)* (LSIGP+LSIGP+ONE)) DO 10 J = 1,IHSH NBAR(J) = NOSH1(J) - IDEL(J,IRHO) - IDEL(J,ISIG) 10 CONTINUE IDELP = 0 IF (M1.NE.0) THEN K1 = IRHO + 1 DO 20 J = K1,ISIG IDELP = IDELP + NBAR(J) 20 CONTINUE ENDIF C IF (M2.NE.0) THEN K1 = IRHOP + 1 DO 30 J = K1,ISIGP IDELP = IDELP + NBAR(J) 30 CONTINUE ENDIF C XMULT = XMULT* (-ONE)**IDELP C C --- DETERMINES RANGES OF K FOR DIRECT AND EXCHANGE INTEGRALS C TRIANGULAR RELATIONS LIMIT (K+1) VALUES TO LIE BETWEEN KD1 AND KD2 C FOR =DIRECT= INTEGRALS, KE1 AND KE2 FOR =EXCHANGE= INTEGRALS C K1 = ABS(LSIG-LSIGP) K2 = LSIG + LSIGP K3 = ABS(LRHO-LRHOP) K4 = LRHO + LRHOP KD1 = MAX(K1,K3) + 1 KD2 = MIN(K2,K4) + 1 K1 = ABS(LRHOP-LSIG) K2 = LRHOP + LSIG K3 = ABS(LRHO-LSIGP) K4 = LRHO + LSIGP KE1 = MAX(K1,K3) + 1 KE2 = MIN(K2,K4) + 1 CDW COLD IF(KE1.GT.LNOEX)KE1=KE2+1 if(nrho.eq.999.or.nrhop.eq.999)stop 'cont' IF(NSIG.EQ.999.AND.NSIGP.EQ.999)KE2=MIN(KE2,LNOEX+1) CDW IF (IBUG1.GT.1) THEN WRITE (IWRITE,3050) KD1,KD2,KE1,KE2 IF (KD1.GT.KD2 .AND. KE1.GT.KE2) WRITE (IWRITE,3000) ENDIF C IF (KD1.GT.KD2 .AND. KE1.GT.KE2) GOTO 390 C C --- ZEROIZE MULTIPLYING FACTORS FOR ALLOWED K-VALUES. THE LOWEST C VALUES KD1 AND KD2 ARE ALWAYS ALLOWED (UNLESS THEY ARE C GREATER THAN KD2 AND KE2 RESPECTIVELY). ALLOWED K-VALUES THEN C STEP BY 2 TO SATISFY THE EVEN CONDITION OF THE REDUCED MATRIX C ELEMENTS, WHICH ARE THEN CALCULATED AND STORED C IF (KD1.GT.KD2) GOTO 60 IF (MOD(LRHO+LRHOP+KD1,2).NE.1 .OR. A MOD(LSIG+LSIGP+KD1,2).NE.1) THEN DO 40 JK1 = KD1,KD2,2 AMULT(JK1) = 0.0D0 RMEDIR(JK1) = 0.0D0 40 CONTINUE C ELSE LPOINT = KPOINT(LRHO+1,LRHOP+1) + (KD1-1-ABS(LRHO-LRHOP))/2 MPOINT = KPOINT(LSIG+1,LSIGP+1) + (KD1-1-ABS(LSIG-LSIGP))/2 DO 50 JK1 = KD1,KD2,2 LPOINT = LPOINT + 1 MPOINT = MPOINT + 1 AMULT(JK1) = 0.0D0 RMEDIR(JK1) = CTABLE(LPOINT)*CTABLE(MPOINT) 50 CONTINUE ENDIF C 60 CONTINUE IF (KE1.GT.KE2) GOTO 90 IF (MOD(LRHO+LSIGP+KE1,2).NE.1 .OR. A MOD(LSIG+LRHOP+KE1,2).NE.1) THEN DO 70 JK1 = KE1,KE2,2 BMULT(JK1) = 0.0D0 RMEEX(JK1) = 0.0D0 70 CONTINUE C ELSE LPOINT = KPOINT(LRHO+1,LSIGP+1) + (KE1-1-ABS(LRHO-LSIGP))/2 MPOINT = KPOINT(LSIG+1,LRHOP+1) + (KE1-1-ABS(LSIG-LRHOP))/2 DO 80 JK1 = KE1,KE2,2 LPOINT = LPOINT + 1 MPOINT = MPOINT + 1 BMULT(JK1) = 0.0D0 RMEEX(JK1) = CTABLE(LPOINT)*CTABLE(MPOINT) 80 CONTINUE ENDIF C C --- SET SENIORITY, S AND L VALUES OF SPECTATOR SHELLS C 90 CONTINUE DO 130 J = 1,IHSH IF (IRHO.NE.J .AND. ISIG.NE.J) THEN DO 100 K = 1,3 J1BAR1(J,K) = J1QN1(J,K) 100 CONTINUE ENDIF C IF (IRHOP.EQ.J .OR. ISIGP.EQ.J) GOTO 130 DO 110 K = 1,3 J1BAR2(J,K) = J1QN2(J,K) 110 CONTINUE IF (IRHO.EQ.J .OR. ISIG.EQ.J) GOTO 130 C C CHECK THAT SPECTATOR SHELLS HAVE SAME RESPECTIVE QUANTUM NUMBERS C DO 120 K = 1,3 IF (J1BAR1(J,K).EQ.J1BAR2(J,K)) GOTO 120 IF (IBUG1.GE.2) WRITE (IWRITE,3010) RETURN C 120 CONTINUE 130 CONTINUE C C --- FROM WHICH ROWS OF NTAB DO WE FIND THE QUANTUM NUMBERS WITH BARS C OR TILDES C NELCTS = NOSH1(ISIG) K2 = NTAB1(NELCTS,LSIG+1) IF (M1.EQ.0) THEN NELCTS = NOSH1(ISIG) - 1 C ELSE NELCTS = NOSH1(IRHO) ENDIF C K1 = NTAB1(NELCTS,LRHO+1) NELCTS = NOSH2(ISIGP) K4 = NTAB1(NELCTS,LSIGP+1) IF (M2.EQ.0) THEN NELCTS = NOSH2(ISIGP) - 1 C ELSE NELCTS = NOSH2(IRHOP) ENDIF C K3 = NTAB1(NELCTS,LRHOP+1) IF (IBUG1.GE.2) WRITE (IWRITE,3060) K1,K2,K3,K4 KK1 = ITAB(K1) KK2 = ITAB(K2) C C STORE IJK3 AND IJK4 (NTAB POSITIONS) FOR INNER LOOPS (JJ3 AND JJ4) C KK3 = 0 DO 170 JJ4 = 1,ITAB(K4) IN3 = 2*LSIGP IJK4 = 3* (JJ4-1) + JTAB(K4) DO 140 K = 2,3 IN1 = NTAB(IJK4+K) - 1 IN2 = J1QN2(ISIGP,K) - 1 IF (IN1.GT.IN2+IN3 .OR. IN1.LT.ABS(IN2-IN3)) GOTO 170 IN3 = 1 140 CONTINUE C DO 160 JJ3 = 1,ITAB(K3) IN3 = 2*LRHOP IJK3 = 3* (JJ3-1) + JTAB(K3) DO 150 K = 2,3 IN1 = NTAB(IJK3+K) - 1 IF (M2.EQ.0) THEN IN2 = NTAB(IJK4+K) - 1 C ELSE IN2 = J1QN2(IRHOP,K) - 1 ENDIF C IF (IN1.GT.IN2+IN3 .OR. IN1.LT.ABS(IN2-IN3)) GOTO 160 IN3 = 1 150 CONTINUE KK3 = KK3 + 1 IJKST3(KK3) = IJK3 IJKST4(KK3) = IJK4 160 CONTINUE 170 CONTINUE IF (KK3.EQ.0) GOTO 360 C C === SUM OVER QUANTUM NUMBERS WITH BARS OR TILDES C C --- TEST TO SEE WHICH PARENT TERMS ARE ALLOWABLE. ONLY TEST THIS ON C L AND S VALUES AT THIS STAGE, BY MEANS OF TRIANGULAR CONDITIONS C FOR TWICE THE QUANTUM NUMBERS, IN ORDER TO USE ONLY INTEGER C QUANTITIES C DO 350 JJ2 = 1,KK2 IN3 = 2*LSIG IJK2 = 3* (JJ2-1) + JTAB(K2) DO 180 K = 2,3 IN1 = NTAB(IJK2+K) - 1 IN2 = J1QN1(ISIG,K) - 1 IF (IN1.GT.IN2+IN3 .OR. IN1.LT.ABS(IN2-IN3)) GOTO 350 IN3 = 1 180 CONTINUE C DO 340 JJ1 = 1,KK1 IN3 = 2*LRHO IJK1 = 3* (JJ1-1) + JTAB(K1) DO 190 K = 2,3 IN1 = NTAB(IJK1+K) - 1 IF (M1.NE.0) THEN IN2 = J1QN1(IRHO,K) - 1 C ELSE IN2 = NTAB(IJK2+K) - 1 ENDIF C IF (IN1.GT.IN2+IN3 .OR. IN1.LT.ABS(IN2-IN3)) GOTO 340 IN3 = 1 190 CONTINUE C DO 320 JJ3 = 1,KK3 IJK3 = IJKST3(JJ3) IJK4 = IJKST4(JJ3) C C SUMMATIONS NOW PERFORMED OVER ALLOWED QUANTUM NUMBERS C THE TILDES CORRESPOND TO IRHO=ISIG AND/OR IRHOP=ISIGP C C --- SET THE REMAINING QUANTUM NUMBERS WITH BARS OR TILDES C DO 200 K = 1,3 J1BAR1(IRHO,K) = NTAB(IJK1+K) J1BAR2(IRHOP,K) = NTAB(IJK3+K) IF (M1.NE.0) THEN J1BAR1(ISIG,K) = NTAB(IJK2+K) C ELSE J1TLD1(K) = NTAB(IJK2+K) ENDIF C IF (M2.NE.0) THEN J1BAR2(ISIGP,K) = NTAB(IJK4+K) C ELSE J1TLD2(K) = NTAB(IJK4+K) ENDIF C 200 CONTINUE C C --- IS POTENTIAL DIAG. IN BARRED QU. NOS. FOR INTERACTING SHELLS C I = ISIG I5 = 0 210 CONTINUE I5 = I5 + 1 IF (J1BAR1(I,1).NE.J1BAR2(I,1) .OR. A J1BAR1(I,2).NE.J1BAR2(I,2) .OR. B J1BAR1(I,3).NE.J1BAR2(I,3)) GOTO 320 IF (I5.EQ.1) THEN I = IRHO IF (M1.EQ.0) I5 = 2 C ELSE IF (I5.EQ.3) THEN I = IRHOP IF (M2.EQ.0) GOTO 220 ENDIF C IF (I5.EQ.2) I = ISIGP IF (I5.LE.3) GOTO 210 220 CONTINUE PICFP = 1.0D0 C C --- EVALUATE FRACTIONAL PARENTAGE COEFFICIENTS C I = 1 CALL MUMDAD(I,ISIG,IRHO,M1,CFPLHS) PICFP = PICFP*CFPLHS IF (ABS(PICFP).LT.EPS) GOTO 320 I = 2 CALL MUMDAD(I,ISIGP,IRHOP,M2,CFPRHS) PICFP = PICFP*CFPRHS IF (ABS(PICFP).LT.EPS) GOTO 320 C C === SET UP J1,J2,J3 AND EVALUATE RECOUPLING COEFFICIENTS C C --- FIRST OF ALL, THE SPIN COEFFICIENTS C I = 3 CALL SETJ1(I,IRHO,ISIG,IRHOP,ISIGP,ISPIND,ISPINE,KK1,KK2, A ITAB(K3),ITAB(K4)) IF (ISPIND.EQ.0) THEN CALL J23SPN(IRHO,ISIG,IRHOP,ISIGP,JSNDIR) ISPIND = 1 ENDIF C C --- DIRECT SPIN INTEGRAL C IF (KD1.GT.KD2) THEN SPINDT = 0.0D0 GOTO 230 C ENDIF C IF (.NOT.FAILSD) THEN IF (ISPIND.NE.2) THEN CALL NJGRAF(SPINDT,FAILSD) ISPIND = 2 IF (FAILSD) GOTO 230 CALL KNJ(J6CSD,J7CSD,J8CSD,J9CSD,JWCSD,K6SD,K7SD,K8SD, A K9SD,KWSD,JDELSD,LDELSD,MPSD,J6PSD, B J7PSD,J8PSD,J9PSD,JWRDSD,NLSMSD,NBJSD,NB6JSD, C K6CPSD,K7CPSD,K8CPSD,K9CPSD,JSM4SD,JSM5SD, D JSM6SD,IN6JSD) C ELSE CALL GENSUM(J6CSD,J7CSD,J8CSD,J9CSD,JWCSD,K6SD,K7SD, A K8SD,K9SD,KWSD,JDELSD,LDELSD,MPSD, B J6PSD,J7PSD,J8PSD,J9PSD,JWRDSD,NLSMSD,NBJSD, C NB6JSD,K6CPSD,K7CPSD,K8CPSD,K9CPSD,JSM4SD, D JSM5SD,JSM6SD,IN6JSD,SPINDT) C ENDIF C ELSE SPINDT = 0.0D0 ENDIF C 230 CONTINUE IF (IBUG1.GE.2) WRITE (IWRITE,3020) SPINDT C C IEXCHG IS ZERO WHENEVER M1=0=M2 , IN WHICH CASE THE EXCHANGE C INTEGRAL HAS ZERO COEFFICIENT. THERE IS THEN NO POINT IN C CALCULATING THIS INTEGRAL, AND SPINEX IS SET ZERO (AT STATEMENT C 93) AS A MARKER OF THIS SITUATION C IF (IEXCHG.EQ.0) GOTO 250 C C --- MODIFY J2 AND J3 TO CALCULATE THE EXCHANGE SPIN INTEGRAL C IF (ISPINE.EQ.0) THEN I = 1 CALL MODJ23(I) ISPINE = 1 ENDIF C C --- EXCHANGE SPIN INTEGRAL C IF (KE1.LE.KE2) GOTO 260 250 CONTINUE SPINEX = 0.0D0 GOTO 270 C 260 CONTINUE IF (.NOT.FAILSE) THEN IF (ISPINE.NE.2) THEN CALL NJGRAF(SPINEX,FAILSE) ISPINE = 2 IF (FAILSE) GOTO 270 CALL KNJ(J6CSE,J7CSE,J8CSE,J9CSE,JWCSE,K6SE,K7SE,K8SE, A K9SE,KWSE,JDELSE,LDELSE,MPSE,J6PSE, B J7PSE,J8PSE,J9PSE,JWRDSE,NLSMSE,NBJSE,NB6JSE, C K6CPSE,K7CPSE,K8CPSE,K9CPSE,JSM4SE,JSM5SE, D JSM6SE,IN6JSE) C ELSE CALL GENSUM(J6CSE,J7CSE,J8CSE,J9CSE,JWCSE,K6SE,K7SE, A K8SE,K9SE,KWSE,JDELSE,LDELSE,MPSE, B J6PSE,J7PSE,J8PSE,J9PSE,JWRDSE,NLSMSE,NBJSE, C NB6JSE,K6CPSE,K7CPSE,K8CPSE,K9CPSE,JSM4SE, D JSM5SE,JSM6SE,IN6JSE,SPINEX) ENDIF C ELSE SPINEX = 0.0D0 ENDIF C 270 CONTINUE IF (IBUG1.GE.2) WRITE (IWRITE,3030) SPINEX C C --- MULTIPLY SPIN INTEGRALS BY PRODUCT OF FRACTIONAL PARENTAGE COEFFS C BDIRCT = SPINDT*PICFP BEXCHG = SPINEX*PICFP C C --- THE ANGULAR RECOUPLING COEFFICIENTS C SET J1,J2,J3 (COMPARE SPIN INTEGRALS) C C IF BOTH SPIN INTEGRALS ARE ZERO, THERE IS NO PURPOSE IN C CALCULATING THE ANGULAR INTEGRALS C IF (ABS(SPINDT).LT.EPS .AND. ABS(SPINEX).LT.EPS) GOTO 310 C I = 2 CALL SETJ1(I,IRHO,ISIG,IRHOP,ISIGP,IANGD,IANGE,KK1,KK2, A ITAB(K3),ITAB(K4)) IF (IANGD.EQ.0) CALL J23ANG(IRHO,ISIG,IRHOP,ISIGP,JANGDI) C C IF THE DIRECT SPIN RECOUPLING COEFFICIENT IS ZERO, WE NEED NOT C CALCULATE THE CORRESPONDING ORBITAL RECOUPLING COEFFICIENT C IF (ABS(SPINDT).LT.EPS) GOTO 290 IF (IANGD.EQ.0) IANGD = 1 C C --- DIRECT ANGULAR INTEGRAL C C CONSIDER ALL ALLOWED K-VALUES C IF (.NOT.FAILAD) THEN DO 280 JK1 = KD1,KD2,2 J1(NJ1S) = 2*JK1 - 1 IF (IANGD.NE.2) THEN IF (KD2.NE.KD1) THEN FREE(NJ1S) = .TRUE. C ELSE FREE(NJ1S) = .FALSE. ENDIF C CALL NJGRAF(ANGDIR,FAILAD) IANGD = 2 IF (FAILAD) GOTO 290 CALL KNJ(J6CAD,J7CAD,J8CAD,J9CAD,JWCAD,K6AD,K7AD,K8AD, A K9AD,KWAD,JDELAD,LDELAD,MPAD,J6PAD, B J7PAD,J8PAD,J9PAD,JWRDAD,NLSMAD,NBJAD,NB6JAD, C K6CPAD,K7CPAD,K8CPAD,K9CPAD,JSM4AD,JSM5AD, D JSM6AD,IN6JAD) C ELSE CALL GENSUM(J6CAD,J7CAD,J8CAD,J9CAD,JWCAD,K6AD,K7AD, A K8AD,K9AD,KWAD,JDELAD,LDELAD,MPAD, B J6PAD,J7PAD,J8PAD,J9PAD,JWRDAD,NLSMAD, C NBJAD,NB6JAD,K6CPAD,K7CPAD,K8CPAD,K9CPAD, D JSM4AD,JSM5AD,JSM6AD,IN6JAD,ANGDIR) ENDIF C C ADD INTO THE COEFFICIENT OF THE SLATER INTEGRAL C AMULT(JK1) = AMULT(JK1) + ANGDIR*BDIRCT C C MULTD=1 WHEN A DIRECT INTEGRAL COEFFICIENT HAS BEEN CALCULATED - C FOR USE, SEE PRNTWT C MULTD = 1 IF (IBUG1.GE.2) WRITE (IWRITE,3070) ANGDIR 280 CONTINUE ENDIF C C IF THE EXCHANGE SPIN RECOUPLING COEFFICIENT IS ZERO, WE NEED NOT C CALCULATE THE CORRESPONDING ORBITAL RECOUPLING COEFFICIENT C 290 CONTINUE IF (ABS(SPINEX).LT.EPS) GOTO 310 C C --- EXCHANGE ANGULAR INTEGRAL C C --- MODIFY J2 AND J3 ARRAYS TO CALCULATE EXCHANGE TERMS C IF (IANGE.EQ.0) THEN I = 2 CALL MODJ23(I) IANGE = 1 ENDIF C C CONSIDER ALL ALLOWED K-VALUES C IF (.NOT.FAILAE) THEN DO 300 JK1 = KE1,KE2,2 J1(NJ1S) = 2*JK1 - 1 IF (IANGE.NE.2) THEN IF (KE2.NE.KE1) THEN FREE(NJ1S) = .TRUE. C ELSE FREE(NJ1S) = .FALSE. ENDIF C CALL NJGRAF(ANGEX,FAILAE) IANGE = 2 IF (FAILAE) GOTO 310 CALL KNJ(J6CAE,J7CAE,J8CAE,J9CAE,JWCAE,K6AE,K7AE,K8AE, A K9AE,KWAE,JDELAE,LDELAE,MPAE,J6PAE, B J7PAE,J8PAE,J9PAE,JWRDAE,NLSMAE,NBJAE,NB6JAE, C K6CPAE,K7CPAE,K8CPAE,K9CPAE,JSM4AE,JSM5AE, D JSM6AE,IN6JAE) C ELSE CALL GENSUM(J6CAE,J7CAE,J8CAE,J9CAE,JWCAE,K6AE,K7AE, A K8AE,K9AE,KWAE,JDELAE,LDELAE,MPAE, B J6PAE,J7PAE,J8PAE,J9PAE,JWRDAE,NLSMAE, C NBJAE,NB6JAE,K6CPAE,K7CPAE,K8CPAE,K9CPAE, D JSM4AE,JSM5AE,JSM6AE,IN6JAE,ANGEX) ENDIF C BMULT(JK1) = BMULT(JK1) - ANGEX*BEXCHG C C MULTE=1 WHEN AN EXCHANGE INTEGRAL COEFFICIENT HAS BEEN CALCULATED C MULTE = 1 IF (IBUG1.GE.2) WRITE (IWRITE,3090) ANGEX 300 CONTINUE ENDIF C C IF FAILSD OR FAILAD IS TRUE, AND FAILSE OR FAILAE IS TRUE, C THEN BOTH THE DIRECT AND EXCHANGE COEFFICIENTS WILL BE ZERO. C 310 CONTINUE IF ((FAILSD.OR.FAILAD) .AND. (FAILSE.OR.FAILAE)) GOTO 390 C 320 CONTINUE 340 CONTINUE 350 CONTINUE C C === INCLUDE MULTIPLICATIVE FACTORS COMMON TO ALL TERMS WITHIN THIS C FOUR-FOLD SUMMATION C 360 CONTINUE IF (MULTD.NE.0) THEN DO 370 JK1 = KD1,KD2,2 AMULT(JK1) = AMULT(JK1)*XMULT*RMEDIR(JK1)*ADIRCT 370 CONTINUE ENDIF C IF (MULTE.NE.0) THEN DO 380 JK1 = KE1,KE2,2 BMULT(JK1) = BMULT(JK1)*XMULT*RMEEX(JK1)*AEXCHG 380 CONTINUE ENDIF C C --- PRINT OUT THE VALUES OF THE COEFFICIENTS OF THE SLATER INTEGRALS C (THE SUBROUTINE PRNTWT IS CALLED FROM RKWTS) C 390 CONTINUE C 3000 FORMAT (' NO POSSIBLE K-VALUES') 3010 FORMAT ( A' SPECTATOR QUANTUM NUMBERS NOT DIAGONAL FOR NON-INTERACTING SHELL BS') 3020 FORMAT (' DIRECT SPIN INTEGRAL =',F10.6) 3030 FORMAT (' EXCHANGE SPIN INTEGRAL =',F10.6) 3040 FORMAT (' NJ,LJ',4 (I8,I4)) 3050 FORMAT (' KD1 =',I4,' KD2 =',I4,' KE1 =',I4,' KE2 =',I4) 3060 FORMAT (' ROWS OF TERM TABLE CONTAINING PARENTS ARE, RESPECTIVELY' A ,2 (I6,I3)) 3070 FORMAT (' DIRECT ANGULAR INTEGRAL =',F10.6) 3090 FORMAT (' EXCHANGE ANGULAR INTEGRAL =',F10.6) END C C C SUBROUTINE GENSUM(J6C,J7C,J8C,J9C,JWC,J6,J7,J8,J9,JW,JDEL,LDEL, A MP,J6P,J7P,J8P,J9P,JWORD,NLSUM,NBJ,NB6J, B K6CP,K7CP,K8CP,K9CP,JSUM4,JSUM5,JSUM6,INV6J, C RECUP) C*********************************************************************** C C GENSUM C C*********************************************************************** C C CARRIES OUT THE SUMMATION OVER COEFFICIENTS DEFINED BY THE ARRAYS C J6,J7,J8,LDEL AND JW TO GIVE RECUP C THE ENTRY IS EITHER MADE FROM NJGRAF OR DIRECTLY ASSUMING THAT THE C ARRAYS J6,...,JW HAVE ALREADY BEEN DETERMINED BY A PREVIOUS C ENTRY TO NJGRAF AND THAT THE SUMMATION IS REQUIRED FOR ANOTHER SET C OF J VALUES DEFINED BY THE ARRAY J1. TIDIED UP MEUDON'89AOUT - WE. C C RECUP IS THE RECOUPLING COEFFICIENT C C SUBROUTINE CALLED: DRACAH C IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) C PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20,KFLS=12, A KFLN=10,KFLV=40) C DIMENSION J6(KFL6),J7(KFL7),J8(KFL8),J9(KFL9),JW(6,KFLW), A LDEL(KFLW,2) C DIMENSION J6P(KFLV),J7P(KFLV),J8P(KFLV),J9P(KFLV),JWORD(6,KFLW), A NBJ(KFLN),NB6J(KFLN),K6CP(KFLN),K7CP(KFLN),K8CP(KFLN), B K9CP(KFLN),JSUM6(KFLS),JSUM4(KFLS,KFLW), C JSUM5(KFLS,KFLW),INV6J(KFLW) C LOGICAL LDIAG,NOEL C DIMENSION MAT(KFLS,KFLS),JMNP(5),JMXP(5),NOEL(KFLS),MAXLP(KFLS), A IST(6),JSUM2(KFLS),JSUM3(KFLS),JSUM(2,KFLW), B JWTEST(KFLW),WSTOR(KFLW),IPAIR(2,2),LDIAG(KFLS) DIMENSION J12(4,KFLS,KFLS) DIMENSION XJ1(KFL1) C LOGICAL FREE C COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 C DATA ZERO,ONE/0.0D0,1.0D0/, A EPSIL/1.D-10/,MXCSVR/4/ C C 1. C EVALUATES ALL TERMS IN J6,J7,J8,J9,LDEL,AND JW WHICH DO NOT INVOL C A SUMMATION.THE RESULT IS STORED IN RECUP AND IASTOR C IF (IBUG3.NE.1) GOTO 20 DO 10 I = 1,M XJ1(I) = (J1(I)-1)*0.5D0 10 CONTINUE PRINT 3080, (XJ1(I),I=1,M) PRINT 3040,NLSUM PRINT 3090 20 CONTINUE MM = M + 1 J1(MM) = 1 C C TEST DELTA FUNCTIONS C J1(MM) = 1 IF (JDEL.LE.0) GOTO 50 DO 40 I = 1,JDEL I1 = LDEL(I,1) I2 = LDEL(I,2) IF (I1.LE.MM .AND. I2.LE.MM) GOTO 30 IF (I1.GT.MM) J1(I1) = J1(I2) IF (I2.GT.MM) J1(I2) = J1(I1) GOTO 40 C 30 CONTINUE IF (J1(I1).EQ.J1(I2)) GOTO 40 RECUP = ZERO GOTO 890 C 40 CONTINUE 50 CONTINUE RECUP = ONE IF (JWC.EQ.0) GOTO 90 C C MULTIPLY RECUP BY ALL RACAH COEFFICIENTS WHICH DO NOT INVOLVE A C SUMMATION C IF (IBUG3.EQ.1) PRINT 3070 DO 80 I = 1,JWC IF (INV6J(I).GT.0) GOTO 80 DO 70 J = 1,6 I1 = JW(J,I) IST(J) = J1(I1) - 1 70 CONTINUE CALL DRACAH(IST(1),IST(2),IST(3),IST(4),IST(5),IST(6),X1) IF (IBUG3.EQ.1) PRINT 3030, (XJ1(K),K=1,6),X1 RECUP = RECUP*X1 80 CONTINUE 90 CONTINUE SQR = 1.0D0 IF (J6C.EQ.0) GOTO 120 C DO 110 I = 1,J6C I1 = J6(I) SQR = SQR*J1(I1) 110 CONTINUE 120 CONTINUE SPR = 1.0D0 IF (J9C.EQ.0) GOTO 140 DO 130 I = 1,J9C I1 = J9(I) SPR = SPR*J1(I1) 130 CONTINUE 140 CONTINUE RECUP = RECUP*SQRT(SQR/SPR) IF (ABS(RECUP).LT.EPSIL) GOTO 900 IASTOR = 0 IF (J7C.EQ.0) GOTO 170 C DO 160 I = 1,J7C I1 = J7(I) IASTOR = IASTOR + J1(I1) - 1 160 CONTINUE 170 CONTINUE IF (J8C.EQ.0) GOTO 200 C DO 190 I = 1,J8C I1 = J8(I) IASTOR = IASTOR + 2* (J1(I1)-1) 190 CONTINUE 200 CONTINUE IF (NLSUM.GT.0) GOTO 210 C IASTOR = IASTOR/2 C NO SUMMATION INVOLVED. END OF COMPUTATION C STOR1 = ONE STOR = ONE IF (MOD(IASTOR,2).EQ.1) RECUP = -RECUP IF (IBUG3.EQ.1) PRINT 3010,RECUP GOTO 890 C C C C 2. C EVALUATION OF THE PART INVOLVING SUMMATIONS. C C 210 CONTINUE NFS = 0 JWR = 0 J6F = 0 J7F = 0 J8F = 0 J9F = 0 NPS = 0 220 CONTINUE NPS = NPS + 1 IF (IBUG3.EQ.1) PRINT 3000,NPS C C 2.0 LOOP ON THE DISCONNECTED SUMMATIONS C C IAS = 0 NSUM = NBJ(NPS) - NFS JWRD = NB6J(NPS) - JWR J6CP = K6CP(NPS) J7CP = K7CP(NPS) J8CP = K8CP(NPS) J9CP = K9CP(NPS) C C 2.1 THE RANGE OF VALUES OF EACH SUMMATION VARIABLE IS C DEFINED BY ESTABLISHING A MATRIX OF THE LINKS BETWEEN C VARIABLES.MAT(I,J) CONTAINS: C I=J ,NUMBER OF POSSIBLE VALUES OF I DUE TO TRIANGULAR C RELATIONS WITH NON-VARIABLES,I.E. CONSTANTS. C I.GT.J NUMBER OF LINKS BETWEEN I AND J THROUGH CONSTANTS C I.LT.J VALUE OF THE CONSTANT,IF THE ABOVE IS 1.IF NOT, C THESE VALUES ARE SRORED IN J12(L,I,J) WHERE THERE C IS ROOM FOR MXCSVR SUCH VALUES (L.LE.4) C C DO 240 I = 1,NSUM DO 230 J = 1,NSUM MAT(I,J) = 0 230 CONTINUE 240 CONTINUE DO 420 I1 = 1,NSUM I1T = I1 + NFS I2 = JSUM6(I1T) DO 410 I3 = 1,I2 I = JSUM5(I1T,I3) J = JSUM4(I1T,I3) GOTO (250,260,270,280,290,300),J C C THE ROWS OF THE IPAIR ARRAYS GIVE LIMITS OF SUMMATION IMPOSED C BY THE TRIANGULAR CONDITION C 250 CONTINUE IPAIR(1,1) = JWORD(2,I) IPAIR(1,2) = JWORD(5,I) IPAIR(2,1) = JWORD(3,I) IPAIR(2,2) = JWORD(6,I) GOTO 310 C 260 CONTINUE IPAIR(1,1) = JWORD(1,I) IPAIR(1,2) = JWORD(5,I) IPAIR(2,1) = JWORD(4,I) IPAIR(2,2) = JWORD(6,I) GOTO 310 C 270 CONTINUE IPAIR(1,1) = JWORD(1,I) IPAIR(1,2) = JWORD(6,I) IPAIR(2,1) = JWORD(4,I) IPAIR(2,2) = JWORD(5,I) GOTO 310 C 280 CONTINUE IPAIR(1,1) = JWORD(2,I) IPAIR(1,2) = JWORD(6,I) IPAIR(2,1) = JWORD(3,I) IPAIR(2,2) = JWORD(5,I) GOTO 310 C 290 CONTINUE IPAIR(1,1) = JWORD(1,I) IPAIR(1,2) = JWORD(2,I) IPAIR(2,1) = JWORD(3,I) IPAIR(2,2) = JWORD(4,I) GOTO 310 C 300 CONTINUE IPAIR(1,1) = JWORD(1,I) IPAIR(1,2) = JWORD(3,I) IPAIR(2,1) = JWORD(2,I) IPAIR(2,2) = JWORD(4,I) 310 CONTINUE DO 400 I4 = 1,2 KM = 0 IF (IPAIR(I4,1).GT.MP) KM = KM + 1 IF (IPAIR(I4,2).GT.MP) KM = KM + 1 JJ1 = IPAIR(I4,1) JJ2 = IPAIR(I4,2) IF (KM-1) 330,370,400 C C ONE VARIABLE LINKED TO TWO CONSTANTS.FIX THE DIAGONAL MAT(I,I) C 330 CONTINUE JT1 = J1(JJ1) - 1 JT2 = J1(JJ2) - 1 JMIN = ABS(JT1-JT2) JMAX = JT1 + JT2 IF (MAT(I1,I1)-1) 340,400,350 C C FIRST TIME C 340 CONTINUE MAT(I1,I1) = (JMAX-JMIN)/2 + 1 JSUM(1,I1) = JMIN JSUM(2,I1) = JMAX GOTO 400 C C IF THERE ARE SEVERAL COUPLES OF CONSTANTS, TAKE THE MORE C STRINGENT COMBINATION C 350 CONTINUE JMIN = MAX(JMIN,JSUM(1,I1)) JMAX = MIN(JMAX,JSUM(2,I1)) IF (JMAX.GE.JMIN) GOTO 360 RECUP = ZERO GOTO 890 C 360 CONTINUE JSUM(1,I1) = JMIN JSUM(2,I1) = JMAX MAT(I1,I1) = (JMAX-JMIN)/2 + 1 GOTO 400 C C ONE VARIABLE LINKED TO ONE CONSTANT AND ONE VARIABLE NON DIAGONAL C ELEMENT C 370 CONTINUE JT1 = MIN(JJ1,JJ2) JT2 = MAX(JJ1,JJ2) - MP IF (JT2.GT.I1) GOTO 400 JT4 = J1(JT1) - 1 K = MAT(I1,JT2) IF (K.EQ.0) GOTO 390 DO 380 LL = 1,K IF (JT4.EQ.J12(LL,JT2,I1)) GOTO 400 380 CONTINUE 390 CONTINUE K = K + 1 IF (K.GT.MXCSVR) GOTO 400 MAT(I1,JT2) = K J12(K,JT2,I1) = JT4 400 CONTINUE 410 CONTINUE 420 CONTINUE C C REDUCE THE DIAGONAL ELEMENTS BY TAKING INTO ACCOUNT THE NON C DIAGONAL ELEMENTS, AND KEEP THE LATTER ONLY IF NEEDED C LOOP STRUCTURE MODIFIED AVOIDING JUMPS BACK IN WE'89AOUT24 C 430 CONTINUE ICHAN = 0 DO 560 I = 1,NSUM NOEL(I) = .TRUE. I2 = 1 I1 = I - 1 IF (I1.NE.0) GOTO 450 440 CONTINUE IF (I.EQ.NSUM) GOTO 560 I1 = NSUM I2 = I + 1 450 CONTINUE DO 550 J = I2,I1 IF (MAT(J,J).EQ.0) GOTO 550 COLD IF (I.EQ.1) GOTO 460 IF (I.LT.J) GOTO 460 !KAB 5/11/03 IF (MAT(I,J).EQ.0) GOTO 550 IK1 = I IK2 = J NOEL(I) = .FALSE. GOTO 470 C 460 CONTINUE IF (MAT(J,I).EQ.0) GOTO 550 IK1 = J IK2 = I C 470 CONTINUE JMIN1 = 0 JMAX1 = 1000 K = MAT(IK1,IK2) DO 490 L1 = 1,K L3 = MAT(J,J) JJ1 = JSUM(1,J) JND = J12(L1,IK2,IK1) JMIN = 1000 JMAX = 0 JMNP(L1) = 0 JMXP(L1) = 1000 DO 480 L2 = 1,L3 JMN = ABS(JND-JJ1) JMX = JND + JJ1 JMIN = MIN(JMN,JMIN) JMAX = MAX(JMX,JMAX) JMNP(L1) = MAX(JMN,JMNP(L1)) JMXP(L1) = MIN(JMX,JMXP(L1)) JJ1 = JJ1 + 2 480 CONTINUE JMIN1 = MAX(JMIN1,JMIN) JMAX1 = MIN(JMAX1,JMAX) 490 CONTINUE IF (MAT(I,I).NE.0) GOTO 500 JSUM(1,I) = JMIN1 JSUM(2,I) = JMAX1 MAT(I,I) = (JMAX1-JMIN1)/2 + 1 ICHAN = ICHAN + 1 GOTO 520 C 500 CONTINUE IF (JSUM(1,I).GE.JMIN1) GOTO 510 JSUM(1,I) = JMIN1 ICHAN = ICHAN + 1 510 CONTINUE IF (JSUM(2,I).LE.JMAX1) GOTO 520 JSUM(2,I) = JMAX1 ICHAN = ICHAN + 1 520 CONTINUE K1 = 0 DO 530 L1 = 1,K IF (JMNP(L1).LE.JSUM(1,I) .AND. A JMXP(L1).GE.JSUM(2,I)) GOTO 530 K1 = K1 + 1 J12(K1,IK2,IK1) = J12(L1,IK2,IK1) 530 CONTINUE IF (K1.EQ.K) GOTO 540 MAT(IK1,IK2) = K1 ICHAN = ICHAN + 1 540 CONTINUE MAT(IK2,IK1) = J12(1,IK2,IK1) C OUT GOTO JKM(171,172) -- SECTION REWRITTEN AT CECAM MEUDON 1989. W.E. C 550 CONTINUE IF (I1.NE.NSUM) GOTO 440 560 CONTINUE IF (ICHAN.NE.0) GOTO 430 C C 2.3 C C CARRY OUT THE SUMMATIONS. C DO 570 I = 1,NSUM LDIAG(I) = MAT(I,I) .EQ. 1 JSUM3(I) = 1 570 CONTINUE DO 580 I = 1,JWRD JWTEST(I) = 1 580 CONTINUE STOR = ZERO STOR1 = ONE NOLP = 0 IP = 1 590 CONTINUE NOLP = NOLP + 1 C 2.3.1 C FIND THE RANGE OF JSUM2(NOLP) C NOLP IS THE INDEX OF THE SUMMATION VARIABLE C JMIN = JSUM(1,NOLP) JMAX = JSUM(2,NOLP) IF (NOEL(NOLP)) GOTO 640 NO1 = NOLP - 1 DO 630 NJ = 1,NO1 IF (MAT(NOLP,NJ)-1) 630,600,610 600 CONTINUE JJ1 = MAT(NJ,NOLP) JJ2 = JSUM2(NJ) JMIN = MAX(JMIN,ABS(JJ2-JJ1)) JMAX = MIN(JMAX,JJ1+JJ2) GOTO 630 C 610 CONTINUE K = MAT(NOLP,NJ) JJ2 = JSUM2(NJ) DO 620 I = 1,K JJ1 = J12(I,NJ,NOLP) JMIN = MAX(JMIN,ABS(JJ2-JJ1)) JMAX = MIN(JMAX,JJ1+JJ2) 620 CONTINUE 630 CONTINUE 640 CONTINUE JSUM2(NOLP) = JMIN MAXLP(NOLP) = JMAX IF (LDIAG(NOLP)) JSUM3(NOLP) = 0 IF (NOLP.LT.NSUM) GOTO 590 DO 850 JJ = JMIN,JMAX,2 JSUM2(NSUM) = JJ C C 2.3.2 C DETERMINE WHICH RACAH COEFFICIENTS NEED RE-EVALUATING AND C SET JWTEST APPROPRIATELY C DO 670 J = IP,NSUM IF (JSUM3(J).LE.0) GOTO 670 I2 = JSUM6(J) DO 660 I1 = 1,I2 I3 = JSUM5(J,I1) JWTEST(I3) = 1 660 CONTINUE 670 CONTINUE DO 700 J = 1,JWRD IF (JWTEST(J).EQ.0) GOTO 700 JWJ = J + JWR DO 680 I = 1,6 I1 = JWORD(I,JWJ) IF (I1.LE.MP) THEN IST(I) = J1(I1) - 1 C ELSE I1 = I1 - MP - NFS IST(I) = JSUM2(I1) ENDIF C 680 CONTINUE CALL DRACAH(IST(1),IST(2),IST(3),IST(4),IST(5),IST(6),X1) WSTOR(J) = X1 IF (IBUG3.NE.1) GOTO 700 DO 690 I = 1,6 XJ1(I) = IST(I)*0.5D0 690 CONTINUE PRINT 3030, (XJ1(I),I=1,6),X1 700 CONTINUE C C 2.3.3 C FORM PRODUCT OF RACAH COEFFICIENTS,(2J+1) FACTORS AND (-1) C FACTORS IN STOR1 C DO 710 I = 1,JWRD STOR1 = STOR1*WSTOR(I) 710 CONTINUE C C IASTOR CONTAINS THE POWER OF (-1) WHICH IS COMMON TO ALL TERMS C IX2 = 0 C IJ6CP = 1 !*4 overflwo nrb DIJ6CP = 1.0D0 IF (J6CP.EQ.J6F) GOTO 740 JB = J6F + 1 DO 730 I = JB,J6CP I1 = J6P(I) - NFS DIJ6CP = DIJ6CP*(JSUM2(I1)+1) !IJ6CP->DIJ6CP nrb 730 CONTINUE 740 CONTINUE IF (J9CP.EQ.J9F) GOTO 760 JB = J9F + 1 DO 750 I = JB,J9CP I1 = J9P(I) - NFS DIJ6CP = DIJ6CP/(JSUM2(I1)+1) !IJ6CP->DIJ6CP nrb 750 CONTINUE 760 CONTINUE C STOR1 = STOR1*SQRT(DBLE(IJ6CP)) !*4 overflwo nrb STOR1 = STOR1*SQRT(DIJ6CP) IF (J7CP.EQ.J7F) GOTO 790 JB = J7F + 1 DO 780 I = JB,J7CP I1 = J7P(I) - NFS IX2 = IX2 + JSUM2(I1) 780 CONTINUE 790 CONTINUE IF (J8CP.EQ.J8F) GOTO 810 JB = J8F + 1 DO 800 I = JB,J8CP I1 = J8P(I) - NFS IX2 = IX2 + 2* JSUM2(I1) 800 CONTINUE 810 CONTINUE IF (MOD(IX2,2).NE.1) GOTO 820 IAS = -1 IX2 = IX2 + 1 820 CONTINUE IX2 = IX2/2 C C 2.3.4 C ADD TERM INTO STOR AND RESET STOR1 TO 1 READY FOR NEXT TERM C IF (MOD(IX2,2).EQ.1) STOR1 = -STOR1 STOR = STOR + STOR1 STOR1 = ONE NSUM1 = NSUM - 1 IF (NSUM1.EQ.0) GOTO 850 DO 830 IK = 1,NSUM1 JSUM3(IK) = 0 830 CONTINUE DO 840 IK = 1,JWRD JWTEST(IK) = 0 840 CONTINUE 850 CONTINUE 860 CONTINUE NOLP = NOLP - 1 IF (NOLP.EQ.0) GOTO 870 IF (LDIAG(NOLP)) GOTO 860 JSUM3(NOLP) = 1 JSUM2(NOLP) = JSUM2(NOLP) + 2 IF (JSUM2(NOLP).GT.MAXLP(NOLP)) GOTO 860 IP = NOLP C CC 2.3.5 PROCEED TO NEXT VARIABLE C GOTO 590 C 870 CONTINUE RECUP = RECUP*STOR IF (IBUG3.EQ.1) PRINT 3050,NPS,STOR,RECUP IF (ABS(RECUP).LT.EPSIL) GOTO 900 JWR = JWRD + JWR NFS = NSUM + NFS J6F = J6CP J7F = J7CP J8F = J8CP J9F = J9CP IASTOR = IASTOR + IAS C C 2.4 PROCEED TO NEXT SUM C IF (NPS.LT.NLSUM) GOTO 220 IASTOR = IASTOR/2 IF (MOD(IASTOR,2).NE.0) RECUP = -RECUP IF (IBUG3.EQ.1) PRINT 3020,RECUP 890 CONTINUE RETURN C NO SUMMATIONS. CHECK THAT THERE ARE NO INCONSISTENCIES. THEN C MULTIPLY BY (-1) FACTOR AND EXIT C 900 CONTINUE RECUP = ZERO GOTO 890 C 3000 FORMAT (' SUM NR.',I3) 3010 FORMAT (' NO SUMMATION. RECOUPLING COEFFICIENT=',G15.8) 3020 FORMAT (' RECOUPLING COEFFICIENT=',G15.8) 3030 FORMAT (6F5.1,10X,G15.8) 3040 FORMAT (' NUMBER OF INDEPENDENT SUMS:',I3) 3050 FORMAT (' SUM NR.',I2,' SUM VALUE=',G15.8,' RECUP=',G15.8) 3070 FORMAT (' NOT INVOLVING SUMMATION VARIABLE') 3080 FORMAT (/'0PRINTOUT FROM SUBROUTINE GENSUM'/ A '0VALUES OF ANGULAR MOMENTA IN *REAL* FORMAT'/ (14F5.1)) 3090 FORMAT (/' RACAH W FUNCTIONS(6J)'/' ARGUMENTS IN *REAL* FORMAT', A 18X,'VALUE') END C C C SUBROUTINE H0WTS(ISIG,ISIGP,TIMES,ICAL) C*********************************************************************** C C H0WTS C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) C PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL23=KFL2*3) C LOGICAL FAIL,FREE C COMMON /CONSTS/ZERO,TENTH,HALF,ONE,TWO,THREE,FOUR,SEVEN,ELEVEN,EPS COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /HOLD/J2SPIN(KFL23),J3SPIN(KFL23),J2ANG(KFL23),J3ANG(KFL23) COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 C C --- OFF-DIAGONAL ONE-ELECTRON HAMILTONIAN MATRIX ELEMENT C ICAL = 0 CALL SETM ICOUNT = 0 C C TEST THAT FINAL ANGULAR MOMENTA ARE EQUAL C DO 10 K = 2,3 IF (J1QN1(M6,K).NE.J1QN2(M6,K)) GOTO 460 10 CONTINUE C C --- DETERMINE INTERACTING SHELLS, ISIG ON L.H.S., ISIGP ON R.H.S., C FOR NON-ZERO ONE-ELECTRON MATRIX ELEMENT, N-1 ELECTRONS MUST BE C COMMON TO BOTH SIDES. THUS THE SUM OF N(I) = NOSH1(I)-NOSH2(I), C I=1,IHSH MUST BE EQUAL TO 0 OR 2 . THUS AT NO STAGE CAN C N(I) BE GREATER THAN 1 . IF THIS SUM IS ZERO, THE TWO C CONFIGURATIONS ARE MADE UP FROM THE SAME ELECTRONS, WITH TWO C DIFFERENT COUPLING SCHEMES. SINCE THE SPHERICAL HARMONICS ARE C EIGENFUNCTIONS OF DEL**2 , THE ORTHOGONALITY OF THE TWO COUPLING C SCHEMES WILL BE MAINTAINED AND ORTHOGONALITY GIVES A ZERO RESULT. C DO 50 I = 1,IHSH N = NOSH1(I) - NOSH2(I) IF (ABS(N).GT.1) GOTO 460 IF (N) 20,50,30 20 CONTINUE ISIGP = I GOTO 40 C 30 CONTINUE ISIG = I 40 CONTINUE ICOUNT = ICOUNT + 1 50 CONTINUE IF (ICOUNT.NE.2) GOTO 460 IF (IBUG4.LT.2) GOTO 60 WRITE (IWRITE,3000) ISIG,ISIGP 60 CONTINUE LSIG = LJ(ISIG) LSIGP = LJ(ISIGP) C C THE ANGULAR MOMENTUM OF THE INTERACTING ELECTRONS MUST BE EQUAL C IF (LSIG-LSIGP) 460,70,460 C C THE SPECTATOR SHELLS MUST HAVE MATCHING QUANTUM NUMBERS C 70 CONTINUE DO 90 J = 1,IHSH IF (J.EQ.ISIG .OR. J.EQ.ISIGP) GOTO 90 DO 80 K = 1,3 IF (J1QN1(J,K).NE.J1QN2(J,K)) GOTO 460 80 CONTINUE 90 CONTINUE C C --- TEST ON TRIANGULAR CONDITIONS C IN3 = 2*LSIG DO 100 K = 2,3 IN1 = J1QN1(ISIG,K) IN2 = J1QN2(ISIG,K) IN4 = J1QN1(ISIGP,K) IN5 = J1QN2(ISIGP,K) IF (IN1.GT. (IN2+IN3) .OR. IN1.LT.ABS(IN2-IN3)) GOTO 460 IF (IN4.GT. (IN5+IN3) .OR. IN4.LT.ABS(IN5-IN3)) GOTO 460 IN3 = 1 100 CONTINUE C C --- CALCULATE FRACTIONAL PARENTAGE COEFFICIENTS C TIMES = ONE IF (LSIG.EQ.0) GOTO 110 N = NOSH1(ISIG) IVI = J1QN1(ISIG,1) ILI = (J1QN1(ISIG,2)-1)/2 ISI = J1QN1(ISIG,3) IVJ = J1QN2(ISIG,1) ILJ = (J1QN2(ISIG,2)-1)/2 ISJ = J1QN2(ISIG,3) CALL CFP(LSIG,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) TIMES = TIMES*COEFP IF (ABS(TIMES).LT.EPS) GOTO 460 N = NOSH2(ISIGP) IVI = J1QN2(ISIGP,1) ILI = (J1QN2(ISIGP,2)-1)/2 ISI = J1QN2(ISIGP,3) IVJ = J1QN1(ISIGP,1) ILJ = (J1QN1(ISIGP,2)-1)/2 ISJ = J1QN1(ISIGP,3) CALL CFP(LSIG,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) TIMES = TIMES*COEFP IF (ABS(TIMES).LT.EPS) GOTO 460 C C --- SET UP J2 AND J3 ARRAYS C 110 CONTINUE M1 = IHSH - 2 M2 = M6 - 2 J2(1,1) = ISIG J2(1,2) = M11 J2(1,3) = M9 J3(1,1) = ISIGP J3(1,2) = M11 J3(1,3) = M10 IF (ISIG.EQ.1) GOTO 120 J2(2,1) = 1 GOTO 130 C 120 CONTINUE J2(2,1) = M9 130 CONTINUE IF (ISIG.EQ.2) GOTO 140 J2(2,2) = 2 GOTO 150 C 140 CONTINUE J2(2,2) = M9 150 CONTINUE J2(2,3) = M4 IF (ISIGP.EQ.1) GOTO 160 J3(2,1) = 1 GOTO 170 C 160 CONTINUE J3(2,1) = M10 170 CONTINUE IF (ISIGP.EQ.2) GOTO 180 J3(2,2) = 2 GOTO 190 C 180 CONTINUE J3(2,2) = M10 190 CONTINUE J3(2,3) = M7 IF (IHSH.LT.3) GOTO 260 DO 250 J = 3,IHSH J2(J,1) = M1 + J J2(J,3) = M1 + J + 1 J3(J,1) = M2 + J IF (J.EQ.IHSH) GOTO 200 J3(J,3) = M2 + J + 1 GOTO 210 C 200 CONTINUE J3(J,3) = M1 + J + 1 210 CONTINUE IF (J.EQ.ISIG) GOTO 220 J2(J,2) = J GOTO 230 C 220 CONTINUE J2(J,2) = M9 230 CONTINUE IF (J.EQ.ISIGP) GOTO 240 J3(J,2) = J GOTO 250 C 240 CONTINUE J3(J,2) = M10 250 CONTINUE C C --- STORE J2 AND J3 ARRAYS FOR USE IN SPIN RECOUPLING COEFFICIENT C 260 CONTINUE I1 = 0 DO 280 J = 1,IHSH DO 270 K = 1,3 I1 = I1 + 1 J2SPIN(I1) = J2(J,K) J3SPIN(I1) = J3(J,K) 270 CONTINUE 280 CONTINUE C C --- ORBITAL RECOUPLING COEFFICIENT C J1(M11) = LSIG + LSIG + 1 K = 2 C C --- SET J1 ARRAY C 290 CONTINUE DO 310 J = 1,IHSH IF (ISIG.EQ.J) GOTO 300 J1(J) = J1QN1(J,K) GOTO 310 C 300 CONTINUE J1(J) = J1QN2(ISIG,K) 310 CONTINUE DO 320 J = M4,M6 J1(J) = J1QN1(J,K) 320 CONTINUE DO 330 J = M7,M8 J1(J) = J1QN2(J-M3,K) 330 CONTINUE J1(M9) = J1QN1(ISIG,K) J1(M10) = J1QN2(ISIGP,K) NJ1S = M11 NJ23S = M4 DO 340 J = 1,M11 FREE(J) = .FALSE. 340 CONTINUE IF (IBUG4.LT.2 .OR. IBUG3.EQ.1) GOTO 360 WRITE (IWRITE,3010) (J1(J),J=1,M11) WRITE (IWRITE,3020) DO 350 J = 1,IHSH WRITE (IWRITE,3030) (J2(J,KL),KL=1,3), (J3(J,KL),KL=1,3) 350 CONTINUE C C --- CALCULATE RECOUPLING COEFFICIENT C 360 CONTINUE CALL NJGRAF(RECUP,FAIL) TIMES = TIMES*RECUP IF (IBUG4.LT.2) GOTO 370 WRITE (IWRITE,3040) TIMES,RECUP 370 CONTINUE IF (ABS(TIMES).LT.EPS) GOTO 460 C C --- SPIN RECOUPLING COEFFICIENT C IF (K.EQ.3) GOTO 400 J1(M11) = 2 K = 3 I1 = 0 DO 390 J = 1,IHSH DO 380 KK = 1,3 I1 = I1 + 1 J2(J,KK) = J2SPIN(I1) J3(J,KK) = J3SPIN(I1) 380 CONTINUE 390 CONTINUE GOTO 290 C C --- INCLUDE MULTIPLICATIVE FACTORS C 400 CONTINUE IDELP = 0 IF (ISIG-ISIGP) 410,450,430 410 CONTINUE JSIG = ISIG + 1 DO 420 J = JSIG,ISIGP IDELP = IDELP + NOSH1(J) 420 CONTINUE GOTO 450 C 430 CONTINUE JSIGP = ISIGP + 1 DO 440 J = JSIGP,ISIG IDELP = IDELP + NOSH2(J) 440 CONTINUE 450 CONTINUE TIMES = TIMES* (-ONE)**IDELP* A SQRT(DBLE(NOSH1(ISIG)*NOSH2(ISIGP))) ICAL = 1 IF (IBUG4.LT.2) GOTO 460 WRITE (IWRITE,3050) TIMES,IDELP 460 CONTINUE C 3000 FORMAT (' ISIG =',I3,' ISIGP =',I3) 3010 FORMAT (' J1 =',I8,36I3) 3020 FORMAT (' J2',18X,'J3') 3030 FORMAT (3I5,I10,2I5) 3040 FORMAT (' TIMES =',F10.6,' RECUP =',F10.6) 3050 FORMAT (' TIMES =',F10.6,' IDELP =',I3) END C C C SUBROUTINE HOUSE(N,A,LENGTH,P,R,ARRAY1) C C ACCEPTS THE UPPER TRIANGLE OF ELEMENTS OF A SYMMETRIC MATRIX, C STORED IN THE LINEAR ARRAY A, AND USING THE HOUSEHOLDER METHOD C REDUCES THIS TO TRI-DIAGONAL FORM, STORING THE NEW MAIN C DIAGONAL ELEMENTS IN POSITION IN A AND ALSO IN THE ARRAY R, C AND THE SUPER-DIAGONAL ELEMENTS IN THE ARRAY P. DETAILS C OF THE TRANSFORMING MATRICES ARE OVERWRITTEN IN THE REDUNDANT C SPACE OF A. N AND LENGTH ARE AS DEFINED IN SUBROUTINE HSLDR. C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION A(LENGTH),P(N),R(N),ARRAY1(N) DATA EPSI/1.0D-16/,ZERO,HALF/0.0D0,0.5D0/ C N1 = N - 1 N2 = N - 2 C C EACH TIME THROUGH THIS LOOP ONE MORE ROW OF THE MATRIX IS C TRANSFORMED TO TRI-DIAGONAL FORM C DO 150 K = 1,N2 DO 10 J = K,N P(J) = 0.0D0 10 CONTINUE C C KB=KC-1 IS THE NUMBER OF ELEMENTS IN THE FIRST (K-1) ROWS OF THE C UPPER TRIANGLE C NJ = N - K KC = ((K-1)* (NJ+N+2))/2 + 1 C C THE SQUARE ROOT OF THE SUM OF THE SQUARES OF THE REMAINING C OFF DIAGONAL ELEMENTS IN ROW K IS FOUND AND STORED IN SUM C SUM = 0.0D0 DO 20 J = 1,NJ SUM = A(KC+J)*A(KC+J) + SUM 20 CONTINUE SUM = SQRT(SUM) IF (SUM.LT.EPSI) GOTO 130 C C SUM IS GIVEN THE SAME SIGN AS THE SUPER DIAGONAL ELEMENT IN ROW K C IF (A(KC+1).LT.ZERO) SUM = -SUM BKR = SUM* (SUM+A(KC+1)) C C THE FIRST ELEMENT OF THE VECTOR FROM WHICH THE TRANSFORMING C MATRIX IS DERIVED IS OVERWRITTEN ON THE OLD SUPER DIAGONAL C ELEMENT IN ROW K. THE REMAINING ELEMENTS ARE ALREADY IN POSITION C IN ROW K. C A(KC+1) = A(KC+1) + SUM C C THE SUPER DIAGONAL ELEMENT IN ROW K OF THE NEW TRI-DIAGONAL C MATRIX IS STORED IN ARRAY P C P(K) = -SUM C C KD IS THE NUMBER OF ELEMENTS IN THE FIRST K ROWS OF THE UPPER C TRIANGLE C KD = (K* (NJ+N+1))/2 C C THE TRANSFORMATION DERIVED FROM THE NJ=(N-K) VECTOR ELEMENTS C STORED IN A IS NOW APPLIED TO THE LAST NJ ROWS OF THE MATRIX C (THE LAST NJ SPACES OF THE ARRAY P ARE SUCCESSIVELY C OVERWRITTEN IN THE PROCESS). C LOL = KC - K + 1 DO 70 M = K,N1 NM = N - M MO = ((NM+N+1)*M)/2 M1 = M - 1 LO = LOL + M1 SUM = P(1+M) DO 30 L = 1,NM SUM = SUM + A(MO+L)*A(LO+L) 30 CONTINUE IF (M.EQ.K) GOTO 60 JM1 = KD + M - K + 1 DO 40 L = K,M1 ARRAY1(L) = A(JM1) JM1 = JM1 + N1 - L 40 CONTINUE DO 50 L = K,M1 SUM = SUM + ARRAY1(L)*A(LOL+L) 50 CONTINUE 60 P(1+M) = SUM/BKR 70 CONTINUE C SUM = 0.0D0 DO 80 L = 1,NJ SUM = SUM + A(KC+L)*P(K+L) 80 CONTINUE SUM = SUM/BKR DO 90 L = 1,NJ P(K+L) = P(K+L) - A(KC+L)*SUM*HALF 90 CONTINUE DO 120 I = 1,NJ AIM = A(KC+I) PIN = P(K+I) DO 100 J = I,NJ ARRAY1(J) = -AIM*P(K+J) - PIN*A(KC+J) 100 CONTINUE KE = KD + (I-1)*NJ - (I* (I-1))/2 DO 110 J = I,NJ A(KE+J) = A(KE+J) + ARRAY1(J) 110 CONTINUE 120 CONTINUE GOTO 150 C 130 CONTINUE P(K) = -EPSI DO 140 I = 1,NJ A(I+KC) = 0.0D0 140 CONTINUE 150 CONTINUE C C THE LAST SUPER DIAGONAL ELEMENT IS ENTERED INTO THE ARRAY P. C ILK = ((N+1)*N)/2 - 1 P(N1) = A(ILK) C C THE MAIN DIAGONAL ELEMENTS PICKED OUT FROM THE ARRAY A ARE C STORED IN THE ARRAY R. C DO 160 I = 1,N ILK = ((I-1)* (2*N-I+2))/2 + 1 R(I) = A(ILK) 160 CONTINUE C END C C C SUBROUTINE HSLDR(N,A,LENGTH,EPSI,EIG,X,NO,P,MSV) C*********************************************************************** C C HSLDR C C*********************************************************************** C C THIS SUBROUTINE ACCEPTS THE UPPER TRIANGLE OF AN N*N SYMMETRIC C MATRIX AND ON THE FIRST CALL DETERMINES ALL THE EIGENVALUES AND C THE FIRST EIGENVECTOR. ON EACH FURTHER CALL ONE MORE OF THE C REMAINING EIGENVECTORS IS CALCULATED. C C DEFINITION OF THE ARGUMENTS. C N........... THE DEGREE OF THE SYMMETRIC MATRIX TO BE C DIAGONALISED. C A........... THE LINEAR ARRAY CONTAINING THE UPPER TRIANGLE OF C THE ORIGINAL MATRIX, OVERWRITTEN ON RETURN BY THE C MAIN DIAGONAL ELEMENTS OF THE TRI-DIAGONAL MATRIX C AND DETAILS OF THE TRANSFORMING MATRICES. C LENGTH...... =(N*(N+1))/2 , THE SIZE OF THE ARRAY A. C EPSI........ THE ACCURACY TO WHICH THE EIGENVALUES ARE TO BE C DETERMINED. C EIG......... THIS ARRAY CONTAINS THE EIGENVALUES ON RETURN. C X........... THIS ARRAY CONTAINS ONE EIGENVECTOR ON RETURN. C NO.......... THIS RUNS FROM 1 TO N AND SPECIFIES WHICH C EIGENVECTOR IS STORED IN X ON RETURN C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION A(LENGTH),EIG(N),X(N),P(MSV,9) C C THE SIZE OF THE MATRIX WAS LIMITED BY THE SIZE OF THE INTERNAL C ARRAYS IN THE SAVE1, SAVE2 AND SAVE3 COMMON BLOCKS. C IF (N.LE.1) THEN X(1) = 1.0D0 EIG(1) = A(1) P(1,2) = A(1) GOTO 60 C ELSE IF (N.EQ.2) THEN P(1,1) = A(2) P(1,2) = A(1) P(2,2) = A(3) GOTO 10 C ENDIF C C IF THE FIRST EIGENVECTOR HAS ALREADY BEEN FOUND THE C TRI-DIAGONALISING AND EIGENVALUE SUBROUTINES ARE SKIPPED ROUND. C IF (NO.NE.1) GOTO 20 C C THE TRI-DIAGONALISING SUBROUTINE IS ENTERED. C CALL HOUSE(N,A,LENGTH,P(1,1),P(1,2),P(1,9)) C C THE ELEMENTS OF THE TRI-DIAGONAL MATRIX ARE USED TO DETERMINE C THE EIGENVALUES. C 10 CONTINUE CALL EIGEN(N,EIG,EPSI,P(1,1),P(1,2),P(1,3),P(1,4)) C C THE EIGENVECTOR OF THE TRI-DIAGONAL MATRIX CORRESPONDING TO A C PARTICULAR EIGENVALUE IS DETERMINED. C 20 CONTINUE CALL VECTOR(N,EIG,X,NO,P(1,1),P(1,2),P(1,3),P(1,4),P(1,5),P(1,6), A P(1,7),P(1,8),P(1,9)) IF (N.EQ.2) GOTO 30 C C THE CORRESPONDING EIGENVECTOR OF THE ORIGINAL MATRIX IS FOUND. C CALL EIGVEC(N,A,LENGTH,X,P) C C NORMALIZING THE EIGENVECTOR. C 30 CONTINUE APP = 0.0D0 DO 40 I = 1,N APP = APP + X(I)*X(I) 40 CONTINUE APP = SQRT(APP) DO 50 I = 1,N X(I) = X(I)/APP 50 CONTINUE C 60 CONTINUE C END C C C SUBROUTINE INTAB C C THIS SUBROUTINE CALLED AT THE END OF DIAGRM,FIXES THE ARRAYS IH C AND IL-SO TO SPEAK HARDWARE AND LOGICAL ADDRESSES OF TRIADS IN C JDIAG.ALSO DETERMINES THE NUMBER OF FREE ENDS NFREE AND THEIR C LOCATION ITFREE. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C COMMON /BUILD/IAL(KFL2B),IF1,IF2,NODE C LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C DO 10 I = 1,M IAL(I) = 1 10 CONTINUE DO 20 I = IFIRST,ILAST J = JDIAG(I,1) K = IAL(J) TAB1(J,K) = I IAL(J) = K + 1 20 CONTINUE IFR = IFIRST - 1 DO 30 I = IFIRST,ILAST IT = I - IFR IL(I) = IT IH(IT) = I 30 CONTINUE J = JDIAG(IFIRST,3) K = IAL(J) IF (K-1) 50,50,40 40 CONTINUE TAB1(J,2) = TAB1(J,1) 50 CONTINUE TAB1(J,1) = IFIRST IAL(J) = 3 J = JDIAG(ILAST,2) TAB1(J,2) = ILAST IAL(J) = 3 NFREE = 0 DO 60 I = IFIRST,ILAST J = JDIAG(I,1) IF (IAL(J).EQ.3) GOTO 60 NFREE = NFREE + 1 ITT = ILAST + NFREE TAB1(J,2) = ITT IL(ITT) = NFREE*1000 ITFREE(NFREE) = I 60 CONTINUE C END C C C SUBROUTINE INTACT(L,LP,IEQUIV) C*********************************************************************** C C INTACT C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /ENAV/COEFCT(5),NINTS,KVALUE(5) C C * THIS SUBROUTINE CONTAINS REAL CONSTANTS IN THE TEXT * C C THIS SUBROUTINE GIVES THE INTERACTION ENERGY BETWEEN TWO SHELLS, C ONE WITH ORBITAL ANGULAR MOMENTUM L , THE OTHER WITH ORBITAL C ANGULAR MOMENTUM LP . NOTICE THAT THE FIRST TERM OF THIS C INTERACTION ENERGY IS ALWAYS F0(L,LP) AND THIS IS NOT GIVEN C IN THIS SUBROUTINE. THUS ONLY THE EXTRA TERMS ARE HERE PRODUCED. C FOR EQUIVALENT ELECTRONS (IEQUIV = 1) , THERE WILL BE FK C INTEGRALS ONLY. FOR NON-EQUIVALENT ELECTRONS (IEQUIV = 2) , C THERE WILL BE GK INTEGRALS ONLY. C C THE EXPRESSIONS FOR THE INTERACTION ENERGIES INVOLVING SHELLS WITH C L.LE.3 ARE GIVEN BY J.C. SLATER , QUANTUM THEORY OF ATOMIC C STRUCTURE, VOL. I, EQUATIONS (14.20) AND (14.22). IN THE LAST C OF HIS EQUATIONS (14.22) , A TERM -1/14 G0(F,FP) IS OMITTED. C THIS IS INCLUDED BELOW. THE INTERACTION ENERGIES FOR G-ELECTRON C SHELLS MAY BE EVALUATED USING HIS EQUATIONS (13.12), (13.17), C (14.19), AND (14,21) C IF (L.GT.4 .OR. LP.GT.4) GOTO 310 GOTO (10,90),IEQUIV C C === EQUIVALENT ELECTRONS C 10 CONTINUE IF (L.NE.LP) GOTO 20 GOTO 30 C 20 CONTINUE WRITE (IWRITE,3000) STOP 30 CONTINUE IF (L.GT.0) GOTO 40 C C S ELECTRONS C NINTS = 0 RETURN C 40 CONTINUE GOTO (50,60,70,80),L C C P ELECTRONS C 50 CONTINUE NINTS = 1 KVALUE(1) = 2 COEFCT(1) = -0.08D0 RETURN C C D ELECTRONS C 60 CONTINUE NINTS = 2 KVALUE(1) = 2 KVALUE(2) = 4 COEFCT(1) = -0.031746032D0 COEFCT(2) = -0.031746032D0 RETURN C C F ELECTRONS C 70 CONTINUE NINTS = 3 KVALUE(1) = 2 KVALUE(2) = 4 KVALUE(3) = 6 COEFCT(1) = -0.020512821D0 COEFCT(2) = -0.013986014D0 COEFCT(3) = -0.017930787D0 RETURN C C G ELECTRONS C 80 CONTINUE NINTS = 4 KVALUE(1) = 2 KVALUE(2) = 4 KVALUE(3) = 6 KVALUE(4) = 8 COEFCT(1) = -0.015278839D0 COEFCT(2) = -0.009519892D0 COEFCT(3) = -0.008227067D0 COEFCT(4) = -0.011856655D0 RETURN C C --- NON-EQUIVALENT ELECTRONS C 90 CONTINUE IF (L.GT.LP) GOTO 100 L1 = L L2 = LP GOTO 110 C 100 CONTINUE L1 = LP L2 = L 110 CONTINUE L1D = L1 + 1 L2D = L2 + 1 GOTO (120,180,230,270,300),L1D C 120 CONTINUE NINTS = 1 KVALUE(1) = L2 GOTO (130,140,150,160,170),L2D C C S - S INTERACTION C 130 CONTINUE COEFCT(1) = -0.5D0 RETURN C C S - P INTERACTION C 140 CONTINUE COEFCT(1) = -0.166666667D0 RETURN C C S - D INTERACTION C 150 CONTINUE COEFCT(1) = -0.1D0 RETURN C C S - F INTERACTION C 160 CONTINUE COEFCT(1) = -0.071428571D0 RETURN C C S - G INTERACTION C 170 CONTINUE COEFCT(1) = -0.055555556D0 RETURN C 180 CONTINUE NINTS = 2 KVALUE(1) = L2 - 1 KVALUE(2) = L2D GOTO (190,200,210,220),L2 C C P - P INTERACTION C 190 CONTINUE COEFCT(1) = -0.166666667D0 COEFCT(2) = -0.066666667D0 RETURN C C P - D INTERACTION C 200 CONTINUE COEFCT(1) = -0.066666667D0 COEFCT(2) = -0.042857143D0 RETURN C C P - F INTERACTION C 210 CONTINUE COEFCT(1) = -0.042857143D0 COEFCT(2) = -0.031746032D0 RETURN C C P - G INTERACTION C 220 CONTINUE COEFCT(1) = -0.031746032D0 COEFCT(2) = -0.025252525D0 RETURN C 230 CONTINUE NINTS = 3 KVALUE(1) = L2 - 2 KVALUE(2) = L2 KVALUE(3) = L2 + 2 L3 = L2 - 1 GOTO (240,250,260),L3 C C D - D INTERACTION C 240 CONTINUE COEFCT(1) = -0.1D0 COEFCT(2) = -0.028571429D0 COEFCT(3) = -0.028571429D0 RETURN C C D - F INTERACTION C 250 CONTINUE COEFCT(1) = -0.042857143D0 COEFCT(2) = -0.019047619D0 COEFCT(3) = -0.021645022D0 RETURN C C D - G INTERACTION C 260 CONTINUE COEFCT(1) = -0.028571429D0 COEFCT(2) = -0.014430014D0 COEFCT(3) = -0.017482517D0 RETURN C 270 CONTINUE NINTS = 4 KVALUE(1) = L2 - 3 KVALUE(2) = L2 - 1 KVALUE(3) = L2 + 1 KVALUE(4) = L2 + 3 L3 = L2 - 2 GOTO (280,290),L3 C C F - F INTERACTION C 280 CONTINUE COEFCT(1) = -0.071428571D0 COEFCT(2) = -0.019047619D0 COEFCT(3) = -0.012987013D0 COEFCT(4) = -0.016650017D0 RETURN C C F - G INTERACTION C 290 CONTINUE COEFCT(1) = -0.031746032D0 COEFCT(2) = -0.012987013D0 COEFCT(3) = -0.009990010D0 COEFCT(4) = -0.013597514D0 RETURN C C G - G INTERACTION C 300 CONTINUE NINTS = 5 KVALUE(1) = 0 KVALUE(2) = 2 KVALUE(3) = 4 KVALUE(4) = 6 KVALUE(5) = 8 COEFCT(1) = -0.055555556D0 COEFCT(2) = -0.014430014D0 COEFCT(3) = -0.008991009D0 COEFCT(4) = -0.007770008D0 COEFCT(5) = -0.011197952D0 RETURN C C --- IF ANGULAR MOMENTUM VALUES ARE TOO LARGE ----- C 310 CONTINUE WRITE (IWRITE,3010) L,LP STOP C 3000 FORMAT ( A' ERROR IN INTACT - EQUIVALENT ELECTRONS CALLED FOR ?????RING ORBI BTAL ANGULAR MOMENTUM') 3010 FORMAT (//' THE ORBITAL ANGULAR MOMENTUM VALUES, WHICH ARE',2I5, A 6X,'ARE TOO LARGE FOR THE CODING OF INTACT'//) END C C C SUBROUTINE J23ANG(IRHO,ISIG,IRHOP,ISIGP,JANGDI) C C === SETS UP J2 AND J3 ARRAYS FOR DIRECT ANGULAR INTEGRAL CALL OF NJSYM C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) C PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL23=KFL2*3) C LOGICAL FREE C COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) COMMON /HOLD/J2SPIN(KFL23),J3SPIN(KFL23),J2ANG(KFL23),J3ANG(KFL23) C C HAVE THE J2 AND J3 ARRAYS ALREADY BEEN SET. IF NOT, THEN GOTO 2 C IF (JANGDI) 10,10,170 C C --- ROWS 3 TO M4 OF SPIN J2 AND J3 ARE SAME AS ROWS 4 TO (M4+1) OF C ANGULAR J2 AND J3 C 10 CONTINUE I1 = 6 DO 30 J = 3,M4 JP1 = J + 1 DO 20 K = 1,3 I1 = I1 + 1 J2(JP1,K) = J2SPIN(I1) J3(JP1,K) = J3SPIN(I1) 20 CONTINUE 30 CONTINUE C C --- SET ROWS 1, 2 AND 3 C IF (M1) 40,50,40 40 CONTINUE J2(3,1) = ISIG GOTO 60 C 50 CONTINUE J2(3,1) = M9 60 CONTINUE IF (M2) 70,80,70 70 CONTINUE J3(3,1) = ISIGP GOTO 90 C 80 CONTINUE J3(3,1) = M11 90 CONTINUE J2(2,3) = M9 J2(2,1) = IRHO J2(2,2) = M13 J2(1,3) = M14 J2(3,2) = M14 J2(3,3) = M10 J2(1,1) = M16 J2(1,2) = M17 J3(3,2) = M16 J3(3,3) = M12 J3(1,2) = M13 J3(1,1) = M17 J3(1,3) = M15 J3(2,3) = M11 J3(2,1) = IRHOP J3(2,2) = M15 C C --- STORE J2 AND J3 FOR USE IN CALCULATING THE EXCHANGE TERM C I1 = 0 DO 110 J = 1,M5 DO 100 K = 1,3 I1 = I1 + 1 J2ANG(I1) = J2(J,K) J3ANG(I1) = J3(J,K) 100 CONTINUE 110 CONTINUE JANGDI = 1 120 CONTINUE IF (IBUG1-1) 160,160,130 C C PRINT-OUT OF VALUES IN NJSYM IF IBUG3=1 C 130 CONTINUE IF (IBUG3-1) 140,160,140 140 CONTINUE WRITE (IWRITE,3000) DO 150 J = 1,M5 WRITE (IWRITE,3010) (J2(J,K),K=1,3), (J3(J,K),K=1,3) 150 CONTINUE 160 CONTINUE RETURN C C --- SET J2 AND J3 ARRAYS FROM STORE OF PREVIOUS CALCULATIONS C 170 CONTINUE I1 = 0 DO 190 J = 1,M5 DO 180 K = 1,3 I1 = I1 + 1 J2(J,K) = J2ANG(I1) J3(J,K) = J3ANG(I1) 180 CONTINUE 190 CONTINUE GOTO 120 C 3000 FORMAT (' J2',18X,'J3') 3010 FORMAT (3I5,I10,2I5) END C C C SUBROUTINE J23SPN(IRHO,ISIG,IRHOP,ISIGP,JSNDIR) C C === SET UP THE J2 AND J3 ARRAYS FOR THE DIRECT SPIN INTEGRAL CALL C OF NJSYM C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) C PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL23=KFL2*3) C LOGICAL FREE C COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) COMMON /HOLD/J2SPIN(KFL23),J3SPIN(KFL23),J2ANG(KFL23),J3ANG(KFL23) C C HAVE THE J2 AND J3 ARRAYS ALREADY BEEN SET. IF NOT, THEN GOTO 2 C IF (JSNDIR) 10,10,540 C C --- SET THIRD ROW OF J2 AND J3 C 10 CONTINUE IF (IRHO-1) 20,30,20 20 CONTINUE J2(3,1) = 1 GOTO 60 C 30 CONTINUE IF (M1) 50,40,50 40 CONTINUE J2(3,1) = M10 GOTO 110 C 50 CONTINUE J2(3,1) = M9 GOTO 110 C 60 CONTINUE IF (IRHO-2) 70,80,70 70 CONTINUE J2(3,2) = 2 GOTO 120 C 80 CONTINUE IF (M1) 90,100,90 90 CONTINUE J2(3,2) = M9 GOTO 120 C 100 CONTINUE J2(3,2) = M10 GOTO 120 C 110 CONTINUE IF (ISIG-2) 70,100,70 120 CONTINUE J2(3,3) = M4 IF (IRHOP-1) 130,140,130 130 CONTINUE J3(3,1) = 1 GOTO 170 C 140 CONTINUE IF (M2) 150,160,150 150 CONTINUE J3(3,1) = M11 GOTO 220 C 160 CONTINUE J3(3,1) = M12 GOTO 220 C 170 CONTINUE IF (IRHOP-2) 180,190,180 180 CONTINUE J3(3,2) = 2 GOTO 230 C 190 CONTINUE IF (M2) 210,200,210 200 CONTINUE J3(3,2) = M12 GOTO 230 C 210 CONTINUE J3(3,2) = M11 GOTO 230 C 220 CONTINUE IF (ISIGP-2) 180,200,180 230 CONTINUE J3(3,3) = M7 C C --- SET ROWS 4,5,.. ETC. C IF (M4-4) 400,240,240 240 CONTINUE DO 390 J = 4,M4 J2(J,1) = M4 + J - 4 J2(J,3) = M4 + J - 3 IF (ISIG+1-J) 250,280,250 250 CONTINUE IF (M1) 260,270,260 260 CONTINUE IF (IRHO+1-J) 270,290,270 270 CONTINUE J2(J,2) = J - 1 GOTO 300 C 280 CONTINUE J2(J,2) = M10 GOTO 300 C 290 CONTINUE J2(J,2) = M9 300 CONTINUE J3(J,1) = M7 + J - 4 IF (J-M4) 320,310,320 310 CONTINUE J3(J,3) = J2(J,3) GOTO 330 C 320 CONTINUE J3(J,3) = M7 + J - 3 330 CONTINUE IF (ISIGP+1-J) 340,370,340 340 CONTINUE IF (M2) 350,360,350 350 CONTINUE IF (IRHOP+1-J) 360,380,360 360 CONTINUE J3(J,2) = J - 1 GOTO 390 C 370 CONTINUE J3(J,2) = M12 GOTO 390 C 380 CONTINUE J3(J,2) = M11 390 CONTINUE C C --- SET FIRST TWO ROWS, CORRESPONDING TO COUPLING OF INTERACTING C ELECTRONS WITHIN THEIR SHELLS C 400 CONTINUE J2(2,3) = M10 J2(1,2) = M13 J2(2,2) = M14 J2(1,3) = M9 IF (M1) 410,420,410 410 CONTINUE J2(1,1) = IRHO J2(2,1) = ISIG GOTO 430 C 420 CONTINUE J2(1,1) = ISIG J2(2,1) = M9 430 CONTINUE J3(2,3) = M12 J3(1,2) = M13 J3(2,2) = M14 J3(1,3) = M11 IF (M2) 440,450,440 440 CONTINUE J3(1,1) = IRHOP J3(2,1) = ISIGP GOTO 460 C 450 CONTINUE J3(1,1) = ISIGP J3(2,1) = M11 C C --- STORE J2,J3 ARRAYS FOR USE IN CALCULATING EXCHANGE INTEGRAL C 460 CONTINUE I1 = 0 DO 480 J = 1,M4 DO 470 K = 1,3 I1 = I1 + 1 J2SPIN(I1) = J2(J,K) J3SPIN(I1) = J3(J,K) 470 CONTINUE 480 CONTINUE JSNDIR = 1 490 CONTINUE IF (IBUG1-1) 530,530,500 C C PRINT-OUT OF VALUES IN NJSYM IF IBUG3=1 C 500 CONTINUE IF (IBUG3-1) 510,530,510 510 CONTINUE WRITE (IWRITE,3000) DO 520 J = 1,M4 WRITE (IWRITE,3010) (J2(J,K),K=1,3), (J3(J,K),K=1,3) 520 CONTINUE 530 CONTINUE RETURN C C --- SET J2 AND J3 ARRAYS FROM STORE OF PREVIOUS CALCULATIONS C 540 CONTINUE I1 = 0 DO 560 J = 1,M4 DO 550 K = 1,3 I1 = I1 + 1 J2(J,K) = J2SPIN(I1) J3(J,K) = J3SPIN(I1) 550 CONTINUE 560 CONTINUE GOTO 490 C 3000 FORMAT (' J2',18X,'J3') 3010 FORMAT (3I5,I10,2I5) END C C C SUBROUTINE KNJ(JD6C,JD7C,JD8C,JD9C,JDWC,JD6,JD7,JD8,JD9,JDW,JDDEL, A LDDEL,MDP,JD6P,JD7P,JD8P,JD9P,JDWORD,NDLSUM, B NDBJ,NDB6J,KD6CP,KD7CP,KD8CP,KD9CP,JDSUM4,JDSUM5, C JDSUM6,INVD6J) C PARAMETER (KFL1=100) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20,KFLS=12, A KFLN=10,KFLV=40) C DIMENSION JD6(KFL6),JD7(KFL7),JD8(KFL8),JD9(KFL9),JDW(6,KFLW), A LDDEL(KFLW,2) C DIMENSION JD6P(KFLV),JD7P(KFLV),JD8P(KFLV),JD9P(KFLV), A JDWORD(6,KFLW),NDBJ(KFLN),NDB6J(KFLN),KD6CP(KFLN), B KD7CP(KFLN),KD8CP(KFLN),KD9CP(KFLN),JDSUM6(KFLS), C JDSUM4(KFLS,KFLW),JDSUM5(KFLS,KFLW),INVD6J(KFLW) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C COMMON /SUMARG/J6P(KFLV),J7P(KFLV),J8P(KFLV),J9P(KFLV), A JWORD(6,KFLW),NLSUM,NBJ(KFLN),NB6J(KFLN),K6CP(KFLN), B K7CP(KFLN),K8CP(KFLN),K9CP(KFLN),JSUM6(KFLS), C JSUM4(KFLS,KFLW),JSUM5(KFLS,KFLW),INV6J(KFLW) C C THIS ROUTINE STORES DATA FOR FUTURE CALLS TO GENSUM C JD6C = J6C JD7C = J7C JD8C = J8C JD9C = J9C JDWC = JWC JDDEL = JDEL MDP = MP NDLSUM = NLSUM IF (J6C.NE.0) THEN DO 10 I = 1,J6C JD6(I) = J6(I) 10 CONTINUE ENDIF C IF (J7C.NE.0) THEN DO 20 I = 1,J7C JD7(I) = J7(I) 20 CONTINUE ENDIF C IF (J8C.NE.0) THEN DO 30 I = 1,J8C JD8(I) = J8(I) 30 CONTINUE ENDIF C IF (J9C.NE.0) THEN DO 40 I = 1,J9C JD9(I) = J9(I) 40 CONTINUE ENDIF C IF (JWC.NE.0) THEN DO 60 I = 1,6 DO 50 J = 1,JWC JDW(I,J) = JW(I,J) 50 CONTINUE 60 CONTINUE DO 70 I = 1,JWC INVD6J(I) = INV6J(I) 70 CONTINUE ENDIF C IF (JDEL.NE.0) THEN DO 90 I = 1,2 DO 80 J = 1,JDEL LDDEL(J,I) = LDEL(J,I) 80 CONTINUE 90 CONTINUE ENDIF C IF (NLSUM.NE.0) THEN DO 110 I = 1,NLSUM NDBJ(I) = NBJ(I) NDB6J(I) = NB6J(I) KD6CP(I) = K6CP(I) KD7CP(I) = K7CP(I) KD8CP(I) = K8CP(I) KD9CP(I) = K9CP(I) 110 CONTINUE ENDIF C DO 120 I = 1,KFLV JD6P(I) = J6P(I) JD7P(I) = J7P(I) JD8P(I) = J8P(I) JD9P(I) = J9P(I) 120 CONTINUE DO 140 I = 1,KFLS JDSUM6(I) = JSUM6(I) DO 130 J = 1,KFLW JDSUM4(I,J) = JSUM4(I,J) JDSUM5(I,J) = JSUM5(I,J) 130 CONTINUE 140 CONTINUE DO 160 I = 1,6 DO 150 J = 1,KFLW JDWORD(I,J) = JWORD(I,J) 150 CONTINUE 160 CONTINUE C END C C C SUBROUTINE LOLPOP(FAIL) C C REDUCES A LOOP WITH ONE LINE AND ONE NODE IN THE FLAT GRAPH. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C DIMENSION KP(3),KS(3) C CHARACTER*6 NAME,NAMSUB COMMON /NAM/NAMSUB DATA NAME/'LOLPOP'/ C DATA KP/2,3,1/ DATA KS/0,1,-1/ C NAMSUB = NAME I1 = NPOINT(1) K3 = 2 IF (I1.EQ.ILAST) K3 = 3 L = JDIAG(I1,K3) CALL DELTA(L,M,FAIL) IF (FAIL) GOTO 20 K = KP(K3) IF (ARR(I1,K).LT.0) THEN C C PHASE2(JDIAG(I1,K)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JDIAG(I1,K) C ENDIF C K1 = KS(K3) IL1 = IL(I1) + K1 I2 = IH(IL1) L1 = JDIAG(I2,1) CALL DELTA(L1,JDIAG(I2,K3),FAIL) IF (FAIL) GOTO 20 IF (ARR(I2,K3).EQ.K1) THEN C C PHASE2(L1): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = L1 C ENDIF C IL2 = IL(I2) + K1 I3 = IH(IL2) K2 = K3 + K1 JDIAG(I3,K2) = L1 ARR(I3,K2) = ARR(I2,1) J9C = J9C + 1 J9(J9C) = L1 J6C = J6C + 1 J6(J6C) = JDIAG(I1,1) IF (K3.EQ.3) GOTO 20 DO 10 I = 3,NBNODE IT = IH(I) ILP = I - 2 IL(IT) = ILP IH(ILP) = IT 10 CONTINUE 20 CONTINUE C END C C C SUBROUTINE MEKEST(KEY,IRHO,ISIG,IRHOP,ISIGP) C*********************************************************************** C C MEKEST C C*********************************************************************** C C KEEP (KEY=2) OR RESTORE (KEY=1) THE COMMON BLOCK MEDEFN, C AND IRHO,ISIG,IRHOP,ISIGP. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) PARAMETER (MXMED=1+5*MXORB2+6*MXORB3) C COMMON /MEDEFN/J(MXMED) C DIMENSION I(MXMED) C SAVE I,I1,I2,I3,I4 C IF (KEY.EQ.2) THEN DO 10 K = 1,MXMED I(K) = J(K) 10 CONTINUE I1 = IRHO I2 = ISIG I3 = IRHOP I4 = ISIGP C ELSE DO 20 K = 1,MXMED J(K) = I(K) 20 CONTINUE IRHO = I1 ISIG = I2 IRHOP = I3 ISIGP = I4 ENDIF C END C C C SUBROUTINE MODJ23(K) C C === MODIFIES THE DIRECT J2 AND J3 ARRAYS FOR EXCHANGE CALL OF NJSYM C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) C PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL23=KFL2*3) C LOGICAL FREE C COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) COMMON /HOLD/J2SPIN(KFL23),J3SPIN(KFL23),J2ANG(KFL23),J3ANG(KFL23) C GOTO (10,40),K C C --- K=1 - SPIN INTEGRALS C 10 CONTINUE MK = M4 I1 = 0 DO 30 J = 1,MK DO 20 K = 1,3 I1 = I1 + 1 J2(J,K) = J2SPIN(I1) J3(J,K) = J3SPIN(I1) 20 CONTINUE 30 CONTINUE J3(1,2) = M14 J3(2,2) = M13 GOTO 70 C C --- K=2 - ANGULAR INTEGRALS C 40 CONTINUE MK = M5 I1 = 0 DO 60 J = 1,MK DO 50 K = 1,3 I1 = I1 + 1 J2(J,K) = J2ANG(I1) J3(J,K) = J3ANG(I1) 50 CONTINUE 60 CONTINUE J2(1,1) = M15 J3(1,3) = M16 70 CONTINUE IF (IBUG1-1) 110,110,80 C C PRINT-OUT OF VALUES IN NJSYM IF IBUG3=1 C 80 CONTINUE IF (IBUG3-1) 90,110,90 90 CONTINUE WRITE (IWRITE,3000) DO 100 J = 1,MK WRITE (IWRITE,3010) (J2(J,K),K=1,3), (J3(J,K),K=1,3) 100 CONTINUE 110 CONTINUE C 3000 FORMAT (' J2',18X,'J3') 3010 FORMAT (3I5,I10,2I5) END C C C SUBROUTINE MUMDAD(II,IJ,IK,M,X) IMPLICIT REAL*8 (A-H,O-Z) C C --- CALLS AND EVALUATES FRACTIONAL PARENTAGE COEFFICIENTS C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) C COMMON /CONSTS/ZERO,TENTH,HALF,ONE,TWO,THREE,FOUR,SEVEN,ELEVEN,EPS COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH(MXORB2,2), A J1QN(MXORB3,3,2),IJFUL(MXORB2) COMMON /INTERM/J1B(MXORB2,3,2),J1T(3,2) C X = ONE LIJ = LJ(IJ) IF (LIJ) 10,10,20 10 CONTINUE IF (M) 60,110,60 20 CONTINUE N = NOSH(IJ,II) IVI = J1QN(IJ,1,II) ILI = (J1QN(IJ,2,II)-1)/2 ISI = J1QN(IJ,3,II) C C IF M=0 THERE ARE QUANTUM NUMBERS WITH TILDES TO CONSIDER C IF (M) 30,40,30 30 CONTINUE IVJ = J1B(IJ,1,II) ILJ = (J1B(IJ,2,II)-1)/2 ISJ = J1B(IJ,3,II) GOTO 50 C 40 CONTINUE IVJ = J1T(1,II) ILJ = (J1T(2,II)-1)/2 ISJ = J1T(3,II) 50 CONTINUE CALL CFP(LIJ,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) X = X*COEFP IF (ABS(X).LT.EPS) GOTO 110 60 CONTINUE LIJ = LJ(IK) IF (LIJ) 110,110,70 70 CONTINUE IF (M) 80,90,80 80 CONTINUE N = NOSH(IK,II) IVI = J1QN(IK,1,II) ILI = (J1QN(IK,2,II)-1)/2 ISI = J1QN(IK,3,II) IVJ = J1B(IK,1,II) ILJ = (J1B(IK,2,II)-1)/2 ISJ = J1B(IK,3,II) GOTO 100 C 90 CONTINUE N = NOSH(IJ,II) - 1 IVI = IVJ ILI = ILJ ISI = ISJ IVJ = J1B(IJ,1,II) ILJ = (J1B(IJ,2,II)-1)/2 ISJ = J1B(IJ,3,II) 100 CONTINUE CALL CFP(LIJ,N,IVI,ILI,ISI,IVJ,ILJ,ISJ,COEFP) X = X*COEFP 110 CONTINUE C END C C C SUBROUTINE NJGRAF(RECUP,FAIL) C*********************************************************************** C C NJGRAF C C*********************************************************************** C*********************************************************************** C C WRITTEN BY: C A. BAR-SHALOM AND M. KLAPISCH C RACAH INSTITUTE OF PHYSICS C HEBREW UNIVERSITY C 91904 JERUSALEM C ISRAEL. C C*********************************************************************** C C-----THIS IS THE MAIN PROGRAM.IT HANDLES ALL THE ANALYSIS OF THE C RECOUPLING COEFFICIENT WITHOUT REFERRING EXPLICITLY TO THE VALUES C OF ANGULAR MOMENTA WHICH ARE IN J1(J),EXCEPT FOR ZERO IN CASE FREE C =.FALSE. .LIKE NJSYM IT PREPARES ARRAYS OF ARGUMENTS FOR PHASE C FACTORS,(2*J+1) FACTORS AND 6J COEFFICIENTS TO BE COMPUTED IN C GENSUM,WHICH CAN ALSO BE CALLED SEPARATELY WHEN ONLY THE NUMERICAL C VALUES OF ANGULAR MOMENTA CHANGE.THESE VARIABLE ANGULAR MOMENTA C SHOULD BE DECLARED FREE(J)=.TRUE.,SO THAT THE FORMULA PREPARED FOR C GENSUM SHOULD BE CORRECT WHEN J1 IS NOT ZERO. C FAIL WILL BE TRUE WHEN THE RECOUPLING COEFFICIENT IS ZERO BECAUSE C OF UNSATISFIED DELTA OR OTHER SIMILAR CAUSES. C C*********************************************************************** C C THIS VERSION HAS BEEN MODIFIED TO HOLD THE ARRAY DIMENSIONS IN C PARAMETER STATEMENTS. THE DIMENSIONS ARE LABELLED: C C KFL1 - DIMENSION OF THE J1 AND FREE ARRAYS IN /COUPLE/, AND THE C FIRST DIMENSION OF THE LINE AND LCOL ARRAYS IN /TREE/. C ALSO THE DIMENSION OF THE SUMVAR ARRAY IN /ARGU/, AND C OF THE INVER ARRAY IN ROUTINE SPRATE. IT IS TESTED FOR M C ON ENTRY TO NJGRAF, AND FOR MP IN ROUTINE SPRATE. C C KFL2 - DIMENSION OF THE J2 AND J3 ARRAYS IN /COUPLE/. THE C DIMENSIONS OF THESE ARRAYS ARE CHECKED ON ENTRY TO NJGRAF. C C KFL2A = (2*KFL2) - DIMENSION OF THE J23, ARROW AND TABS ARRAYS IN C /TREE/. C KFL2B = (4*KFL2) - DIMENSION OF THE JDIAG,ARR,IL AND IH ARRAYS C IN /GRAPH/, AND OF THE IAL ARRAY IN /BUILD/ C KFL2C = (2*KFL2+2) - DIMENSION OF THE NPOINT ARRAY IN /GRAPH/. C C KFL6 - DIMENSION OF THE J6 ARRAY IN /ARGU/. TESTED IN SPRATE. C C KFL7 - DIMENSION OF THE J7 ARRAY IN /ARGU/. TESTED IN SPRATE. C C KFL8 - DIMENSION OF THE J8 ARRAY IN /ARGU/. TESTED IN SPRATE. C C KFL9 - DIMENSION OF THE J9 ARRAY IN /ARGU/. TESTED IN SPRATE. C C KFLW - DIMENSION OF THE JW AND LDEL ARRAYS IN /ARGU/, AND C OF THE JWORD AND INV6J ARRAYS IN /SUMARG/. ALSO THE SECOND C DIMENSION OF THE JSUM4 AND JSUM5 ARRAYS IN /SUMARG/. C IN ADDITION IT GIVES THE DIMENSIONS OF A NUMBER OF C TEMPORARY WORKING ARRAYS IN ROUTINES SPRATE AND GENSUM. C KFLW IS TESTED IN SPRATE. C C KFLS - DIMENSION OF THE JSUM6 ARRAY AND FIRST DIMENSION OF THE C JSUM4 AND JSUM5 ARRAYS IN /SUMARG/. ALSO GIVES THE C DIMENSIONS OF SOME TEMPORARY WORKING ARRAYS IN SPRATE C AND GENSUM. KFLS IS THE MAXIMUM NUMBER OF SUMMATION C VARIABLES IN A PARTICULAR SUM, AND IS TESTED IN SPRATE. C C KFLN - DIMENSION OF THE NBJ, NB6J, K6CP, K7CP, K8CP AND K9CP C ARRAYS IN /SUMARG/. KFLN IS THE MAXIMUM NUMBER OF C SUMS ALLOWED, AND IS TESTED IN ROUTINE SPRATE. C C KFLV - DIMENSION OF THE J6P, J7P, J8P AND J9P ARRAYS IN C /SUMARG/, AND OF THE JNS ARRAY IN ROUTINE VAR. C KFLV IS TESTED IN VAR. C C KFLZ - DIMENSION OF THE JZERO ARRAY IN /ZER/. KFLZ IS TESTED IN C ROUTINE ZERO. C C****************************************************************** C C OTHER CHANGES: C C (1) THIS VERSION HAS BEEN ALTERED TO USE THE RACAH C COEFFICIENT ROUTINE WRITTEN BY STAN SCOTT AND ALAN HIBBERT C (CPC 28 189-200 1982). SUBROUTINE FACTT MUST BE CALLED TO C SET UP THE ARRAY GAM WHICH HOLDS LN(N!) BEFORE NJGRAF IS C CALLED. C (2) THE SUBROUTINES ORDER AND SETDIM HAVE BEEN RENAMED AS C ORDTRI AND SETDM, TO AVOID CONFLICT IN THE QUB CODES. C (3) COMMON BLOCKS /ARGU/ AND /SUMARG/ HAVE BEEN REMOVED FROM C SUBROUTINE GENSUM, AND THEIR CONTENTS ARE TRANSFERRED VIA C THE ARGUMENT LIST INSTEAD. /SUMARG/ HAS THUS BEEN INCLUDED C IN SUBROUTINE NJGRAF TO ACHIEVE THIS. C (4) THE TIMING ROUTINES HAVE BEEN REMOVED FROM THE NJGRAF C PACKAGE AS THEY WERE CDC SPECIFIC AND SEEMED UNNECESSARY. C C**************************************************************** C IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20,KFLS=12, A KFLN=10,KFLV=40) C LOGICAL FAIL,FIND C INTEGER ARROW LOGICAL TABS COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR LOGICAL CUT COMMON /CUTDIG/CUT INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C COMMON /SUMARG/J6P(KFLV),J7P(KFLV),J8P(KFLV),J9P(KFLV), A JWORD(6,KFLW),NLSUM,NBJ(KFLN),NB6J(KFLN),K6CP(KFLN), B K7CP(KFLN),K8CP(KFLN),K9CP(KFLN),JSUM6(KFLS), C JSUM4(KFLS,KFLW),JSUM5(KFLS,KFLW),INV6J(KFLW) C CHARACTER*6 NAME,NAMSUB COMMON /NAM/NAMSUB DATA NAME/'NJGRAF'/ C C TEST THE DIMENSION OF THE J1 ARRAY C IF (M+1.GT.KFL1) THEN PRINT 3000,M + 1,KFL1 STOP C ENDIF C C TEST THE DIMENSIONS OF THE J2 AND J3 ARRAYS C IF (N-1.GT.KFL2) THEN PRINT 3010,N - 1,KFL2 STOP C ENDIF C C 1.0 C BUILDING UP OF THE UNSTRUCTURED GRAPH C FAIL = .FALSE. J6C = 0 J7C = 0 J8C = 0 J9C = 0 JWC = 0 JDEL = 0 CALL SETDM NFIN = 0 CUT = .FALSE. CALL SETTAB(FAIL) M = M + 1 IF (FAIL) GOTO 130 M = M - 1 JF = 0 JF1 = 0 C C 1.1 C LOCATING AND HANDLING OF ZEROS C CALL ZERO(JF1,JF,FAIL) IF (FAIL) GOTO 130 MP = M IF (NBTR.EQ.0) GOTO 140 JUMP = 1 C C 2. C BUILDING OF A FLAT DIAGRAM OUT OF THE UNSTRUCTURED GRAPH. C THERE MAY BE SEVERAL FLAT DIAGRAMS OUT OF THE ORIGINAL C GRAPH,IN CASE OF POSSIBLE CUTS.THEN THE FLAT DIAGRAMS C WILL HAVE FREE ENDS. C 10 CONTINUE CALL DIAGRM(JUMP) NFIN = MAX(0,NFREE-2) IF (NFIN.EQ.0) GOTO 20 JUMP = 3 C C 2.1 C HANDLING OF FREE ENDS IF A CUT WAS FOUND C CALL CUTNL(FAIL) IF (FAIL) GOTO 130 GOTO 50 C 20 CONTINUE JUMP = 2 IF (NFREE-1) 50,30,40 30 CONTINUE CALL CUT1L(FAIL) IF (FAIL) GOTO 130 GOTO 50 C 40 CONTINUE CALL CUT2L(FAIL) IF (FAIL) GOTO 130 50 CONTINUE NBTR = NBTR + NFIN IF (NBTR.NE.0) CUT = .TRUE. C C 3. C ANALYSIS OF THE FLAT DIAGRAM. C CLOSED CIRCUITS OF INCREASING ORDER NC ARE SEARCHED,ANALYSED,AND C TAKEN OUT OF THE FLAT DIAGRAM,THUS REDUCING THE NUMBER OF NODES, C NBNODE. C NC = 0 60 CONTINUE NC = NC + 1 CALL SEARCH(FIND) IF (.NOT.FIND) GOTO 60 NCP = NC - 2 JPOL = 0 IF (M.EQ.MP .AND. NC.GT.3) CALL SETDM IF (IPARTL.GT.2) CALL POLYGN(JPOL) GOTO (70,80,90,100),NC C 70 CONTINUE CALL LOLPOP(FAIL) IF (FAIL) GOTO 130 GOTO 110 C 80 CONTINUE CALL BUBBLE(JPOL,FAIL) IF (FAIL) GOTO 130 GOTO 110 C 90 CONTINUE CALL TRIANG(FAIL) IF (FAIL) GOTO 130 GOTO 110 C 100 CONTINUE CALL SQUARE 110 CONTINUE NBNODE = NBNODE - 2 IF (NBNODE.EQ.0) GOTO 120 IFIRST = IH(1) ILAST = IH(NBNODE) C PRINTJ IS AN ALL PURPOSE PRINTING SUBROUTINE CALLED FROM MANY PLA CALL PRINTJ(NAMSUB,8) IF (NBNODE.EQ.NFIN) GOTO 120 NC = NCP C C PROCEED TO OTHER CIRCUITS OF ORDER NC-1 C GOTO 60 C 120 CONTINUE IF (NBTR.EQ.0) GOTO 140 IF (JUMP.EQ.3) CALL ORDTRI C C AT THIS STAGE,THE FLAT DIAGRAM HAS BEEN REDUCED TO NODES C INVOLVING FREE ENDS.PROCEED TO BUILD OTHER FLAT DIAGRAMS C IF NECESSARY. C GOTO 10 C C ALL PARTS OF THE ORIGINAL GRAPH HAVE BEEN REDUCED. C 130 CONTINUE RECUP = 0.0D0 M = M - 1 RETURN C 140 CONTINUE CALL PRINTJ(NAME,0) C C 4.PREPARATION OF THE RESULTS,AND SEPARATION IN SEVERAL SUMS C IF CUTS HAVE BEEN DETECTED,ALSO IN THE FLAT DIAGRAM ITSELF C CALL SPRATE(M) M = M - 1 C C 5. GENSUM COMPUTES THE NUMERICAL VALUE OF THE RECOUPLING COEFFICIEN C CALL GENSUM(J6C,J7C,J8C,J9C,JWC,J6,J7,J8,J9,JW,JDEL,LDEL, A MP,J6P,J7P,J8P,J9P,JWORD,NLSUM,NBJ,NB6J,K6CP,K7CP, B K8CP,K9CP,JSUM4,JSUM5,JSUM6,INV6J,RECUP) C 3000 FORMAT (1X,'DIMENSION ERROR IN NJGRAF. M+1=',I5,' KFL1=',I5) 3010 FORMAT (' DIMENSION ERROR IN NJGRAF. N-1=',I5,' KFL2=',I5) END C C C SUBROUTINE NORM(N,X) C C NORMALISES THE VECTOR X OF DIMENSION N C SUCH THAT THE LARGEST COMPONENT IS UNITY. C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION X(N) C G = ABS(X(1)) DO 10 I = 2,N GA = ABS(X(I)) IF (GA.GT.G) G = GA 10 CONTINUE DO 20 I = 1,N X(I) = X(I)/G 20 CONTINUE C END C C C FUNCTION NTAB1(NELCTS,K) C*********************************************************************** C C NTAB1 C C*********************************************************************** C C THIS SUBROUTINE CALCULATES THE ROW OF NTAB CORRESPONDING TO THE C PARENTS WHICH MAY GIVE RISE TO THE TERM ASSOCIATED WITH SHELL C LAMBDA . E.G. IF WE SEEK THE ROW OF NTAB CONTAINING THE PARENTS C OF ONE OF THE P**3 TERMS, THE ROW = VALUE OF NTAB1 IS THAT C CONTAINING THE P**2 TERMS C C USE IS MADE OF THE FACT THAT THE LIST OF POSSIBLE PARENTS (SEE C WHITE - ATOMIC SPECTRA - APPENDIX) IS SYMMETRICAL ABOUT THE C CONFIGURATION L**(2L+1) C C C --- FOR ONE ELECTRON IN A TERM, THE PARENT IS ALWAYS A SINGLET S TERM C COMMON /INFORM/IREAD,IWRITE,IPUNCH C IF (NELCTS-1) 20,10,20 10 CONTINUE NTAB1 = 2 RETURN C C OTHERWISE THE VALUE OF NTAB1 DEPENDS ON THE L VALUE OF THE C ELECTRONS C 20 CONTINUE GOTO (30,40,70,100,250,120,140,160,180, A 200,220),K C C --- FOR S ELECTRONS, THE ONLY OTHER POSSIBILITY IS THAT NELCTS=2 C 30 CONTINUE NTAB1 = 1 RETURN C C --- P ELECTRONS - ARE WE BEYOND P**4 C 40 CONTINUE IF (NELCTS-4) 60,60,50 50 CONTINUE NELCTS = 8 - NELCTS 60 CONTINUE NTAB1 = 1 + NELCTS RETURN C C --- D ELECTRONS - ARE WE BEYOND D**6 C 70 CONTINUE IF (NELCTS-6) 90,90,80 80 CONTINUE NELCTS = 12 - NELCTS 90 CONTINUE NTAB1 = 4 + NELCTS RETURN C C --- F ELECTRONS - ARE THERE MORE THAN TWO. IF SO, THE PROGRAMME NEEDS C AN F-SHELL COEFFICIENT-OF-FRACTIONAL-PARENTAGE ROUTINE, AND THE C ARRAYS IN /TERMS/ NEED EXTENDING C 100 CONTINUE IF (NELCTS-2) 10,110,270 110 CONTINUE NTAB1 = 11 RETURN C C C L=3 ELECTRONS SHOULD NOT BE MORE THAN TWO C C 120 CONTINUE IF (NELCTS-2) 10,130,240 130 CONTINUE NTAB1 = 13 RETURN C 140 CONTINUE IF (NELCTS-2) 10,150,240 150 CONTINUE NTAB1 = 14 RETURN C 160 CONTINUE IF (NELCTS-2) 10,170,240 170 CONTINUE NTAB1 = 15 RETURN C 180 CONTINUE IF (NELCTS-2) 10,190,240 190 CONTINUE NTAB1 = 16 RETURN C 200 CONTINUE IF (NELCTS-2) 10,210,240 210 CONTINUE NTAB1 = 17 RETURN C 220 CONTINUE IF (NELCTS-2) 10,230,240 230 CONTINUE NTAB1 = 18 240 CONTINUE RETURN C C --- G ELECTRONS - ARE THERE MORE THAN TWO. IF SO, THE PROGRAMME C NEEDS A G-SHELL COEFFICIENT-OF-FRACTIONAL-PARENTAGE ROUTINE, AND C THE ARRAYS IN /TERMS/ NEED EXTENDING C 250 CONTINUE IF (NELCTS-2) 10,260,280 260 CONTINUE NTAB1 = 12 RETURN C 270 CONTINUE WRITE (IWRITE,3000) GOTO 290 C 280 CONTINUE WRITE (IWRITE,3010) 290 CONTINUE STOP C 3000 FORMAT (//// A' STOP AND EXTEND THE NTAB AND ITAB ARRAYS TO ALLOW CORE F-ELECTRO BNS'/ C' YOU WILL ALSO REQUIRE A COMPLETE FRACTIONAL PARENTAGE ROUTINE FO DR F-ELECTRONS'//) 3010 FORMAT (//// A' STOP AND EXTEND THE NTAB AND ITAB ARRAYS TO ALLOW CORE G-ELECTRO BNS'/ C' YOU WILL ALSO REQUIRE A COMPLETE FRACTIONAL PARENTAGE ROUTINE FO DR G-ELECTRONS'//) END SUBROUTINE ORDTRI C C THIS SUBROUTINE ORDERS THE TRIADS WHICH WERE LEFT WITH FREE ENDS C AS CONSEQUENCE OF CUTTING,SO THAT THE NEW GRAPH WILL START THERE. C C Corrections by SMW Bailey, DAMTP QUB, 29 Apr 1998. The value of C NBT1 has been changed at two places to put a triad in the correct C position, and a small loop has been moved to avoid overwriting a C triad at position NM. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP INTEGER ARROW LOGICAL TABS COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C COMMON /BUILD/IAL(KFL2B),IF1,IF2,NODE C C COMMON /KEEP/JKP(2,3),JARR(2,3),IT2,IT3,IT5 C CHARACTER*6 NAME DATA NAME/'ORDTRI'/ C DO 10 I = 1,MP IAL(I) = 0 10 CONTINUE IF (NFIN.NE.0) GOTO 70 NF = NBTR - ITFREE(1) IF (IT5.NE.0) GOTO 20 NBT1 = NBTR - 1 N0 = 0 NFT = NFREE ASSIGN 60 TO ISW GOTO 160 C 20 CONTINUE NFT = IT5 - IT2 NM = NFT + NBTR + 1 NBT1 = NBTR N0 = 0 ASSIGN 40 TO ISW GOTO 160 C 40 CONTINUE N0 = NFT NBT1 = NBT1 + 1 NFT = IT3 - IT5 ASSIGN 50 TO ISW GOTO 160 C 50 CONTINUE NBT1 = K - NFT DO 30 J = 1,3 JDIAG(NBTR,J) = JKP(1,J) ARR(NBTR,J) = JARR(1,J) JDIAG(NM,J) = JKP(2,J) ARR(NM,J) = JARR(2,J) 30 CONTINUE C 60 CONTINUE NODE = NBT1 + NFT CALL CHANGE(NODE,2) GOTO 130 C 70 CONTINUE NBT1 = NBTR - 1 NBT = NBT1 + NFIN NBTT = NBT + 1 NB = 0 80 CONTINUE DO 120 I = 1,NBNODE I1 = IH(I) IF (IL(I1).GT.ILAST) GOTO 120 I2 = NBT1 + I IF (I1.GT.NBTT) GOTO 90 IF (I1.EQ.I2) GOTO 110 IF (IL(I2).LE.NBNODE) GOTO 120 90 CONTINUE DO 100 J = 1,3 JDIAG(I2,J) = JDIAG(I1,J) ARR(I2,J) = ARR(I1,J) 100 CONTINUE IL(I1) = ILAST + I 110 CONTINUE NB = NB + 1 IL(I2) = 0 120 CONTINUE IF (NB.NE.NFIN) GOTO 80 NODE = NBT 130 CONTINUE IF1 = JDIAG(NBTR,1) IF2 = JDIAG(NBTR,3) DO 150 I = NBTR,NODE DO 140 K = 1,3 J = JDIAG(I,K) IAL(J) = IAL(J) + 1 140 CONTINUE 150 CONTINUE ILAST = NODE CALL PRINTJ(NAME,8) RETURN C 160 CONTINUE IF (NF.GT.0) GOTO 170 NFR = N0 I1 = 1 GOTO 180 C 170 CONTINUE NFR = NFT + 1 I1 = -1 180 CONTINUE DO 200 I = 1,NFT IK = NFR + I1*I IT = ITFREE(IK) K = NBT1 + IK DO 190 J = 1,3 JDIAG(K,J) = JDIAG(IT,J) ARR(K,J) = ARR(IT,J) 190 CONTINUE 200 CONTINUE GOTO ISW(40,60,50) C END C C C SUBROUTINE ORTHOG(LET) C*********************************************************************** C C ORTHOG C C*********************************************************************** C C THIS SUBROUTINE CHECKS FOR POSSIBLE ORTHOGONALITY DUE TO C COUPLING DIFFERENCES OR UNEVEN PARITY C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) C COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) C C --- DO PSI AND PSIP CONTAIN THE SAME NUMBERS OF ELECTRONS C DO PSI AND PSIP HAVE THE SAME TOTAL PARITY C N5 = 0 N6 = 0 N7 = 0 DO 10 I = 1,IHSH L1 = LJ(I) L2 = NOSH1(I) L3 = NOSH2(I) N5 = N5 + L2 N6 = N6 + L3 N7 = N7 + L1* (L2-L3) 10 CONTINUE C C CHECK ON NUMBER OF ELECTRONS C IF (N5-N6) 20,40,20 20 CONTINUE IF (IBUG1-1) 110,30,30 30 CONTINUE WRITE (IWRITE,3020) GOTO 110 C C CHECK ON PARITY C 40 CONTINUE IF (N7-N7/2*2) 50,70,50 50 CONTINUE IF (IBUG1-1) 110,60,60 60 CONTINUE WRITE (IWRITE,3030) GOTO 110 C 70 CONTINUE N1 = 2*IHSH - 1 C C --- IS THE FINAL STATE THE SAME FOR PSI AND PSIP C DO 80 K = 2,3 IF (J1QN1(N1,K)-J1QN2(N1,K)) 90,80,90 80 CONTINUE GOTO 120 C 90 CONTINUE IF (IBUG1-1) 110,100,100 100 CONTINUE WRITE (IWRITE,3000) C C --- THE TWO CONFIGURATIONS WILL HAVE ZERO HAMILTONIAN MATRIX ELEMENT C 110 CONTINUE LET = 0 RETURN C 120 CONTINUE C C --- NO OBVIOUS ANGULAR MOMENTUM ORTHOGONALITY C LET = 1 C 3000 FORMAT (' DIFFERING RESULTANT ANGULAR MOMENTUM') 3020 FORMAT ( A ' THE TWO CONFIGURATIONS HAVE DIFFERING NUMBERS OF NS)???ONS' B ) 3030 FORMAT (' THE TWO CONFIGURATIONS HAVE DIFFERING TOTAL PARITY') END C C C SUBROUTINE POLYGN(JPOL) C C THIS ROUTINE REDUCES A CIRCUIT OF ARBITRARY ORDER NC.IT EXCHANGES C NODES ON THE FLAT DIAGRAM UNTIL THE DISTANCE ON THE AXIS BETWEEN C NODES EQUEALS ONE.EACH EXCHANGE INTRODUCES A SUMMATION VARIABLE C AND A 6J SYMBOL.THE CIRCUIT HAS A MAXIMUM OF NPART=2 DISCONNECTED C PARTS ON THE AXIS. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C CHARACTER*6 NAME DATA NAME/'POLYGN'/ C NC1 = NC + 1 NC2 = NC NBC = IPARTL - 2 10 CONTINUE DO 90 I = 1,NBC IT2 = NPOINT(NC1-I) IT1 = NPOINT(NC2-I) JB = JDIAG(IT1,1) JC = JDIAG(IT2,1) JDIAG(IT1,1) = JC JDIAG(IT2,1) = JB JAR = ARR(IT1,1) ARR(IT1,1) = ARR(IT2,1) ARR(IT2,1) = JAR JE = JDIAG(IT1,2) MP = MP + 1 SUMVAR(MP) = .TRUE. JDIAG(IT1,2) = MP JDIAG(IT2,3) = MP IF (TAB1(JB,1)-IT1) 30,20,30 20 CONTINUE TAB1(JB,1) = IT2 GOTO 40 C 30 CONTINUE TAB1(JB,2) = IT2 40 CONTINUE IF (TAB1(JC,1)-IT2) 60,50,60 50 CONTINUE TAB1(JC,1) = IT1 GOTO 70 C 60 CONTINUE TAB1(JC,2) = IT1 70 CONTINUE IF (ARR(IT1,2).GT.0) GOTO 80 C C PHASE2(JE): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JE C ARR(IT1,2) = 1 ARR(IT2,3) = -1 80 CONTINUE JWC = JWC + 1 JW(1,JWC) = JB JW(2,JWC) = MP JW(3,JWC) = JE JW(4,JWC) = JC JW(5,JWC) = JDIAG(IT2,2) JW(6,JWC) = JDIAG(IT1,3) J6(J6C+1) = MP J6C = J6C + 2 J6(J6C) = MP 90 CONTINUE NC = NC - NBC IF (NC.LE.4) GOTO 100 NBC = IPARTS - 2 NC1 = IPARTS + 1 NC2 = IPARTS GOTO 10 C 100 CONTINUE IF (NPART.EQ.1) GOTO 110 NPOINT(3) = NPOINT(NC1) NPOINT(4) = NPOINT(NC1+1) 110 CONTINUE IF (NC.EQ.2) JPOL = 1 CALL PRINTJ(NAME,10) C END C C C SUBROUTINE PRINTJ(NAMES,JP) C C THIS SUBROUTINE PRINTS INTERMEDIATE RESULTS IN STANDARD FORM FROM C WHEREVER IT IS CALLED. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2,KFLZ=99) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C CHARACTER IM,IP,IS(3) CHARACTER*4 I6,I7,I8,I9,IJ1 CHARACTER*6 NAMES,NSETTB CHARACTER*8 IBLANK,IFREE,IFR C DIMENSION IX(7),JTAB(KFL1,3) INTEGER ARROW LOGICAL TABS COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP INTEGER ARR,TAB1 LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC C OUT COMMON/CONST/I6C,I7C,I8C,I9C,IDEL,IWC -- OUT WE'90MAR16: EQUIV++ COMMON /ZER/NZERO,JZERO(KFLZ) COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 C EQUIVALENCE (I6C,IX(1)), (I7C,IX(2)), (I8C,IX(3)), (I9C,IX(4)), A (IDEL,IX(5)), (IWC,IX(6)) C DATA IBLANK,IFREE,IP,IM/' ','FREE END','+','-'/ DATA NSETTB/'SETTAB'/ DATA I6,I7,I8,I9,IJ1/'I6= ','I7= ','I8= ','I9= ','J1= '/ DATA I6C,I7C,I8C,I9C,IDEL,IWC/6*1/ C IF (IBUG3.NE.1) RETURN PRINT 3170,NAMES JUMP = JP IF (JUMP.NE.0) GOTO 20 DO 10 I = 1,7 IX(I) = 1 10 CONTINUE PRINT 3120,IJ1, (J1(I),I=1,M) 20 CONTINUE IF (JUMP.LT.8) GOTO 90 PRINT 3000,NBNODE,NBTR,NFIN,IFIRST,ILAST,NFREE JUMP = JUMP - 8 PRINT 3010 K = 0 DO 50 I = 1,NBNODE IT = IH(I) IFR = IBLANK JT = JDIAG(IT,1) IF (TAB1(JT,2).EQ.IT .AND. JT.NE.JDIAG(IFIRST,3)) GOTO 30 K = K + 1 JTAB(K,1) = JT JTAB(K,2) = TAB1(JT,1) JTAB(K,3) = TAB1(JT,2) 30 CONTINUE IF (TAB1(JT,2).GT.ILAST) IFR = IFREE DO 40 J = 1,3 IS(J) = IP IF (ARR(IT,J).LT.1) IS(J) = IM 40 CONTINUE PRINT 3020, (IS(J),J=1,3) PRINT 3030,IL(IT),IT,IFR, (JDIAG(IT,J),J=1,3) 50 CONTINUE PRINT 3040 NTIME = 0 JT = JDIAG(IFIRST,3) IF (JT.EQ.JDIAG(ILAST,2)) GOTO 60 IF (TAB1(JT,2).GE.1000) GOTO 60 GOTO 70 C 60 CONTINUE K = K + 1 JTAB(K,1) = JT JTAB(K,2) = TAB1(JT,1) JTAB(K,3) = TAB1(JT,2) 70 CONTINUE NTIME = NTIME + 1 IF (NTIME.EQ.2) GOTO 80 JT = JDIAG(ILAST,2) IF (TAB1(JT,2).EQ.1000) GOTO 60 80 CONTINUE PRINT 3050, ((JTAB(I,J),J=1,3),I=1,K) PRINT 3060, (I,SUMVAR(I),I=1,MP) 90 CONTINUE IF (JUMP.LT.4) GOTO 120 JUMP = JUMP - 4 NBTR1 = 2*N - 2 PRINT 3070,NBTR1 K = 0 DO 110 I = 1,NBTR1 IF (TABS(I)) GOTO 110 K = K + 1 DO 100 J = 1,3 IS(J) = IP IF (ARROW(I,J).LT.1) IS(J) = IM 100 CONTINUE PRINT 3080, (IS(J),J=1,3) PRINT 3090,K,I, (J23(I,J),J=1,3) 110 CONTINUE PRINT 3100 MM = M IF (NAMES.NE.NSETTB) MM = M - 1 PRINT 3110, (I, (LINE(I,J),LCOL(I,J),J=1,2),I=1,MM) 120 CONTINUE IF (JUMP.LT.2) GOTO 130 JUMP = JUMP - 2 PRINT 3150,NC,NPART,IPARTL,IPARTS,ICROSS, (NPOINT(I),I=1,NC) 130 CONTINUE IF (JUMP.LT.1) GOTO 140 PRINT 3160,NZERO, (I,JZERO(I),I=1,NZERO) 140 CONTINUE IF (J6C.GE.I6C) PRINT 3120,I6, (J6(I),I=I6C,J6C) IF (J7C.GE.I7C) PRINT 3120,I7, (J7(I),I=I7C,J7C) IF (J8C.GE.I8C) PRINT 3120,I8, (J8(I),I=I8C,J8C) IF (J9C.GE.I9C) PRINT 3120,I9, (J9(I),I=I9C,J9C) IF (JDEL.GE.IDEL) PRINT 3130, ((LDEL(I,J),J=1,2),I=IDEL,JDEL) IF (JWC.GE.IWC) PRINT 3140, ((JW(J,I),J=1,6),I=IWC,JWC) I6C = J6C + 1 I7C = J7C + 1 I8C = J8C + 1 I9C = J9C + 1 IDEL = JDEL + 1 IWC = JWC + 1 C 3000 FORMAT (/10X,'NBNODE=',I3,10X,'NBTR=',I3,10X,'NFIN=',I3,/10X, A 'IFIRST=',I3,10X,'ILAST=',I3,9X,'NFREE=',I3) 3010 FORMAT (//7X,'IL',3X,'IH',14X,'JDIAG'//) 3020 FORMAT (26X,3 (2X,A1)) 3030 FORMAT (7X,I2,3X,I2,2X,A8,2X,3I3/) 3040 FORMAT (/5X,'TAB1'/) 3050 FORMAT (4 (I3,')',2X,I3,I5,5X)) 3060 FORMAT (/2X,'SUMVAR=',15 (I3,L1)) 3070 FORMAT (//10X,'J23',10X,'NBTR1=',I3//) 3080 FORMAT (16X,3 (2X,A1)) 3090 FORMAT (I9,I5,2X,3I3/) 3100 FORMAT (/3X,'J L1 K1 L2 K2') 3110 FORMAT (4 (I4,')',I3,I3,I4,I3)) 3120 FORMAT (/3X,A4,3X,3 (20I3/)) 3130 FORMAT (/3X,'DELTA=',7 (I5,I3)) 3140 FORMAT (/3X,'JW(ARG. OF 6J)',6I3) 3150 FORMAT (//' NC=',I2,4X,'NPART=',I2,4X,'IPARTL=',I2,4X,'IPARTS=', A I2,4X,'ICROSS=',I2,4X,/2X,'NPOINT=',20I3) 3160 FORMAT (//' NZERO=',I2,5X,12 (I4,')',I3)) 3170 FORMAT (///3X,'PRINT OUT AFTER CALLING SUBROUTINE ',A7) END C C C SUBROUTINE REDUCE(IRHO,ISIG,IRHOP,ISIGP,LESSEN) C*********************************************************************** C C REDUCE C C*********************************************************************** C C THIS SUBROUTINE REMOVES SPECTATOR SINGLET S SHELLS WHICH HAVE C NO EFFECT IN ANGULAR OR SPIN INTEGRALS C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) DIMENSION LEAVE(MXORB2) COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) C C LMIN INITIALLY SET LARGE C LMIN = 99 ICOUNT = 0 DO 30 I = 1,IHSH C C NO INTERACTING SHELL MAY BE REMOVED C IF (I.EQ.IRHO .OR. I.EQ.ISIG .OR. I.EQ.IRHOP .OR. A I.EQ.ISIGP) GOTO 10 C C IF A SPECTATOR SHELL HAS SINGLET S COUPLING ON BOTH SIDES OF C THE MATRIX ELEMENT, IT MAY, IN GENERAL, BE REMOVED, AS IT HAS NO C EFFECT IN FANO C IF (J1QN1(I,1).EQ.0 .AND. J1QN2(I,1).EQ.0) GOTO 20 10 CONTINUE ICOUNT = ICOUNT + 1 LEAVE(ICOUNT) = I GOTO 30 C 20 CONTINUE IF (LJ(I).GE.LMIN) GOTO 30 LMIN = LJ(I) ILMIN = I 30 CONTINUE IF (ICOUNT.EQ.IHSH) GOTO 120 C C IF A CHANGE IN THE COMMON BLOCK MEDEFN IS TO BE MADE, C ITS PRESENT SITUATION MUST BE PRESERVED BY A CALL OF MEKEEP C CALL MEKEST(2,IRHO,ISIG,IRHOP,ISIGP) C C IF ONLY ONE SHELL WOULD BE LEFT IN THIS WAY, THE ONE, DESTINED C FOR REMOVAL, WITH THE LOWEST L-VALUE MUST BE RETAINED TO DEFINE A C COUPLING C IF (ICOUNT.EQ.1) GOTO 90 C C --- MODIFY THE COMMON BLOCK MEDEFN C 40 CONTINUE DO 60 I = 1,ICOUNT J = LEAVE(I) IF (J.EQ.IRHO) IRHO = I IF (J.EQ.ISIG) ISIG = I IF (J.EQ.IRHOP) IRHOP = I IF (J.EQ.ISIGP) ISIGP = I NJ(I) = NJ(J) LJ(I) = LJ(J) NOSH1(I) = NOSH1(J) NOSH2(I) = NOSH2(J) IJFUL(I)=IJFUL(J) DO 50 K = 1,3 J1QN1(I,K) = J1QN1(J,K) J1QN2(I,K) = J1QN2(J,K) 50 CONTINUE 60 CONTINUE ISUBH = IHSH - 1 DO 80 I = 2,ICOUNT J = LEAVE(I) II = ICOUNT + I - 1 IJ = ISUBH + J DO 70 K = 1,3 J1QN1(II,K) = J1QN1(IJ,K) J1QN2(II,K) = J1QN2(IJ,K) 70 CONTINUE 80 CONTINUE IHSH = ICOUNT GOTO 100 C C THIS SITUATION ONLY OCCURS IF IRHO=ISIG=IRHOP=ISIGP C 90 CONTINUE J = LEAVE(1) II = MIN(J,ILMIN) LEAVE(1) = II LEAVE(2) = J + ILMIN - II ICOUNT = 2 GOTO 40 C 100 CONTINUE IF (IBUG1.LE.1) GOTO 110 I2HSH = IHSH + IHSH - 1 WRITE (IWRITE,3000) WRITE (IWRITE,3020) ((J1QN1(J,K),K=1,3),J=1,I2HSH) WRITE (IWRITE,3010) WRITE (IWRITE,3020) ((J1QN2(J,K),K=1,3),J=1,I2HSH) C C LESSEN = 0 IF NO CHANGE IN MEDEFN C = 1 OTHERWISE C 110 CONTINUE LESSEN = 1 RETURN C 120 CONTINUE LESSEN = 0 C 3000 FORMAT (/' NEW DEFINITION OF COUPLING SCHEMES'/ A ' FOR THIS SET OF RHO, SIG, RHOP, SIGP'//10X, B ' L.H.S. OF HAMILTONIAN MATRIX ELEMENT DEFINED BY') 3010 FORMAT (10X,' R.H.S. OF HAMILTONIAN MATRIX ELEMENT DEFINED BY') 3020 FORMAT (10X,' J1QN ', (I5,2I3)) END C C C FUNCTION RME(L,LP,K) C*********************************************************************** C C RME C C*********************************************************************** C C --- EVALUATES THE REDUCED MATRIX ELEMENT (L//C(K)//LP) - SEE FANO C AND RACAH, IRREDUCIBLE TENSORIAL SETS, CHAP. 14, P. 81 C IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' PARAMETER (MXFCT=500) C COMMON /CONSTS/ZERO,TENTH,HALF,ONE,TWO,THREE,FOUR,SEVEN,ELEVEN,EPS COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /FACT/GAMMA(MZFAC) COMMON /FACTS/GAM(MXFCT) C IF (K.GT. (L+LP) .OR. K.LT.ABS(L-LP)) GOTO 20 I2G = L + LP + K IG = I2G/2 IF (I2G-2*IG) 10,40,10 10 CONTINUE RME = ZERO RETURN C 20 CONTINUE IF (IBUG1-1) 10,10,30 30 CONTINUE WRITE (IWRITE,3000) L,LP,K GOTO 10 C 40 CONTINUE IF (IG) 20,50,60 50 CONTINUE RME = ONE RETURN C 60 CONTINUE I1 = IG - L I2 = IG - LP I3 = IG - K C NRB IMX=MAX(I1,I2,I3,IG+1) IF(2*IMX+1.GT.MZFAC)THEN AL1 = GAM(I1+1) AL2 = GAM(2*I1+1) ALP1 = GAM(I2+1) ALP2 = GAM(2*I2+1) AK1 = GAM(I3+1) AK2 = GAM(2*I3+1) QUSQRT = (2*L+1)*(2*LP+1)*EXP(AL2+ALP2+AK2-GAM(I2G+2)) RME = SQRT(QUSQRT)*EXP(GAM(IG+1)-AL1-ALP1-AK1) ELSE C NRB AL1 = GAMMA(I1+1) AL2 = GAMMA(2*I1+1) ALP1 = GAMMA(I2+1) ALP2 = GAMMA(2*I2+1) AK1 = GAMMA(I3+1) AK2 = GAMMA(2*I3+1) QUSQRT = (2*L+1)* (2*LP+1)*AL2*ALP2*AK2/GAMMA(I2G+2) RME = SQRT(QUSQRT)*GAMMA(IG+1)/ (AL1*ALP1*AK1) ENDIF C 3000 FORMAT (//' L =',I3,' LP =',I3,' K =',I3, A ' RME SET ZERO SINCE ???ANGLE DOES NOT MATCH') END C C C SUBROUTINE SEARCH(FIND) C LOGICAL FIND C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4,KFL2A=2*KFL2,KFL2B=4*KFL2, A KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC COMMON /INFORM/IREAD,IWRITE,IPUNCH C CHARACTER*6 NAME DATA NAME/'SEARCH'/ C C-----THIS SUBROUTINE LOCATES CIRCUITS OR LOOPS OF ORDER NC.NPOINT(NC) C ARE THE INDICES OF THE POINTS(TRIADS) PERTAINING TO THE FIRST C SUCH LOOP FOUND. C NPART IS THE NUMBER OF SEPARATE PARTS(GROUPS OF CONTIGUOUS POINTS) C ON THE AXIS OF THE FLAT GRAPH.IPARTS IS THE NUMBER OF POINTS IN C THE SMALLEST PART.IPARTL IS THE NUMBER OF POINTS IN THE LARGEST C PART. C THIS SUBROUTINE FINDS ALL THE POSSIBLE LOOPS OF ORDER 3 AND 4.FOR C NC.GE.5,IT LOOKS FOR ONLY THOSE WHO ARE PARTITIONNED IN NPART.LE.2 C , WHICH CAN EVENTUALLY C REDUCE TO A LOOP OF ORDER 4 WITHOUT BREAKING THE BASIC STRUCTURE C OF THE FLAT GRAPH.ICROSS=-1,IF LINES CROSS C-------------------------------------------------------------------- C C INITIALIZATION C FIND = .FALSE. NCM1 = NC - 1 NCM = NC - 2 ICROSS = 0 C C FIRST TREAT TWO CASES THAT DO NOT INVOLVE DO LOOPS: C C 1. ONE ISOLATED POINT,EITHER THE FIRST OR THE LAST C NPART = 1 IPARTL = NC - 1 IPARTS = 1 C C A. FIRST C I1 = IFIRST K3 = 3 K2 = 2 10 CONTINUE JA = JDIAG(I1,1) JC = JDIAG(I1,K3) C IF (JA.EQ.JC) THEN IF (NC.GT.1) THEN WRITE (IWRITE,3000) I1,K3,JA,JC,NC STOP C ENDIF C NPOINT(1) = I1 GOTO 160 C ENDIF C I2 = TAB1(JA,K2) I3 = TAB1(JC,K2) C IF (ABS(IL(I3)-IL(I2))-NCM.LT.0) THEN WRITE (IWRITE,3010) I2,I3,JA,JC,K2,NC STOP C ENDIF C IF (ABS(IL(I3)-IL(I2))-NCM.GT.0) THEN C C B. LAST C IF (I1.NE.IFIRST) GOTO 30 I1 = ILAST K3 = 2 K2 = 1 GOTO 10 C ENDIF C IC = 1 NPOINT(IC) = I1 I20 = MIN(I2,I3) I21 = IL(I20) I31 = I21 + NCM1 C DO 20 II = I21,I31 IC = IC + 1 NPOINT(IC) = IH(II) 20 CONTINUE C IF (NC.LE.2) THEN IF (JDIAG(IFIRST,1).NE.JDIAG(ILAST,1)) THEN C C PHASE(I1,JDIAG,KFL2B): C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C J7(J7C+1) = JDIAG(I1,1) J7(J7C+2) = JDIAG(I1,2) J7C = J7C + 3 J7(J7C) = JDIAG(I1,3) C GOTO 160 C ENDIF C ENDIF C IF (I1.NE.ILAST) THEN IT = I2 JT = JDIAG(ILAST,2) K4 = 2 I4 = ILAST C ELSE IT = I3 JT = JDIAG(IFIRST,3) K4 = 3 I4 = IFIRST ENDIF C IF (IT.EQ.I20) THEN C C PHASE (I1,JDIAG,KFL2B): C PHASE FACTOR ARISING FROM NON-CYCLIC PERMUTATION OF ARGUMENTS IN C TRIAD L. C J7(J7C+1) = JDIAG(I1,1) J7(J7C+2) = JDIAG(I1,2) J7C = J7C + 3 J7(J7C) = JDIAG(I1,3) C ENDIF C IF ((JT.EQ.JA) .OR. (JT.EQ.JC)) CALL CHANGE(I4,K4) GOTO 160 C C 2. TWO ISOLATED POINTS,FIRST AND LAST C 30 CONTINUE IF (NC.EQ.1) RETURN IF (NC.LE.3) GOTO 50 IPARTL = NC - 2 IPARTS = 1 I1 = IFIRST I2 = ILAST JA = JDIAG(I1,1) JB = JDIAG(I1,3) C IF (TAB1(JA,2).NE.I2) THEN JA = JDIAG(I1,3) JB = JDIAG(I1,1) IF (TAB1(JA,2).NE.I2) GOTO 50 ENDIF C IF (JA.EQ.JDIAG(I2,1)) THEN JC = JDIAG(I2,2) C ELSE JC = JDIAG(ILAST,1) ENDIF C I3 = TAB1(JB,2) I4 = TAB1(JC,1) IDIST = IL(I4) - IL(I3) C IF (ABS(IDIST)- (NCM-1).LT.0) THEN WRITE (IWRITE,3020) I3,I4,JB,JC,IDIST,NC STOP C ENDIF C IF (ABS(IDIST)- (NCM-1).EQ.0) THEN NPOINT(1) = ILAST NPOINT(2) = IFIRST ICROSS = SIGN(1,IDIST) IC = 2 I20 = MIN(I3,I4) I21 = IL(I20) I31 = I21 + NCM C DO 40 II = I21,I31 IC = IC + 1 NPOINT(IC) = IH(II) 40 CONTINUE C IF (JA.EQ.JDIAG(IFIRST,1)) CALL CHANGE(IFIRST,3) IF (JA.EQ.JDIAG(ILAST,1)) CALL CHANGE(ILAST,2) GOTO 160 C ENDIF C C FIRST GENERAL CASE: ALL POINTS IN ONE GROUP C 50 CONTINUE NPART = 1 IPARTS = 0 IPARTL = NC K3 = 1 C DO 80 IN = 1,NBNODE I = IH(IN) 60 CONTINUE JA = JDIAG(I,K3) IF (I.NE.TAB1(JA,2)) THEN I2 = TAB1(JA,2) C IF (IL(I2)-IN-NCM1.LT.0) THEN WRITE (IWRITE,3030) IN,I,I2,IL(I2),JA,NC STOP C ENDIF C IF (IL(I2)-IN-NCM1.EQ.0) THEN I21 = IL(I2) IC = 0 C DO 70 II = IN,I21 IC = IC + 1 NPOINT(IC) = IH(II) 70 CONTINUE C IF (JA.EQ.JDIAG(IFIRST,3)) CALL CHANGE(IFIRST,3) IF (JA.EQ.JDIAG(ILAST,2)) CALL CHANGE(ILAST,2) GOTO 160 C ENDIF C ENDIF C IF (IN.EQ.1) THEN IF (K3.NE.3) THEN K3 = 3 GOTO 60 C ELSE K3 = 1 ENDIF C ENDIF C 80 CONTINUE C C SEARCH DID NOT FIND LOOP NC .LE. 3 C IF (NC.LE.3) RETURN C C GENERAL CASE OF LOOP PARTITIONNED IN 2 GROUPS. DO LOOP C ON IPARTS C NPART = 2 NC2 = NC/2 K3 = 1 K2 = 1 C DO 150 IPS = 2,NC2 JPS = IPS - 1 NBN = NBNODE - JPS C DO 140 I1 = 1,NBN I = IH(I1) I2 = IH(I1+JPS) 90 CONTINUE JA = JDIAG(I,K3) JD = JDIAG(I2,K2) C IF (I.EQ.TAB1(JA,1)) THEN II2 = TAB1(JD,2) II1 = TAB1(JA,2) C ELSE II1 = TAB1(JA,1) II2 = TAB1(JD,1) ENDIF C IDIST = IL(II1) - IL(II2) C IF (ABS(IDIST)- (NCM-JPS).LT.0) THEN WRITE (IWRITE,3040) JPS,I1,I,I2,JA,JD,II1,II2,IDIST,NC STOP C ENDIF C IF (ABS(IDIST)- (NCM-JPS).GT.0) GOTO 130 ICROSS = SIGN(1,IDIST) IC = 0 I21 = IL(I2) C DO 110 II = I1,I21 IC = IC + 1 NPOINT(IC) = IH(II) 110 CONTINUE C I20 = MIN(II1,II2) I30 = MAX(II1,II2) I21 = IL(I20) I31 = IL(I30) C DO 120 II = I21,I31 IC = IC + 1 NPOINT(IC) = IH(II) 120 CONTINUE C IPARTS = IPS IPARTL = NC - IPS IF ((JDIAG(IFIRST,3).EQ.JA) .OR. A (JDIAG(IFIRST,3).EQ.JD)) CALL CHANGE(IFIRST,3) IF ((JDIAG(ILAST,2).EQ.JA) .OR. A (JDIAG(ILAST,2).EQ.JD)) CALL CHANGE(ILAST,2) GOTO 160 C 130 CONTINUE IF (I1.EQ.1) THEN IF (K3.EQ.3) THEN K3 = 1 GOTO 140 C ELSE K3 = 3 GOTO 90 C ENDIF C ENDIF C IF (I2.EQ.ILAST) THEN IF (K2.NE.2) THEN K2 = 2 GOTO 90 C ENDIF C ENDIF C 140 CONTINUE 150 CONTINUE C C SEARCH DID NOT FIND CIRCUIT OF ORDER NC C RETURN C C LOOP FOUND C 160 CONTINUE FIND = .TRUE. CALL PRINTJ(NAME,10) C C ERROR PRINTOUT C 3000 FORMAT (' ERROR IN SEARCH. I1,K3,JA,JC,NC = ',5I5) 3010 FORMAT (' ERROR IN SEARCH. I2,I3,JA,JC,K2,NC = ',6I5) 3020 FORMAT (' ERROR IN SEARCH. I3,I4,JB,JC,IDIST,NC = ',6I5) 3030 FORMAT (' ERROR IN SEARCH. IN,I,I2,IL(I2),JA,NC = ',6I5) 3040 FORMAT (' ERROR IN SEARCH. JPS,I1,I,I2,JA,JD,II1,II2,IDIST,NC = ', A 10I5) END C C C SUBROUTINE SETDM C C SET DIMENSIONS OF ARRAYS. C PARAMETER (KFL1=100) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C COMMON /DIM/J6CC,J7CC,J8CC,J9CC,JWCC,JDELC C JWCC = JWC JDELC = JDEL J6CC = J6C J7CC = J7C J8CC = J8C J9CC = J9C C END C C C SUBROUTINE SETJ1(K,IRHO,ISIG,IRHOP,ISIGP,ITST1,ITST2,K1,K2,K3,K4) C C === SETS J1 AND FREE ARRAYS FOR DIRECT INTEGRAL CALLS OF NJGRAF C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) PARAMETER (KFL1=100,KFL2=MXORB+4) C LOGICAL FREE C COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 COMMON /NJLJ/NRHO,LRHO,NSIG,LSIG,NRHOP,LRHOP,NSIGP,LSIGP COMMON /INTERM/J1BAR1(MXORB2,3),J1BAR2(MXORB2,3),J1TLD1(3), A J1TLD2(3) COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) C DO 10 J = 1,IHSH J1(J) = J1BAR1(J,K) 10 CONTINUE DO 20 J = M4,M6 J1(J) = J1QN1(J,K) 20 CONTINUE DO 30 J = M7,M8 J1(J) = J1QN2(J-M3,K) 30 CONTINUE J1(M10) = J1QN1(ISIG,K) J1(M12) = J1QN2(ISIGP,K) IF (M1) 40,50,40 40 CONTINUE J1(M9) = J1QN1(IRHO,K) GOTO 60 C 50 CONTINUE J1(M9) = J1TLD1(K) 60 CONTINUE IF (M2) 70,80,70 70 CONTINUE J1(M11) = J1QN2(IRHOP,K) GOTO 90 C 80 CONTINUE J1(M11) = J1TLD2(K) C C K=2 IMPLIES ANGULAR TERMS , K=3 IMPLIES SPIN TERMS C 90 CONTINUE IF (K-2) 110,110,100 100 CONTINUE J1(M13) = 2 J1(M14) = 2 MLIMIT = M14 NJ1S = M14 NJ23S = M5 GOTO 120 C 110 CONTINUE J1(M13) = 2*LRHO + 1 J1(M14) = 2*LSIG + 1 J1(M15) = 2*LRHOP + 1 J1(M16) = 2*LSIGP + 1 MLIMIT = M16 NJ1S = M17 NJ23S = M18 C C PRINT-OUT OF VALUES IN NJGRAF IF IBUG3=1 C 120 CONTINUE IF (IBUG1.GT.1 .AND. IBUG3.NE.1) WRITE (IWRITE,3000) (J1(J),J=1, A MLIMIT) C C IF ITST1.NE.2 OR ITST2.NE.2 THEN NJGRAF IS BEING CALLED, SO SET C THE ELEMENTS OF THE FREE ARRAY. C C IF ((ITST1.NE.2) .OR. (ITST2.NE.2)) THEN C DO 130 J = 1,MLIMIT FREE(J) = .FALSE. 130 CONTINUE IF (K1.NE.1) FREE(IRHO) = .TRUE. IF (M1.EQ.0) THEN IF (K2.NE.1) FREE(M9) = .TRUE. C ELSE IF (K2.NE.1) FREE(ISIG) = .TRUE. ENDIF C IF (M2.EQ.0 .AND. K4.NE.1) FREE(M11) = .TRUE. C C PRINT-OUT OF VALUES IN NJGRAF IF IBUG3=1 C IF (IBUG1.GT.1 .AND. IBUG3.NE.1) WRITE (IWRITE,3010) (FREE(J), A J=1,MLIMIT) C ENDIF C C 3000 FORMAT (' J1: ',36I3) 3010 FORMAT ('FREE: ',36L3) END C C C SUBROUTINE SETM C*********************************************************************** C C SETM C C*********************************************************************** C C --- SET CONSTANTS USEFUL IN INNER SUBROUTINES C IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) C COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) COMMON /MVALUE/M1,M2,M3,M4,M5,M6,M7,M8,M9,M10,M11,M12,M13,M14,M15, A M16,M17,M18,M19,M20 C M3 = IHSH - 1 M4 = IHSH + 1 M5 = IHSH + 2 M6 = 2*IHSH - 1 M7 = M6 + 1 M8 = M3 + M6 M9 = M8 + 1 M10 = M8 + 2 M11 = M8 + 3 M12 = M8 + 4 M13 = M8 + 5 M14 = M8 + 6 M15 = M8 + 7 M16 = M8 + 8 M17 = M8 + 9 M18 = IHSH + 3 C END C C C SUBROUTINE SETTAB(FAIL) C C BUILDS UP THE UNSTRUCTURED GRAPH C SETS THE ARRAY J23,CONTAINING THE TWO LISTS OF ORIGINAL TRIADS C J2 AND J3,AND THE CORRESPONDING ARROWS ON THE ANGULAR MOMENTA C LINES.ALSO ESTABLISHES THE NUMERICAL AND PHASE FACTORS CONNECTING C RECOUPLING COEFFICIENT AND GRAPHS,ACCORDING TO YUTSIS,LEVINSON AND C VANAGAS.FOR THIS PURPOSE DETERMINES THE TOTAL J C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C LOGICAL TABS INTEGER ARROW COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR C LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C COMMON /BUILD/IAL(KFL2B),IF1,IF2,NODE C CHARACTER*6 NAME DATA NAME/'SETTAB'/ C IPR = N - 1 NBTR = IPR + IPR DO 20 I = 1,IPR DO 10 J = 1,2 J23(I,J) = J2(I,J) ARROW(I,J) = 1 10 CONTINUE TABS(I) = .FALSE. J23(I,3) = J2(I,3) ARROW(I,3) = -1 20 CONTINUE IPR1 = IPR + 1 DO 40 I = IPR1,NBTR II = I - IPR DO 30 J = 1,2 J23(I,J) = J3(II,J) ARROW(I,J) = -1 30 CONTINUE TABS(I) = .FALSE. J23(I,3) = J3(II,3) ARROW(I,3) = 1 40 CONTINUE DO 50 J = 1,NBTR J8(J) = J23(J,1) 50 CONTINUE J8C = NBTR + IPR NB1 = NBTR + 1 DO 60 J = NB1,J8C I = J - IPR J8(J) = J23(I,3) 60 CONTINUE J6C = NBTR DO 70 J = 1,J6C J6(J) = J23(J,3) 70 CONTINUE DO 80 I = 1,M SUMVAR(I) = .FALSE. IAL(I) = 1 80 CONTINUE DO 100 I = 1,NBTR DO 90 J = 1,3 JI = J23(I,J) K = IAL(JI) LINE(JI,K) = I LCOL(JI,K) = J IAL(JI) = K + 1 90 CONTINUE 100 CONTINUE IT = 0 DO 130 I = 1,NBTR JT = J23(I,3) IF (IAL(JT).EQ.3) GOTO 120 IF (IT.EQ.1) GOTO 110 JT1 = JT IT = 1 GOTO 130 C 110 CONTINUE CALL DELTA(JT1,JT,FAIL) IF (FAIL) GOTO 150 K = LINE(JT,1) KC = LCOL(JT,1) LINE(JT1,2) = K LCOL(JT1,2) = KC J23(K,KC) = JT1 IAL(JT) = 1 GOTO 140 C C OTHERJ(I,JT,L,LC,K): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C 120 CONTINUE L = LINE(JT,1) IF (L.EQ.I .OR. TABS(L)) THEN K = 1 L = LINE(JT,2) LC = LCOL(JT,2) C ELSE K = 2 LC = LCOL(JT,1) ENDIF C IF (LC.EQ.3) GOTO 140 130 CONTINUE 140 CONTINUE J9(J9C+1) = JT J9C = J9C + 2 J9(J9C) = JT 150 CONTINUE CALL PRINTJ(NAME,4) C END C C C SUBROUTINE SETUPE(JA,JB,NJCOMP,LJCOMP) C*********************************************************************** C C SETUPE C C*********************************************************************** C C === GENERATES THE ARRAYS NJ,LJ - DEFINING THE QUANTUM NUMBERS OF THE C SHELLS, NOSH - DEFINING THE OCCUPATION OF THE SHELLS, J1QN - C DEFINING THE COUPLING OF THE SHELLS, FOR EACH OF THE TWO C CONFIGURATIONS CONSIDERED. ONLY THOSE SHELLS OCCURRING IN AT C LEAST ONE CONFIGURATION ARE INCLUDED. C AT LEAST TWO SHELLS MUST BE CONSIDERED OCCUPIED. C THUS (1S)**2 HELIUM MUST BE TREATED AS ,E.G., (1S)**2(2S)**0 C THE SIZE OF THE ARRAYS HERE CALCULATED IS ARRANGED TO BE NO C GREATER THAN IS NECESSARY TO INCLUDE ALL ORBITALS WHICH ARE C DEEMED TO BE OCCUPIED IN EITHER OR BOTH OF THE CONFIGURATIONS C JA,JB C INCLUDE 'PARAM' C PARAMETER (MXOCC1=MZOCC+1,MXOCC2=2*MZOCC+1) PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) C DIMENSION NJCOMP(MXORB2),LJCOMP(MXORB2) C COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH(MXORB2,2), A J1QN(MXORB3,3,2),IJFUL(MXORB2) COMMON /MSTATE/MCFG,MOCCSH(MZNC2),MOCORB(MXOCC1,MZNC2), A MELCSH(MXOCC1,MZNC2),M1QNRD(MXOCC2,3,MZNC2),KCFG, B KOCCSH(MZNC2),KOCORB(MXOCC1,MZNC2),KELCSH(MXOCC1,MZNC2), C K1QNRD(MXOCC2,3,MZNC2),MAXOR C C NOTICE THE DIFFERENT NAMES IN THE COMMON BLOCK MEDEFN - WE C STORE NOSH1(I=1,10) AS NOSH((I=1,10),1) AND NOSH2(I=1,10) AS C NOSH((I=1,10),2) AND USE THE FACT THAT NOSH1 AND NOSH2 WILL THEN C BE EQUIVALENT TO THE SINGLE 2-DIMENSIONAL ARRAY NOSH. SIMILARLY C FOR J1QN C C --- INITIALIZE BASIC QUANTITIES - (I1+1) RUNS OVER 1,MAXORB, IHSH IS C THE CURRENT VALUE OF THE HIGHEST OCCUPIED SHELL YET CONSIDERED, C WHILE I2HSH=2*IHSH-1 C I1 = 0 IHSH = 0 I2HSH = -1 IA = MOCCSH(JA) IB = KOCCSH(JB) C C --- TEST ON WHETHER LIMIT OF I1 HAS BEEN REACHED C 10 CONTINUE IF (I1-MAXOR) 20,350,350 C C --- INCREASE BASIC QUANTITIES C 20 CONTINUE I1 = I1 + 1 I3 = IHSH + 1 I5 = I2HSH + I3 C C --- IS THE SHELL I1 OCCUPIED IN JA C DO 30 J = 1,IA IF (I1-MOCORB(J,JA)) 30,40,30 30 CONTINUE NA = 1 GOTO 50 C 40 CONTINUE NA = 2 J1 = J C C --- IS THE SHELL I1 OCCUPIED IN JB C 50 CONTINUE DO 60 J = 1,IB IF (I1-KOCORB(J,JB)) 60,70,60 60 CONTINUE NB = 1 GOTO 80 C 70 CONTINUE NB = 2 J2 = J C C IF THE SHELL I1 IS NOT OCCUPIED IN EITHER JA OR JB, IGNORE THE C SHELL, DO NOT INCREASE IHSH, AND CONSIDER NEXT SHELL BY INCREASING C I1 C 80 CONTINUE IF (NA-1) 90,90,100 90 CONTINUE IF (NB-1) 10,10,100 C C --- IF THE SHELL I1 IS OCCUPIED IN EITHER JA OR JB - C (1) IF IHSH.GT.1, THEN ALREADY AT LEAST TWO SHELLS AND THE C RESULTING COUPLINGS HAVE BEEN STORED. WE MUST THUS MAKE ROOM FOR C THE QUANTUM NUMBERS OF THIS NEW SHELL BETWEEN THE QUANTUM NUMBERS C OF THE PREVIOUS SHELLS AND THE QUANTUM NUMBERS OF THE INTERMEDIATE C COUPLINGS OF THE CONFIGURATIONS. THUS THE LATTER SET ARE =MOVED C ALONG= TO MAKE ROOM FOR THE NEW SHELL C (2) IF IHSH.LE.1, THERE ARE NO INTERMEDIATE COUPLING QUANTUM C NUMBERS, AND SO THERE IS NOTHING TO MOVE C 100 CONTINUE IF (IHSH-1) 150,150,110 110 CONTINUE DO 140 I = 1,2 DO 130 J = I3,I2HSH I4 = I5 - J DO 120 K = 1,3 J1QN(I4+1,K,I) = J1QN(I4,K,I) 120 CONTINUE 130 CONTINUE 140 CONTINUE 150 CONTINUE IHSH = I3 I2HSH = I2HSH + 2 NC = NA I = 1 IC = J1 JC = JA C C --- FIRST CONSIDER THE L.H.S. (I=1) OF THE MATRIX ELEMENT. NC=1 MEANS C UNOCCUPIED, REPRESENTED BY A DUMMY SINGLET S SHELL, AND THE C ADDITIONAL SET OF COUPLING QUANTUM NUMBERS WILL BE THE SAME AS THE C LAST SET OF COUPLING QUANTUM NUMBERS ALREADY OBTAINED. C NC=2 MEANS OCCUPIED. THEN ALL THE NEW QUANTUM NUMBERS (BOTH FOR C THE SHELL AND FOR THE COUPLING OF THIS SHELL TO THE RESULTANT OF C THE PREVIOUS ONES) ARE DEFINED IN THE CORRESPONDING J1QNRD ARRAY. C NOSH - THE NUMBER OF ELECTRONS IN THIS SHELL, IS DEFINED BY THE C APPROPRIATE ENTRY IN NELCSH . THE R.H.S. IS THEN CONSIDERED C SIMILARLY (I=2) C 160 CONTINUE GOTO (170,210),NC C 170 CONTINUE NOSH(IHSH,I) = 0 J1QN(IHSH,1,I) = 0 J1QN(IHSH,2,I) = 1 J1QN(IHSH,3,I) = 1 IF (IHSH-2) 320,180,190 180 CONTINUE J1QN(3,1,I) = 0 J1QN(3,2,I) = J1QN(1,2,I) J1QN(3,3,I) = J1QN(1,3,I) GOTO 320 C 190 CONTINUE DO 200 K = 1,3 J1QN(I2HSH,K,I) = J1QN(I2HSH-1,K,I) 200 CONTINUE GOTO 320 C 210 CONTINUE IF (I.GE.2) GOTO 280 NOSH(IHSH,I) = MELCSH(IC,JC) DO 270 K = 1,3 J1QN(IHSH,K,I) = M1QNRD(IC,K,JC) C C IS THIS THE FIRST OCCUPIED SHELL OF EITHER CONFIGURATION. IF SO, C THEN THERE ARE NO INTERMEDIATE COUPLINGS TO CONSIDER AT THIS STAGE C IF (IHSH-1) 270,270,220 C C IS THIS THE FIRST OCCUPIED SHELL OF THIS CONFIGURATION, THOUGH NOT C THE FIRST OF THE OTHER CONFIGURATION. IF SO, THE INTERMED9ATE C COUPLING FORMED HAS THE SAME L,S VALUES AS THIS OCCUPIED SHELL, C SINCE WE COUPLE THE SHELL TO A DUMMY SINGLET S. C 220 CONTINUE IF (IC-1) 230,230,250 230 CONTINUE I2 = 1 IF (K-1) 260,240,260 C C SENIORITY SET (ARBITRARILY) ZERO FOR INTERMEDIATE COUPLING C 240 CONTINUE J1QN(I2HSH,1,I) = 0 GOTO 270 C 250 CONTINUE I2 = MOCCSH(JC) + IC - 1 260 CONTINUE J1QN(I2HSH,K,I) = M1QNRD(I2,K,JC) 270 CONTINUE GOTO 320 C 280 CONTINUE NOSH(IHSH,I) = KELCSH(ICE,JCE) DO 310 K = 1,3 J1QN(IHSH,K,I) = K1QNRD(ICE,K,JCE) IF (IHSH.LE.1) GOTO 310 C C IS THIS THE FIRST OCCUPIED SHELL OF THIS CONFIGURATION, THOUGH NOT C THE FIRST OF THE OTHER CONFIGURATION. IF SO, THE INTERMEDIATE C COUPLING FORMED HAS THE SAME L,S VALUES AS THIS OCCUPIED SHELL, C SINCE WE COUPLE THE SHELL TO A DUMMY SINGLET S. C IF (ICE.GT.1) GOTO 290 I2 = 1 IF (K.NE.1) GOTO 300 J1QN(I2HSH,1,I) = 0 GOTO 310 C 290 CONTINUE I2 = KOCCSH(JCE) + ICE - 1 300 CONTINUE J1QN(I2HSH,K,I) = K1QNRD(I2,K,JCE) 310 CONTINUE 320 CONTINUE IF (I-2) 330,340,340 330 CONTINUE NC = NB I = 2 ICE = J2 JCE = JB GOTO 160 C C --- SET THE NJ AND LJ VALUES OF THE OCCUPIED SHELLS C 340 CONTINUE NJ(IHSH) = NJCOMP(I1) LJ(IHSH) = LJCOMP(I1) IJFUL(IHSH)=I1 C C --- RETURN TO 1 TO SEE IF MAXORB HAS BEEN REACHED C GOTO 10 C 350 CONTINUE C END C C C SUBROUTINE SPRATE(M) C C THIS SUBROUTINE PREPARES THE INFORMATION TO BE TRANSFERED TO C GENSUM FOR NUMERICAL EVALUATION. C PARAMETER (KFL1=100) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20,KFLS=12, A KFLN=10,KFLV=40) C LOGICAL SUM6J,T6J,JT,JS DIMENSION JTEM4(KFLS,KFLW),JTEM5(KFLS,KFLW),JTEM6(KFLS), A NSUM6J(KFLW),J6SUM(KFLW) DIMENSION SUM6J(KFLW),T6J(KFLW),JT(KFLS),JS(KFLS),INVER(KFL1), A JNSUM(KFLS),JINV(KFLS),N6JN(KFLW),IN6J(KFLW), B JSUMT(KFLW,6) C LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C COMMON /DIM/J6CC,J7CC,J8CC,J9CC,JWCC,JDELC C COMMON /SUMARG/J6P(KFLV),J7P(KFLV),J8P(KFLV),J9P(KFLV), A JWORD(6,KFLW),NLSUM,NBJ(KFLN),NB6J(KFLN),K6CP(KFLN), B K7CP(KFLN),K8CP(KFLN),K9CP(KFLN),JSUM6(KFLS), C JSUM4(KFLS,KFLW),JSUM5(KFLS,KFLW),INV6J(KFLW) C LOGICAL CUT COMMON /CUTDIG/CUT C CHARACTER*4 NAME C C C 1. TEST THAT ARRAY DIMENSIONS HAVE NOT BEEN EXCEEDED. C IF (MP.LE.KFL1) GOTO 10 NMX = KFL1 NPX = MP NAME = 'KFL1' GOTO 60 C 10 CONTINUE IF (JWC.LE.KFLW) GOTO 20 NMX = KFLW NPX = JWC NAME = 'KFLW' GOTO 60 C 20 CONTINUE IF (J6C.LE.KFL6) GOTO 30 NMX = KFL6 NPX = J6C NAME = 'KFL6' GOTO 60 C 30 CONTINUE IF (J7C.LE.KFL7) GOTO 40 NMX = KFL7 NPX = J7C NAME = 'KFL7' GOTO 60 C 40 CONTINUE IF (J8C.LE.KFL8) GOTO 50 NMX = KFL8 NPX = J8C NAME = 'KFL8' GOTO 60 C 50 CONTINUE IF (J9C.LE.KFL9) GOTO 70 NMX = KFL9 NPX = J9C NAME = 'KFL9' 60 CONTINUE PRINT 3000,NAME,NPX,NMX STOP C C 2. DETERMINATION OF EFFECTIVE SUMMATION VARIABLES AND THEIR C RELATIONSHIPS WITH 6J COEFFICIENTS. C 70 CONTINUE DO 80 I = 1,JWC INV6J(I) = 0 SUM6J(I) = .FALSE. 80 CONTINUE NSUM = 0 NLSUM = 0 IF (MP.EQ.M) RETURN M1 = M + 1 DO 90 I = M1,MP IF (.NOT.SUMVAR(I)) GOTO 90 NSUM = NSUM + 1 JSUM6(NSUM) = 0 INVER(I) = NSUM 90 CONTINUE IF (NSUM.EQ.0) RETURN IF (NSUM.LE.KFLS) GOTO 100 NMX = KFLS NPX = NSUM NAME = 'NSUM' GOTO 60 C 100 CONTINUE KT = 0 DO 130 I = 1,JWC DO 120 J = 1,6 IK = JW(J,I) IF (.NOT.SUMVAR(IK)) GOTO 120 IF (SUM6J(I)) GOTO 110 SUM6J(I) = .TRUE. KT = KT + 1 J6SUM(KT) = 0 NSUM6J(KT) = I INV6J(I) = KT 110 CONTINUE ISK = INVER(IK) I2 = JSUM6(ISK) + 1 JSUM6(ISK) = I2 JSUM4(ISK,I2) = J JSUM5(ISK,I2) = KT I3 = J6SUM(KT) + 1 J6SUM(KT) = I3 JSUMT(KT,I3) = ISK 120 CONTINUE 130 CONTINUE CALL VAR(J6,J6P,J6C,J6CP,J6CC,SUMVAR,MP,M,INVER) CALL VAR(J7,J7P,J7C,J7CP,J7CC,SUMVAR,MP,M,INVER) CALL VAR(J8,J8P,J8C,J8CP,J8CC,SUMVAR,MP,M,INVER) CALL VAR(J9,J9P,J9C,J9CP,J9CC,SUMVAR,MP,M,INVER) IF (CUT) GOTO 180 NLSUM = 1 NBJ(1) = NSUM NB6J(1) = KT K6CP(1) = J6CP K7CP(1) = J7CP K8CP(1) = J8CP K9CP(1) = J9CP DO 150 I = 1,KT I1 = NSUM6J(I) DO 140 J = 1,6 JWORD(J,I) = JW(J,I1) 140 CONTINUE 150 CONTINUE DO 170 I = 1,NSUM ISU = JSUM6(I) DO 160 J = 1,ISU I1 = JSUM5(I,J) J1 = JSUM4(I,J) JWORD(J1,I1) = MP + I 160 CONTINUE 170 CONTINUE GOTO 410 C C 3.SEPARATION OF VARIABLES AND SUMS IN CASE A CUT WAS DETECTED. C 180 CONTINUE K6C = 0 K7C = 0 K8C = 0 K9C = 0 NJ = 0 N6J = 0 DO 190 I = 1,KT T6J(I) = .FALSE. 190 CONTINUE DO 200 I = 1,NSUM JT(I) = .FALSE. JS(I) = .FALSE. 200 CONTINUE J = 1 210 CONTINUE NJ = NJ + 1 JNSUM(NJ) = J JINV(J) = NJ JT(J) = .TRUE. 220 CONTINUE JS(J) = .TRUE. JS6 = JSUM6(J) DO 250 I = 1,JS6 I6J = JSUM5(J,I) IF (T6J(I6J)) GOTO 230 T6J(I6J) = .TRUE. N6J = N6J + 1 N6JN(N6J) = NSUM6J(I6J) IN6J(I6J) = N6J 230 CONTINUE J6J = J6SUM(I6J) DO 240 K = 1,J6J JK = JSUMT(I6J,K) IF (JT(JK)) GOTO 240 NJ = NJ + 1 JNSUM(NJ) = JK JINV(JK) = NJ JT(JK) = .TRUE. 240 CONTINUE 250 CONTINUE DO 260 JJ = 1,NSUM J = JJ IF (JS(JJ) .OR. .NOT.JT(JJ)) GOTO 260 GOTO 220 C 260 CONTINUE NLSUM = NLSUM + 1 IF (NLSUM.LE.KFLN) GOTO 280 NMX = KFLN NPX = NLSUM NAME = 'KFLN' GOTO 60 C 280 CONTINUE NBJ(NLSUM) = NJ NB6J(NLSUM) = N6J IF (J6CP.EQ.0) GOTO 290 CALL CHVAR(J6P,J6CP,K6C,JT,JINV,NSUM) 290 CONTINUE K6CP(NLSUM) = K6C IF (J7CP.EQ.0) GOTO 300 CALL CHVAR(J7P,J7CP,K7C,JT,JINV,NSUM) 300 CONTINUE K7CP(NLSUM) = K7C IF (J8CP.EQ.0) GOTO 310 CALL CHVAR(J8P,J8CP,K8C,JT,JINV,NSUM) 310 CONTINUE K8CP(NLSUM) = K8C IF (J9CP.EQ.0) GOTO 320 CALL CHVAR(J9P,J9CP,K9C,JT,JINV,NSUM) 320 CONTINUE K9CP(NLSUM) = K9C IF (NJ.EQ.NSUM) GOTO 340 DO 330 JJ = 1,NSUM J = JJ IF (.NOT.JT(JJ)) GOTO 210 330 CONTINUE 340 CONTINUE DO 360 I = 1,KT I1 = N6JN(I) DO 350 J = 1,6 JWORD(J,I) = JW(J,I1) 350 CONTINUE 360 CONTINUE DO 380 I = 1,NSUM IK = JNSUM(I) I2 = JSUM6(IK) JTEM6(I) = I2 DO 370 J = 1,I2 JTEM4(I,J) = JSUM4(IK,J) K = JSUM5(IK,J) JTEM5(I,J) = IN6J(K) 370 CONTINUE 380 CONTINUE DO 400 I = 1,NSUM I2 = JTEM6(I) JSUM6(I) = I2 DO 390 J = 1,I2 I1 = JTEM5(I,J) J1 = JTEM4(I,J) JSUM4(I,J) = J1 JSUM5(I,J) = I1 JWORD(J1,I1) = I + MP 390 CONTINUE 400 CONTINUE 410 CONTINUE C 3000 FORMAT (' DIMENSION ERROR FOR ',A4,I5, A ' IS OUT OF ALLOWED RANGE',I3) END C C C SUBROUTINE SQUARE C C REDUCES A CIRCUIT OF ORDER 4 IN THE TWO CASES WHICH ARE LEFT C OVER BY POLYGN,NAMELY TWO DISCONNECTED GROUPS OF TWO POINTS C AND ONE GROUP OF TWO POINTS PLUS THE TWO ENDS OF THE AXIS.IN C THE LATTER, THE END OF THE AXIS IS TRANSFERRED TO THE BEGINNING. C IN THIS PROCESS,ONE SUMMATION VARIABLE AND TWO 6J SYMBOLS ARE C INTRODUCED. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C CHARACTER*6 NAME,NAMSUB COMMON /NAM/NAMSUB DATA NAME/'SQUARE'/ C NAMSUB = NAME MP = MP + 1 SUMVAR(MP) = .TRUE. IT1 = NPOINT(1) IT2 = NPOINT(2) C IF (ICROSS.EQ.1) THEN IT3 = NPOINT(3) IT4 = NPOINT(4) K23 = 3 K32 = 2 C ELSE IT3 = NPOINT(4) IT4 = NPOINT(3) K23 = 2 K32 = 3 ENDIF C L4 = JDIAG(IT2,1) C IF (ARR(IT2,1).LE.0) THEN C C PHASE2 (L4): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = L4 C ARR(IT2,1) = 1 ARR(IT3,1) = -1 ENDIF C L2 = JDIAG(IT1,1) IF (ARR(IT1,1).GT.0) THEN C C PHASE2(L2): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = L2 C ENDIF C JWC = JWC + 1 JW(1,JWC) = L4 JW(2,JWC) = L2 JW(3,JWC) = JDIAG(IT2,2) JJ1 = JDIAG(IT1,3) JW(4,JWC) = JJ1 JW(5,JWC) = MP JW(6,JWC) = JDIAG(IT1,2) IF (ARR(IT1,2).LT.0) THEN C C PHASE2 (JDIAG(IT1,2)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JDIAG(IT1,2) C ENDIF C JWC = JWC + 1 JW(1,JWC) = L4 JW(2,JWC) = L2 JJ3 = JDIAG(IT3,K23) JJ2 = JDIAG(IT4,K32) JW(3,JWC) = JJ3 JW(4,JWC) = JJ2 JW(5,JWC) = MP JW(6,JWC) = JDIAG(IT3,K32) IF (ARR(IT3,K32).LT.0) THEN C C PHASE2 (JDIAG(IT3,K32)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JDIAG(IT3,K32) C ENDIF C J6(J6C+1) = MP J6C = J6C + 2 J6(J6C) = MP C IF (NPART.EQ.1) THEN ITMIN = IT2 ITMAX = IT3 C ELSE ITMIN = MIN(IT2,IT3) ITMAX = MAX(IT2,IT3) ENDIF C ITMN = MIN(IT1,IT4) ITMX = MAX(IT1,IT4) C TAB1(MP,1) = ITMIN TAB1(MP,2) = ITMAX JDIAG(IT2,1) = MP JDIAG(IT3,1) = MP JDIAG(IT2,3) = JJ1 ARR(IT2,3) = ARR(IT1,3) JDIAG(IT3,K32) = JJ2 ARR(IT3,K32) = ARR(IT4,K32) C IF (ICROSS.EQ.1) THEN J7(J7C+1) = L2 J7(J7C+2) = L4 C C PHASE2 (L4): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = L4 C J7C = J7C + 3 J7(J7C) = MP C ELSE C C PHASE2 (JJ2): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JJ2 C ENDIF C ITLL = IL(ITMN) ITHL = IL(ITMX) C DO 10 I = ITLL + 1,ITHL - 1 IT = IH(I) ILP = I - 1 IL(IT) = ILP IH(ILP) = IT 10 CONTINUE IF (ITHL.NE.NBNODE) THEN DO 20 I = ITHL + 1,NBNODE IT = IH(I) ILP = I - 2 IL(IT) = ILP IH(ILP) = IT 20 CONTINUE ENDIF C IF (NPART.NE.2) THEN TAB1(JJ1,1) = IH(1) TAB1(JJ1,2) = IH(NBNODE-2) ENDIF C END C C C SUBROUTINE TENSOR(KA,ISPIN,IRHO,ISIG,VSHELL) C*********************************************************************** C C TENSOR C C*********************************************************************** C C REDUCED MATRIX ELEMENT OF SUMMATIONS OF ONE PARTICLE TENSOR C OPERATORS BY W D ROBB CATALOGUE NUMBER AAKF (CPC 6 (1973) 132) C AND CORRECTION DECK CATALOGUE NUMBER AAKF000A C C*********************************************************************** C C W. D. ROBB - NOVEMBER 1971 C MODIFIED TO USE THE NJGRAF PACKAGE--DECEMBER 86 C C*********************************************************************** C C A ROUTINE FOR THE EVALUATION OF ANGULAR AND SPIN FACTORS IN THE C REDUCED MATRIX ELEMENT OF ANY ONE-ELECTRON TENSOR OPERATOR BETWEEN C ARBITRARILY COUPLED L-S CONFIGURATIONS C C*********************************************************************** C C ** NOTE THAT THE DEFINITIONS OF TENSOR OPERATORS USED ARE THOSE C OF FANO AND RACAH, IRREDUCIBLE TENSORIAL SETS, ACADEMIC PRESS 1959 C C*********************************************************************** IMPLICIT REAL*8 (A-H,O-Z) C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (MXORB2=MXORB+2,MXORB3=2*MXORB+3) C PARAMETER (KFL1=100,KFL2=MXORB+4) C PARAMETER (ZERO=0.0D0,ONE=1.0D0, A HALF=0.5D0) C DIMENSION J2STO(KFL2,3),J3STO(KFL2,3),JMEM(16),VSHELL(MXORB2) C LOGICAL FAIL,FREE C COMMON /COUPLE/NJ1S,NJ23S,J1(KFL1),J2(KFL2,3),J3(KFL2,3), A FREE(KFL1) COMMON /INFORM/IREAD,IWRITE,IPUNCH COMMON /MEDEFN/IHSH,NJ(MXORB2),LJ(MXORB2),NOSH1(MXORB2), A NOSH2(MXORB2),J1QN1(MXORB3,3),J1QN2(MXORB3,3),IJFUL(MXORB2) COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /TERMS/NROWS,ITAB(18),JTAB(18),NTAB(189) C NBUG6 = IBUG6 AJF = ONE RML = ZERO RPL = ZERO NTOT = 0 DO 10 IS = 1,IHSH VSHELL(IS) = ZERO 10 CONTINUE IHSHP1 = IHSH + 1 I2HSH = IHSH*2 - 1 C C PRINT OUT THE OCCUPATION AND COUPLING ARRAYS C IF (NBUG6-1) 30,20,30 20 CONTINUE WRITE (IWRITE,3000) (NJ(I),LJ(I),I=1,IHSH) WRITE (IWRITE,3010) (NOSH1(J),J=1,IHSH) WRITE (IWRITE,3010) (NOSH2(J),J=1,IHSH) WRITE (IWRITE,3020) ((J1QN1(J,K),K=1,3),J=1,I2HSH) WRITE (IWRITE,3020) ((J1QN2(J,K),K=1,3),J=1,I2HSH) C C TEST FOR AT MOST ONE ELECTRON DIFFERENCE IN CONFIGURATIONS C 30 CONTINUE NOSHUM = 0 DO 40 K = 1,IHSH NOSHUM = NOSHUM + ABS(NOSH1(K)-NOSH2(K)) 40 CONTINUE IF (NOSHUM-2) 50,50,1020 C C TEST FOR TRIANGLE RELATION BETWEEN KA AND TOTAL ANGULAR MOMENTA C 50 CONTINUE IF (ISPIN.EQ.0) GOTO 60 K = 3 IF (ISPIN.EQ.2) GOTO 70 IF (J1QN1(I2HSH,2).NE.J1QN2(I2HSH,2)) GOTO 1020 GOTO 70 C 60 CONTINUE K = 2 IF (J1QN1(I2HSH,3).NE.J1QN2(I2HSH,3)) GOTO 1020 70 CONTINUE LB = J1QN1(I2HSH,K) - 1 NB = J1QN2(I2HSH,K) - 1 MB = KA + KA BTST = TRITST(MB,LB,NB) IF (BTST) 1020,80,1020 C C DETERMINE IRHO AND ISIGMA, THE NUMBERS OF THE OCCUPIED SHELLS C 80 CONTINUE IRHO = 0 ISIG = 0 DO 110 J = 1,IHSH NX = NOSH1(J) - NOSH2(J) + 2 GOTO (90,110,100),NX C 90 CONTINUE ISIG = J GOTO 110 C 100 CONTINUE IRHO = J 110 CONTINUE IF (IRHO.NE.0) GOTO 120 IRHO = 1 ISIG = 1 120 CONTINUE C C THE BEGINNING OF THE LOOP OVER ALL SHELLS C 130 CONTINUE IF (IRHO.NE.ISIG) GOTO 150 IF (NBUG6-1) 150,140,150 140 CONTINUE WRITE (IWRITE,3140) IRHO 150 CONTINUE NTOT = NTOT + 1 L1 = LJ(IRHO) + 1 L2 = LJ(ISIG) + 1 AJF = A DBLE(J1QN1(I2HSH,2)) / B DBLE(2*LJ(IRHO)+1) IF (ISPIN.EQ.1) AJF = DBLE(J1QN1(I2HSH,3))*HALF IF (ISPIN.EQ.2) AJF = AJF*DBLE(J1QN1(I2HSH,3))*HALF IF (IRHO-ISIG) 240,160,240 C C FIND THE PARENT TERMS GIVEN BY ALLOWED J VALUES IN NTAB WITH IRHO C 160 CONTINUE NELCTS = NOSH1(IRHO) K1 = NTAB1(NELCTS,L1) KK1 = ITAB(K1) DO 230 JJ1 = 1,KK1 IJK1 = 3* (JJ1-1) + JTAB(K1) DO 220 K = 2,3 IJKK = IJK1 + K IF (K.EQ.3) GOTO 170 LA = NTAB(IJKK) MA = 2*LJ(IRHO) + 1 NA = J1QN1(IRHO,K) GOTO 180 C 170 CONTINUE LA = NTAB(IJKK) - 1 MA = 1 NA = J1QN1(IRHO,K) - 1 180 CONTINUE ATST = TRITST(LA,MA,NA) IF (ATST) 200,190,200 190 CONTINUE IF (K-3) 220,210,220 200 CONTINUE JMEM(JJ1) = 0 GOTO 230 C 210 CONTINUE JMEM(JJ1) = 1 220 CONTINUE 230 CONTINUE C C PARENTAGE CHECK C 240 CONTINUE IF (IRHO-ISIG) 250,330,250 250 CONTINUE NELCTS = NOSH1(IRHO) K1 = NTAB1(NELCTS,L1) NELCTS = NOSH2(ISIG) K2 = NTAB1(NELCTS,L2) KK1 = ITAB(K1) KK2 = ITAB(K2) DO 270 JJ1 = 1,KK1 IJK1 = 3* (JJ1-1) + JTAB(K1) DO 260 K = 2,3 IJKK = IJK1 + K MSAM1 = NTAB(IJKK) - J1QN2(IRHO,K) IF (MSAM1.NE.0) GOTO 270 IF (K.EQ.3) GOTO 290 260 CONTINUE 270 CONTINUE IF (NBUG6-1) 1120,280,1120 280 CONTINUE WRITE (IWRITE,3040) GOTO 1120 C 290 CONTINUE DO 310 JJ2 = 1,KK2 IJK2 = 3* (JJ2-1) + JTAB(K2) DO 300 K = 2,3 IJKK = IJK2 + K MSAM2 = NTAB(IJKK) - J1QN1(ISIG,K) IF (MSAM2.NE.0) GOTO 310 IF (K.EQ.3) GOTO 330 300 CONTINUE 310 CONTINUE IF (NBUG6-1) 1120,320,1120 320 CONTINUE WRITE (IWRITE,3040) GOTO 1120 C C SET J2 AND J3 . SAME FOR L AND S C 330 CONTINUE M1 = IHSH - 2 M2 = 2*M1 + 1 M3 = 3*IHSH - 1 M4 = M3 + 1 M5 = M3 + 2 M10 = M5 + 1 NJ1S = M10 + 1 J2(1,1) = M10 J2(1,2) = NJ1S J2(1,3) = M5 J2(2,1) = IRHO J2(2,2) = M5 J2(2,3) = M3 J3(1,1) = ISIG J3(1,2) = M10 J3(1,3) = M4 IF (IRHO-1) 350,340,350 340 CONTINUE J2(3,1) = M3 GOTO 360 C 350 CONTINUE J2(3,1) = 1 360 CONTINUE IF (IRHO-2) 380,370,380 370 CONTINUE J2(3,2) = M3 GOTO 390 C 380 CONTINUE J2(3,2) = 2 390 CONTINUE J2(3,3) = IHSHP1 IF (ISIG-1) 410,400,410 400 CONTINUE J3(2,1) = M4 GOTO 420 C 410 CONTINUE J3(2,1) = 1 420 CONTINUE IF (ISIG-2) 440,430,440 430 CONTINUE J3(2,2) = M4 GOTO 450 C 440 CONTINUE J3(2,2) = 2 450 CONTINUE J3(2,3) = 2*IHSH IF (IHSH-3) 540,460,460 460 CONTINUE DO 530 J = 4,IHSHP1 L = J - 1 J2(J,1) = M1 + L J2(J,3) = M1 + J J3(L,1) = M2 + L J3(L,3) = M2 + J IF (IRHO-L) 490,480,490 480 CONTINUE J2(J,2) = M3 GOTO 500 C 490 CONTINUE J2(J,2) = L 500 CONTINUE IF (ISIG-L) 520,510,520 510 CONTINUE J3(L,2) = M4 GOTO 530 C 520 CONTINUE J3(L,2) = L 530 CONTINUE 540 CONTINUE M6 = IHSHP1 J3(M6,1) = M3 - 1 J3(M6,2) = NJ1S J3(M6,3) = I2HSH IF (IHSH-1) 560,550,560 550 CONTINUE J3(M6,1) = M4 J3(M6,3) = M3 560 CONTINUE DO 580 J = 1,IHSHP1 DO 570 K = 1,3 J2STO(J,K) = J2(J,K) J3STO(J,K) = J3(J,K) 570 CONTINUE 580 CONTINUE C C RECOUPLING COEFFICIENTS C JMEM1 = J1QN1(IRHO,1) JMEM2 = J1QN1(IRHO,2) JMEM3 = J1QN1(IRHO,3) JMEM4 = J1QN2(ISIG,1) JMEM5 = J1QN2(ISIG,2) JMEM6 = J1QN2(ISIG,3) IF (IRHO-ISIG) 650,590,650 C C BEGINNING OF LOOP OVER ALL PARENT TERMS C 590 CONTINUE JJ1 = 1 600 CONTINUE IF (NBUG6-1) 620,610,620 610 CONTINUE WRITE (IWRITE,3120) JJ1 620 CONTINUE IF (JMEM(JJ1).EQ.1) GOTO 640 IF (NBUG6-1) 1050,630,1050 630 CONTINUE WRITE (IWRITE,3100) GOTO 1050 C 640 CONTINUE IJK1 = 3* (JJ1-1) + JTAB(K1) NI1 = NTAB(IJK1+1) NI2 = NTAB(IJK1+2) NI3 = NTAB(IJK1+3) J1QN2(IRHO,1) = NI1 J1QN1(ISIG,1) = NI1 J1QN2(IRHO,2) = NI2 J1QN1(ISIG,2) = NI2 J1QN2(IRHO,3) = NI3 J1QN1(ISIG,3) = NI3 650 CONTINUE K = 2 M7 = M3 - IHSH M9 = M7 + 1 M11 = M3 - 1 M12 = IHSH - 1 RECUPS = ONE NJ23S = M6 + 1 C C SET UP THE J1 ARRAY FOR THE ANGULAR AND SPIN RECOUPLING C COEFFICIENTS C 660 CONTINUE DO 670 J = 1,IHSH IF (J.EQ.IRHO .OR. J.EQ.ISIG) GOTO 670 DO 665 KK = 1,3 IF (J1QN1(J,KK).NE.J1QN2(J,KK)) GOTO 680 665 CONTINUE 670 CONTINUE GOTO 690 C 680 CONTINUE RECUPS = ZERO IF (IRHO.NE.ISIG) GOTO 1020 GOTO 1040 C 690 CONTINUE IF (K-3) 700,710,700 700 CONTINUE J1(M5) = 2*LJ(IRHO) + 1 J1(M10) = 2*LJ(ISIG) + 1 J1(NJ1S) = 2*KA + 1 IF (ISPIN.EQ.1) J1(NJ1S) = 1 J1(M3) = JMEM2 J1(M4) = JMEM5 IF (IRHO.EQ.ISIG) GOTO 720 J1(M3) = J1QN1(IRHO,K) J1(M4) = J1QN2(ISIG,K) GOTO 720 C 710 CONTINUE J1(M5) = 2 J1(M10) = 2 J1(NJ1S) = 1 IF (ISPIN.NE.0) J1(NJ1S) = 2*KA + 1 J1(M3) = JMEM3 J1(M4) = JMEM6 IF (IRHO.EQ.ISIG) GOTO 720 J1(M3) = J1QN1(IRHO,K) J1(M4) = J1QN2(ISIG,K) 720 CONTINUE DO 750 J = 1,IHSH IF (IRHO-J) 740,730,740 730 CONTINUE J1(J) = J1QN2(IRHO,K) GOTO 750 C 740 CONTINUE J1(J) = J1QN1(J,K) 750 CONTINUE IF (IHSH.EQ.1) GOTO 780 DO 760 J = M6,M7 J1(J) = J1QN1(J,K) 760 CONTINUE DO 770 J = M9,M11 JM12 = J - M12 J1(J) = J1QN2(JM12,K) 770 CONTINUE C C PRINT OUT THE J1,J2 AND J3 ARRAYS C C 10 STATEMENTS REPLACED BY THE FOLLOWING 4 -- WE'90MAR15: 780 CONTINUE IF (K.GE.3) GOTO 790 IF (NBUG6.NE.1) GOTO 790 WRITE (IWRITE,3050) (J1(J),J=1,NJ1S) WRITE (IWRITE,3060) ((J2(I,J),J=1,3), (J3(I,J),J=1,3),I=1,IHSHP1) C C EVALUATE ORBITAL AND SPIN RECOUPLING COEFFICIENTS C NJGRAF IS CALLED EACH TIME SO SET THE ELEMENTS OF THE FREE ARRAY C TO .FALSE. C 790 CONTINUE DO 800 I = 1,NJ1S FREE(I) = .FALSE. 800 CONTINUE C CALL NJGRAF(RECUP,FAIL) C RECUPS = RECUPS*RECUP IF (K-3) 810,840,840 810 CONTINUE IF (NBUG6-1) 830,820,830 820 CONTINUE WRITE (IWRITE,3070) RECUP 830 CONTINUE 840 CONTINUE K = K + 1 DO 860 J = 1,IHSHP1 DO 850 KK = 1,3 J2(J,KK) = J2STO(J,KK) J3(J,KK) = J3STO(J,KK) 850 CONTINUE 860 CONTINUE IF (K.EQ.3) GOTO 660 IF (NBUG6-1) 880,870,880 870 CONTINUE WRITE (IWRITE,3080) RECUP C C FIRST FRACTIONAL PARENTAGE COEFFICIENT C 880 CONTINUE LIJ = LJ(IRHO) COEFP = ONE IF (LIJ) 890,900,890 890 CONTINUE N = NOSH1(IRHO) IV1 = JMEM1 IL1 = (JMEM2-1)/2 IS1 = JMEM3 IV2 = J1QN2(IRHO,1) IL2 = (J1QN2(IRHO,2)-1)/2 IS2 = J1QN2(IRHO,3) CALL CFP(LIJ,N,IV1,IL1,IS1,IV2,IL2,IS2,COEFP) RECUPS = RECUPS*COEFP 900 CONTINUE IF (IRHO-ISIG) 910,920,910 910 CONTINUE IF (ABS(RECUPS).LT.1.0D-5) GOTO 1020 C C SECOND FRACTIONAL PARENTAGE COEFFICIENT C 920 CONTINUE LIJ = LJ(ISIG) COEFP = ONE IF (LIJ) 940,940,930 930 CONTINUE N = NOSH2(ISIG) IV1 = JMEM4 IL1 = (JMEM5-1)/2 IS1 = JMEM6 IV2 = J1QN1(ISIG,1) IL2 = (J1QN1(ISIG,2)-1)/2 IS2 = J1QN1(ISIG,3) CALL CFP(LIJ,N,IV1,IL1,IS1,IV2,IL2,IS2,COEFP) 940 CONTINUE RECUPS = RECUPS*COEFP IF (ABS(RECUPS).LT.1.0D-8 .AND. IRHO.NE.ISIG) GOTO 1020 C C PERMUTATION FACTOR C IDELP = 2 IF (IRHO-ISIG) 960,1000,980 960 CONTINUE JRHO = IRHO + 1 DO 970 J = JRHO,ISIG IDELP = IDELP + NOSH1(J) 970 CONTINUE GOTO 1000 C 980 CONTINUE JSIG = ISIG + 1 DO 990 J = JSIG,IRHO IDELP = IDELP + NOSH2(J) 990 CONTINUE 1000 CONTINUE MINUS = (-1)**IDELP C C MULTIPLICATIVE FACTOR C IF (IRHO-ISIG) 1010,1040,1010 1010 CONTINUE SQRN = SQRT(DBLE(NOSH1(IRHO)*NOSH2(ISIG))) VALML = SQRN*RECUPS*DBLE(MINUS) GOTO 1030 C 1020 CONTINUE VALML = ZERO 1030 CONTINUE RML = RML + VALML C RESULT STORED IN VSHELL IF (NTOT.EQ.0) NTOT = 1 VSHELL(NTOT) = RML*SQRT(AJF) GOTO 1100 C 1040 CONTINUE VALUML = RECUPS IF (NBUG6.NE.0) WRITE (IWRITE,3130) JJ1,VALUML RPL = RPL + VALUML 1050 CONTINUE IF (IRHO.NE.ISIG) GOTO 1060 JJ1 = JJ1 + 1 IF (JJ1.LE.KK1) GOTO 600 1060 CONTINUE J1QN1(IRHO,1) = JMEM1 J1QN1(IRHO,2) = JMEM2 J1QN1(IRHO,3) = JMEM3 J1QN2(ISIG,1) = JMEM4 J1QN2(ISIG,2) = JMEM5 J1QN2(ISIG,3) = JMEM6 ANL = DBLE(NOSH1(IRHO))*RPL C C RESULTS STORED IN VSHELL C IF (NTOT.EQ.0) NTOT = 1 VSHELL(NTOT) = ANL*SQRT(AJF) IF (NBUG6-1) 1090,1080,1090 1080 CONTINUE WRITE (IWRITE,3090) IRHO,ANL 1090 CONTINUE IRHO = IRHO + 1 ISIG = ISIG + 1 RPL = ZERO IF (IRHO-IHSH) 130,130,1100 1100 CONTINUE IF (NBUG6-1) 1120,1110,1120 1110 CONTINUE WRITE (IWRITE,3110) (VSHELL(N),N=1,NTOT) 1120 CONTINUE C 3000 FORMAT (//' NJ,LJ ', (I6,I3)) 3010 FORMAT (//' NOSH ',10I4) 3020 FORMAT (//' J1QN ',30I3) 3040 FORMAT (/'0PARENT TERMS NOT FOUND'//) 3050 FORMAT (/'0J1'/ (24I5)) 3060 FORMAT (' J2',19X,'J3'/ (3I5,I10,2I5)) 3070 FORMAT (///' ORBITAL RECOUPLING COEFF=',E20.8) 3080 FORMAT (///' SPIN RECOUPLING COEFF=',E20.8//) 3090 FORMAT (/' THE CONTRIBUTION FROM SHELL',I2,' IS',F15.8) 3100 FORMAT (//' THIS IS NOT A PARENT') 3110 FORMAT (///' VSHELL=',8F15.8) 3120 FORMAT (//' FRACTIONAL PARENT TERMS',I2) 3130 FORMAT (//' THE CONTRIBUTION FROM FRACTIONAL PARENTAGE TERMS',I2, A ' IS',F15.8) 3140 FORMAT (//' SHELL',I2) END C C C SUBROUTINE TRDEL(JJ1,JJ2,JJ3,NBN,FAIL) C C TEST FOR TRIANGULAR DELTA.IF NOT SATISFIED FAIL=.TRUE. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP LOGICAL CUT COMMON /CUTDIG/CUT LOGICAL FREE COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) C IF (SUMVAR(JJ1) .OR. SUMVAR(JJ2) .OR. SUMVAR(JJ3)) GOTO 10 IF (NBN.GT.4) CUT = .TRUE. IF (FREE(JJ1) .OR. FREE(JJ2) .OR. FREE(JJ3)) GOTO 10 I1 = J1(JJ1) I2 = J1(JJ2) I3 = J1(JJ3) IF (I1.GE. (ABS(I2-I3)+1) .AND. I1.LE. (I2+I3-1)) GOTO 10 FAIL = .TRUE. 10 CONTINUE C END C C C SUBROUTINE TRIANG(FAIL) C C REDUCES A TRIANGLE HAVING ONE APEX AT EITHER END OF THE AXIS OF C THE FLAT DIAGRAM. C THIS INTRODUCES ONE 6J SYMBOL AND SOME PHASE FACTORS . C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFL2C=2*KFL2+2) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL C INTEGER ARR,TAB1 COMMON /GRAPH/JDIAG(KFL2B,3),ARR(KFL2B,3),TAB1(KFL1,2),IL(KFL2B), A IH(KFL2B),NPOINT(KFL2C),NBNODE,IFIRST,ILAST,IPARTS,IPARTL, B NPART,ICROSS,NFREE,ITFREE(KFL2A),NFIN,NC LOGICAL SUMVAR COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP C CHARACTER*6 NAME,NAMSUB COMMON /NAM/NAMSUB DATA NAME/'TRIANG'/ C NAMSUB = NAME IT1 = NPOINT(1) IT2 = NPOINT(2) IT3 = NPOINT(3) JWC = JWC + 1 JW(1,JWC) = JDIAG(IT3,2) JW(2,JWC) = JDIAG(IT2,3) JW(3,JWC) = JDIAG(IT3,1) IF (ARR(IT3,1).GT.0) THEN C C PHASE2(JW(3,JWC)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JW(3,JWC) C ENDIF C JW(4,JWC) = JDIAG(IT2,1) IF (ARR(IT2,1).LT.0) THEN C C PHASE2(JW(4,JWC)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JW(4,JWC) C ENDIF C K23 = 3 IF (IT1.EQ.IFIRST) K23 = 2 JW(5,JWC) = JDIAG(IT1,K23) JW(6,JWC) = JDIAG(IT3,3) CALL TRDEL(JW(1,JWC),JW(2,JWC),JW(5,JWC),NBNODE,FAIL) IF (FAIL) GOTO 60 IF (ARR(IT3,3).GT.0) THEN C C PHASE2(JW(6,JWC)): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JW(6,JWC) C ENDIF C JT1 = JW(5,JWC) JDIAG(IT3,1) = JT1 JDIAG(IT3,3) = JW(2,JWC) ARR(IT3,1) = ARR(IT1,K23) ARR(IT3,3) = ARR(IT2,3) IF (IT1.EQ.IFIRST) GOTO 10 TAB1(JT1,1) = IT3 TAB1(JT1,2) = IH(NBNODE-1) K12 = 1 GOTO 20 C 10 CONTINUE TAB1(JT1,1) = IH(2) TAB1(JT1,2) = IT3 K12 = 2 20 CONTINUE IL3 = IL(IT3) IF (IT1.EQ.ILAST) GOTO 40 IL2 = IL(IT2) - 1 DO 30 I = 2,IL2 IT = IH(I) ILP = I - 1 IL(IT) = ILP IH(ILP) = IT 30 CONTINUE 40 CONTINUE DO 50 I = IL3,NBNODE IT = IH(I) ILP = I - K12 IL(IT) = ILP IH(ILP) = IT 50 CONTINUE 60 CONTINUE C END C C C FUNCTION TRITST(L,M,N) C*********************************************************************** C C TRITST C C*********************************************************************** C C IF TRITST=1.0 THE TRIANGLE RELATION IS NOT SATISFIED C IF TRITST=0.0 THE TRIANGLE RELATION IS SATISFIED C IMPLICIT REAL*8 (O,T,Z) PARAMETER (ZERO=0.0D0,ONE=1.0D0) C LMN = ABS(L-M) C IF (N.GE.LMN .AND. N.LE.L+M) THEN TRITST = ZERO ELSE TRITST = ONE ENDIF C END C C C SUBROUTINE VAR(JN,JNS,JNC,JNSC,JBC,SUMVAR,MP,M,INVER) C C TEST FOR VARIABLE CHARACTER AND PUT IN JNS IF YES,AND JN NOW C CONTAINS 0. C PARAMETER (KFLV=40) C LOGICAL SUMVAR(MP) DIMENSION JN(JNC),JNS(KFLV),INVER(MP) C JNSC = 0 IF (JBC.EQ.JNC) GOTO 20 JBBC = JBC + 1 DO 10 I = JBBC,JNC I1 = JN(I) IF (.NOT.SUMVAR(I1)) GOTO 10 JNSC = JNSC + 1 IF (JNSC.GT.KFLV) THEN PRINT 3000,JNSC,KFLV STOP C ENDIF C J = INVER(I1) JNS(JNSC) = J JN(I) = M 10 CONTINUE 20 CONTINUE C 3000 FORMAT (' DIMENSION ERROR IN VAR. JNSC=',I5,' KFLV=',I5) END C C C SUBROUTINE VECTOR(N,EIG,X,NO,P,R,RN,XM, XCH,B,U,V,W) C C TAKES ARRAYS R OF MAIN DIAGONAL ELEMENTS, P OF SUPER DIAGONAL C ELEMENTS, EIG OF EIGENVALUES, OF THE TRI-DIAGONAL MATRIX C AND BY MEANS OF INVERSE ITERATIONS DETERMINES C AN EIGENVECTOR OF THE TRI-DIAGONAL MATRIX. C IMPLICIT REAL*8 (A-H,O-Z) LOGICAL ITER DIMENSION EIG(N),X(N),P(N),R(N),RN(N),XM(N), XCH(N),B(N),U(N), A V(N),W(N) DATA EPSI/1.0D-16/ C P(N)=0.0D0 !INITIALIZE C C THE ARRAY RN HOLDS THE MAIN DIAGONAL ELEMENTS OF A NEW C TRI-DIAGONAL MATRIX. C DO 10 K = 1,N RN(K) = R(K) - EIG(NO) 10 CONTINUE C C BY MEANS OF GAUSSIAN ELIMINATION THE NEW TRI-DIAGONAL MATRIX C IS TRANSFORMED INTO UPPER TRIANGULAR FORM. THE ROW MULTIPLIERS C ARE STORED IN ARRAY XM. IF ROW I IS INTERCHANGED WITH ROW I+1 C LXCH(I)=1 THE MAIN DIAGONAL ELEMENTS OF THE UPPER TRIANGULAR C MATRIX ARE STORED IN THE ARRAY U, THE NEXT DIAGONAL IN THE C ARRAY V AND THE LAST DIAGONAL IN THE ARRAY W. C C LXCH TYPE INTEGER INCOMPATIBLE WITH CALLING BY HSLDR. TRY XCH - NRB. C PA = RN(1) QA = P(1) N1 = N - 1 DO 50 I = 1,N1 C C DETERMINE IF A ROW INTERCHANGE IS NECESSARY. C GA = ABS(P(I)) PPA = ABS(PA) IF (GA.LE.EPSI) GOTO 20 IF (GA.GT.PPA) GOTO 30 C C NO INTERCHANGE. C 20 CONTINUE U(I) = PA V(I) = QA W(I) = 0.0D0 HA = P(I) PA = RN(I+1) QA = P(I+1) XCH(I) = 0.0D0 GOTO 40 C C INTERCHANGE. C 30 CONTINUE U(I) = P(I) V(I) = RN(I+1) W(I) = P(I+1) HA = PA PA = QA QA = 0.0D0 XCH(I) = 1.0D0 C C THE ROW MULTIPLIER IS DETERMINED. C 40 CONTINUE XM(I) = HA/U(I) C C ROW I IS MULTIPLIED BY XM(I) AND SUBTRACTED FROM ROW I+1. C HA = 0.0D0 PA = PA - XM(I)*V(I) QA = QA - XM(I)*W(I) 50 CONTINUE C C THE SINGLE ELEMENT IN THE LAST ROW IS PLACED IN U(N) C U(N) = PA IF (ABS(U(N)).LT.EPSI) U(N) = EPSI C C THE ARRAY X IS RESERVED FOR THE CURRENT ESTIMATE OF THE C EIGENVECTOR. THE SUBROUTINE BAKSUB IS CALLED AND AN ESTIMATE OF C THE EIGENVECTOR OBTAINED. C DO 70 I = 1,N B(I) = 1.0D0 70 CONTINUE CALL BAKSUB(N,X,B,U,V,W) C C THE EIGENVECTOR STORED IN X IS NORMALISED. C CALL NORM(N,X) C C TO SAVE COMPUTING TIME, A FURTHER ITERATION TO YIELD A BETTER C ESTIMATE OF THE VECTOR CAN BE SUPPRESSED BY SETTING ITER=.TRUE. C ITER = .FALSE. IF (ITER) GOTO 130 C C A NEW COLUMN VECTOR RELATED TO X BY EXACTLY THE SAME ROW C INTERCHANGES AND MULTIPLICATIONS BY WHICH THE UPPER TRIANGULAR C MATRIX WAS OBTAINED FROM THE NEW TRI-DIAGONAL MATRIX, IS C CALCULATED AND STORED IN B. C DO 100 I = 1,N B(I) = X(I) 100 CONTINUE N1 = N - 1 DO 120 I = 1,N1 IF ( XCH(I).LT.0.5D0) GOTO 110 AC = B(I) BC = B(I+1) B(I) = BC B(I+1) = AC 110 B(I+1) = B(I+1) - XM(I)*B(I) 120 CONTINUE C C THE NEW VECTOR STORED IN B IS NORMALISED. C CALL NORM(N,B) C C THE SUBROUTINE BAKSUB IS CALLED AND A NEW ESTIMATE OF THE C EIGENVECTOR IS OBTAINED. C CALL BAKSUB(N,X,B,U,V,W) C C THE NEW ESTIMATE OF THE EIGENVECTOR IS NORMALISED. C CALL NORM(N,X) 130 CONTINUE C END C C C SUBROUTINE WAY(L,KA,KB,ICH,NB) C C C TESTS ONE STEP FORWARD IF THE WAY IS FREE.FIRST AND SECOND C ARGUMENTS ARE INTERCHANGED OR NOT ACCORDING TO ICH=-1,OR +1 C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4,KFL2A=2*KFL2,KFL2B=4*KFL2) C LOGICAL TABS INTEGER ARROW COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR C COMMON /BUILD/IAL(KFL2B),IF1,IF2,NODE C K1 = J23(L,KA) K2 = J23(L,KB) NB = IAL(K1) + IAL(K2) - 1 IF (NB) 20,10,60 10 CONTINUE NB1 = IAL(K1) - IAL(K2) IF (NB1) 70,60,60 C C OTHERJ(L,K1,L1,LC1,LA) C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C 20 CONTINUE L1 = LINE(K1,1) IF (L1.EQ.L .OR. TABS(L1)) THEN L1 = LINE(K1,2) LC1 = LCOL(K1,2) ELSE LC1 = LCOL(K1,1) ENDIF C C OTHERJ(L,K2,L2,LC2,LB): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C L2 = LINE(K2,1) IF (L2.EQ.L .OR. TABS(L2)) THEN L2 = LINE(K2,2) LC2 = LCOL(K2,2) ELSE LC2 = LCOL(K2,1) ENDIF C C C NEIBOR(LC1,I1,I2): C GIVES THE POSITIONS OF THE OTHER TWO ARGUMENTS IN THE TRIAD. C IF (LC1.LT.2) THEN I1 = 2 I2 = 3 C ELSE IF (LC1.EQ.2) THEN I1 = 3 I2 = 1 C ELSE I1 = 1 I2 = 2 ENDIF C C NEIBOR(LC2,I3,I4): C GIVES THE POSITIONS OF THE OTHER TWO ARGUMENTS IN THE TRIAD. C IF (LC2.LT.2) THEN I3 = 2 I4 = 3 C ELSE IF (LC2.EQ.2) THEN I3 = 3 I4 = 1 C ELSE I3 = 1 I4 = 2 ENDIF C JI1 = J23(L1,I1) JI2 = J23(L1,I2) JI3 = J23(L2,I3) JI4 = J23(L2,I4) IA = IAL(JI1) + IAL(JI2) IB = IAL(JI3) + IAL(JI4) NBP = IB + IA + 1 NBM = IB - IA GOTO (60,30,40,30,50),NBP C 30 CONTINUE IF (NBM) 70,60,60 40 CONTINUE IF (NBM) 70,50,60 50 CONTINUE IF (JI3.EQ.IF1 .OR. JI3.EQ.IF2 .OR. JI4.EQ.IF1 .OR. A JI4.EQ.IF2) GOTO 70 60 CONTINUE ICH = 1 GOTO 80 C 70 CONTINUE ICH = -1 80 CONTINUE C END C C C SUBROUTINE ZERO(J,JZ,FAIL) C C SUPPRESSES ONE LINE AND TWO NODES OF THE UNSTRUCTURED GRAPH C INTRODUCES ZEROS IN THE TRIADS J23.AS A CONSEQUENCE THE OTHER C TWO ARGUMENTS OF THE TRIAD ARE PUT EQUAL.IF THERE WAS ALREADY C A ZERO IN THE TRIAD WHICH IS CHANGED IT IS A SPECIAL CASE. C INCLUDE 'PARAM' C PARAMETER (MXORB=(MZLR1*(2*MZNR1-MZLR1+1))/2) PARAMETER (KFL1=100,KFL2=MXORB+4) PARAMETER (KFL2A=2*KFL2,KFL2B=4*KFL2,KFLZ=99) PARAMETER (KFL6=120,KFL7=150,KFL8=120,KFL9=40,KFLW=20) C LOGICAL FAIL COMMON /DEBUG/IBUG1,IBUG2,IBUG3,IBUG4,IBUG5,IBUG6,IBUG7,IBUG8, A IBUG9 COMMON /DIM/J6CC,J7CC,J8CC,J9CC,JWCC,JDELC COMMON /ZER/NZERO,JZERO(KFLZ) COMMON /BUILD/IAL(KFL2B),IF1,IF2,NODE C LOGICAL TABS,FREE,SUMVAR,CUT,NOCUT INTEGER ARROW COMMON /TREE/J23(KFL2A,3),ARROW(KFL2A,3),LINE(KFL1,2), A LCOL(KFL1,2),TABS(KFL2A),NBTR COMMON /COUPLE/M,N,J1(KFL1),J2(KFL2,3),J3(KFL2,3),FREE(KFL1) COMMON /ARGU/J6C,J7C,J8C,J9C,JWC,J6(KFL6),J7(KFL7),J8(KFL8), A J9(KFL9),JW(6,KFLW),JDEL,LDEL(KFLW,2),SUMVAR(KFL1),MP COMMON /CUTDIG/CUT C CHARACTER*6 NAME DATA NAME/'ZERO '/ C NOCUT = .FALSE. NZERO = 0 IF (J.LT.1) GOTO 10 C C OTHERJ(0,JZ,LIN,LC,K1): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C LIN = LINE(JZ,1) IF (LIN.EQ.0 .OR. TABS(LIN)) THEN K1 = 1 LIN = LINE(JZ,2) LC = LCOL(JZ,2) C ELSE K1 = 2 LC = LCOL(JZ,1) ENDIF C I = NZERO GOTO 50 C 10 CONTINUE DO 20 I = 1,M IF (J1(I).NE.1 .OR. FREE(I) .OR. IAL(I).LE.1) GOTO 20 NZERO = NZERO + 1 IF (NZERO.GT.KFLZ) THEN PRINT 3000,NZERO,KFLZ STOP C ENDIF C JZERO(NZERO) = I 20 CONTINUE NOCUT = .TRUE. M = M + 1 J1(M) = 1 SUMVAR(M) = .FALSE. FREE(M) = .FALSE. IF (NZERO.EQ.0) GOTO 160 IF (IBUG3.EQ.1) CALL PRINTJ(NAME,1) I = 0 30 CONTINUE I = I + 1 JZ = JZERO(I) J = 0 40 CONTINUE J = J + 1 LIN = LINE(JZ,J) IF (TABS(LIN)) GOTO 150 LC = LCOL(JZ,J) C C NEIBOR(LC,L1,L2): C GIVES THE POSITIONS OF THE OTHER TWO ARGUMENTS IN THE TRIAD. C 50 CONTINUE IF (LC.LT.2) THEN L1 = 2 L2 = 3 C ELSE IF (LC.EQ.2) THEN L1 = 3 L2 = 1 C ELSE L1 = 1 L2 = 2 ENDIF C JJ1 = J23(LIN,L1) JJ2 = J23(LIN,L2) IF (JJ1.EQ.JJ2) THEN J6C = J6C + 1 J6(J6C) = JJ1 GOTO 110 ENDIF C C CALL DELTA(JJ1,JJ2,FAIL) - INLINED C C TEST FOR DELTA(JJ1,JJ2).IF THEY ARE SUMMATION VARIABLES,THE SECOND C IS CHANGED INTO THE FIRST EVERYWHERE.IF THEY ARE FIXED,THEIR C VALUE IS CHECKED,AND FAIL PUT TO .TRUE. IF THEY DIFFER. C IF (IBUG3.EQ.1) PRINT 3001,JJ1,SUMVAR(JJ1),JJ2,SUMVAR(JJ2) IF (SUMVAR(JJ1) .AND. SUMVAR(JJ2)) GOTO 21 IF (FREE(JJ1) .OR. FREE(JJ2)) THEN JDEL = JDEL + 1 LDEL(JDEL,1) = JJ1 LDEL(JDEL,2) = JJ2 SUMVAR(JJ1) = .FALSE. SUMVAR(JJ2) = .FALSE. GOTO 161 ENDIF IF (J1(JJ1).NE.J1(JJ2)) FAIL = .TRUE. CUT = .TRUE. GOTO 161 C 21 CONTINUE IF (J6C.NE.J6CC) THEN DO 31 II = J6CC + 1,J6C IF (J6(II).EQ.JJ2) J6(II) = JJ1 31 CONTINUE ENDIF IF (J7C.NE.J7CC) THEN DO 51 II = J7CC + 1,J7C IF (J7(II).EQ.JJ2) J7(II) = JJ1 51 CONTINUE ENDIF IF (J8C.NE.J8CC) THEN DO 71 II = J8CC + 1,J8C IF (J8(II).EQ.JJ2) J8(II) = JJ1 71 CONTINUE ENDIF IF (J9C.NE.J9CC) THEN DO 91 II = J9CC + 1,J9C IF (J9(II).EQ.JJ2) J9(II) = JJ1 91 CONTINUE ENDIF IF (JWC.NE.JWCC) THEN DO 121 II = JWCC + 1,JWC IF (JW(1,II).EQ.JJ2) JW(1,II) = JJ1 IF (JW(2,II).EQ.JJ2) JW(2,II) = JJ1 IF (JW(3,II).EQ.JJ2) JW(3,II) = JJ1 IF (JW(4,II).EQ.JJ2) JW(4,II) = JJ1 IF (JW(5,II).EQ.JJ2) JW(5,II) = JJ1 IF (JW(6,II).EQ.JJ2) JW(6,II) = JJ1 121 CONTINUE ENDIF IF (JDEL.NE.JDELC) THEN DO 151 II = JDELC + 1,JDEL IF (LDEL(II,1).EQ.JJ2) LDEL(II,1) = JJ1 IF (LDEL(II,2).EQ.JJ2) LDEL(II,2) = JJ1 151 CONTINUE SUMVAR(JJ2) = .FALSE. ENDIF 161 CONTINUE IF (FAIL) GOTO 160 IF ((J1(JJ1).NE.1 .AND. J1(JJ2).NE.1) * .OR. J1(JJ1).LT.J1(JJ2)) GOTO 100 IF (J1(JJ1).GT.J1(JJ2)) GOTO 90 IF (NZERO.EQ.0) GOTO 100 DO 80 JJX = I,NZERO JJZ = JZERO(JJX) IF (JJ1.EQ.JJZ) GOTO 100 IF (JJ2.EQ.JJZ) GOTO 90 80 CONTINUE GOTO 100 C 90 CONTINUE JJZ = JJ2 JJ2 = JJ1 JJ1 = JJZ C C OTHERJ(LIN,JJ1,LO1,LCO1,K1): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C 100 CONTINUE LO1 = LINE(JJ1,1) IF (LO1.EQ.LIN .OR. TABS(LO1)) THEN K1 = 1 LO1 = LINE(JJ1,2) LCO1 = LCOL(JJ1,2) C ELSE K1 = 2 LCO1 = LCOL(JJ1,1) ENDIF C C OTHERJ(LIN,JJ2,LO2,LCO2,K2): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C LO2 = LINE(JJ2,1) IF (LO2.EQ.LIN .OR. TABS(LO2)) THEN LO2 = LINE(JJ2,2) LCO2 = LCOL(JJ2,2) C ELSE LCO2 = LCOL(JJ2,1) ENDIF C J9C = J9C + 1 J9(J9C) = JJ1 J23(LO2,LCO2) = JJ1 LINE(JJ1,K1) = LO2 LCOL(JJ1,K1) = LCO2 110 CONTINUE IF (ARROW(LIN,L1).LT.ARROW(LIN,L2)) THEN C C PHASE2(JJ1): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JJ1 C ELSE IF (ARROW(LIN,L1).EQ.ARROW(LIN,L2)) THEN ARROW(LO1,LCO1) = 1 ARROW(LO2,LCO2) = -1 ENDIF TABS(LIN) = .TRUE. NBTR = NBTR - 1 IF (NBTR.EQ.0) GOTO 160 IF (LO1.NE.LO2) GOTO 150 L = 6 - LCO1 - LCO2 JT = J23(LO1,L) IF (J1(JT).EQ.1 .AND. .NOT.FREE(JT)) GOTO 150 C CALL DELTA(JT,M,FAIL) - INLINED C C TEST FOR DELTA(JT,M).IF THEY ARE SUMMATION VARIABLES,THE SECOND C IS CHANGED INTO THE FIRST EVERYWHERE.IF THEY ARE FIXED,THEIR C VALUE IS CHECKED,AND FAIL PUT TO .TRUE. IF THEY DIFFER. C IF (IBUG3.EQ.1) PRINT 3001,JT,SUMVAR(JT),M,SUMVAR(M) IF (SUMVAR(JT) .AND. SUMVAR(M)) GOTO 22 IF (FREE(JT) .OR. FREE(M)) THEN JDEL = JDEL + 1 LDEL(JDEL,1) = JT LDEL(JDEL,2) = M SUMVAR(JT) = .FALSE. SUMVAR(M) = .FALSE. GOTO 162 ENDIF IF (J1(JT).NE.J1(M)) FAIL = .TRUE. CUT = .TRUE. GOTO 162 C 22 CONTINUE IF (J6C.NE.J6CC) THEN DO 32 II = J6CC + 1,J6C IF (J6(II).EQ.M) J6(II) = JT 32 CONTINUE ENDIF IF (J7C.NE.J7CC) THEN DO 52 II = J7CC + 1,J7C IF (J7(II).EQ.M) J7(II) = JT 52 CONTINUE ENDIF IF (J8C.NE.J8CC) THEN DO 72 II = J8CC + 1,J8C IF (J8(II).EQ.M) J8(II) = JT 72 CONTINUE ENDIF IF (J9C.NE.J9CC) THEN DO 92 II = J9CC + 1,J9C IF (J9(II).EQ.M) J9(II) = JT 92 CONTINUE ENDIF IF (JWC.NE.JWCC) THEN DO 122 II = JWCC + 1,JWC IF (JW(1,II).EQ.M) JW(1,II) = JT IF (JW(2,II).EQ.M) JW(2,II) = JT IF (JW(3,II).EQ.M) JW(3,II) = JT IF (JW(4,II).EQ.M) JW(4,II) = JT IF (JW(5,II).EQ.M) JW(5,II) = JT IF (JW(6,II).EQ.M) JW(6,II) = JT 122 CONTINUE ENDIF IF (JDEL.NE.JDELC) THEN DO 152 II = JDELC + 1,JDEL IF (LDEL(II,1).EQ.M) LDEL(II,1) = JT IF (LDEL(II,2).EQ.M) LDEL(II,2) = JT 152 CONTINUE SUMVAR(M) = .FALSE. ENDIF 162 CONTINUE IF (FAIL) GOTO 160 C C NEIBOR(L,L1,L2): C GIVES THE POSITIONS OF THE OTHER TWO ARGUMENTS IN THE TRIAD. C IF (L.LT.2) THEN L1 = 2 L2 = 3 C ELSE IF (L.EQ.2) THEN L1 = 3 L2 = 1 C ELSE L1 = 1 L2 = 2 ENDIF C JTF = J23(LO1,L1) IF ((ARROW(LO1,L1)-ARROW(LO1,L2)).LT.0) THEN C C PHASE2(JTF): ADDS A PHASE FACTOR (-1)**2J C J8C = J8C + 1 J8(J8C) = JTF C ENDIF C J6C = J6C + 1 J6(J6C) = JTF NBTR = NBTR - 1 TABS(LO1) = .TRUE. C C OTHERJ(LO1,JT,LIN,LC,K): C GIVES THE OTHER TRIAD WHERE A GIVEN J OCCURS AND ITS POSITION. C LIN = LINE(JT,1) IF (LIN.EQ.LO1 .OR. TABS(LIN)) THEN LIN = LINE(JT,2) LC = LCOL(JT,2) C ELSE LC = LCOL(JT,1) ENDIF C GOTO 50 C 150 CONTINUE IF (J.EQ.1) GOTO 40 IF (NBTR.EQ.0) GOTO 160 IF (I.LT.NZERO) GOTO 30 160 CONTINUE IF (IBUG3.EQ.1) CALL PRINTJ(NAME,4) IF (NOCUT) CUT = .FALSE. C 3000 FORMAT (1X,'DIMENSION ERROR IN ZERO. NZERO=',I5,' KFLZ=',I5) 3001 FORMAT (/' FROM DELTA JA=',I2,L2,5X,'JB=',I2,L2) END